ALL 0.9.3
A Loadbalacing Library
Loading...
Searching...
No Matches
ALL_Staggered.hpp
Go to the documentation of this file.
1/*
2Copyright 2018-2020 Rene Halver, Forschungszentrum Juelich GmbH, Germany
3Copyright 2018-2020 Godehard Sutmann, Forschungszentrum Juelich GmbH, Germany
4
5Redistribution and use in source and binary forms, with or without modification,
6are permitted provided that the following conditions are met:
7
81. Redistributions of source code must retain the above copyright notice, this
9 list of conditions and the following disclaimer.
10
112. Redistributions in binary form must reproduce the above copyright notice,
12this list of conditions and the following disclaimer in the documentation and/or
13 other materials provided with the distribution.
14
153. Neither the name of the copyright holder nor the names of its contributors
16may be used to endorse or promote products derived from this software without
17 specific prior written permission.
18
19THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
20ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
21WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
22DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR
23ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
24(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
25LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
26ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
27(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
28SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
29*/
30
31#ifndef ALL_STAGGERED_HEADER_INCLUDED
32#define ALL_STAGGERED_HEADER_INCLUDED
33
35#include "ALL_LB.hpp"
36#include "ALL_Defines.h"
37#include <exception>
38#include <mpi.h>
39#include <vector>
40
41namespace ALL {
42
52template <class T, class W> class Staggered_LB : public LB<T, W> {
53public:
60 Staggered_LB(int d, W w, T g) : LB<T, W>(d, g) {
61 this->setWork(w);
62
63 // array of MPI communicators for each direction (to collect work
64 // on each plane)
65 communicators.resize(d);
66 for (int _d = 0; _d < d; ++_d)
67 communicators.at(_d) = MPI_COMM_NULL;
68
69 nNeighbors.resize(2 * d);
70 }
71
74 ;
75
77 void setup() override;
78
83 void balance(int step) override;
84
88 virtual std::vector<int> &getNeighbors() override;
89
92 virtual void setAdditionalData(known_unused const void *data) override {}
93
97 virtual W getEstimatedEfficiency() override {return (W)-1;};
98
99private:
101 MPI_Datatype MPIDataTypeT;
103 MPI_Datatype MPIDataTypeW;
104
106 std::vector<MPI_Comm> communicators;
107
109 std::vector<int> neighbors;
111 std::vector<int> nNeighbors;
112
114 void find_neighbors();
115};
116
117template <class T, class W> Staggered_LB<T, W>::~Staggered_LB() {
118 int dimension = this->getDimension();
119 for (int i = 1; i < dimension; ++i) {
120 if (communicators.at(i) != MPI_COMM_NULL)
121 MPI_Comm_free(&communicators.at(i));
122 }
123}
124
125// setup routine for the tensor-based load-balancing scheme
126// requires:
127// this->globalComm (int): cartesian MPI communicator, from
128// which separate sub communicators
129// are derived in order to represent
130// each plane of domains in the system
131template <class T, class W> void Staggered_LB<T, W>::setup() {
132 int status;
133 int dimension = this->getDimension();
134
135 // check if Communicator is cartesian
136 MPI_Topo_test(this->globalComm, &status);
137 if (status != MPI_CART) {
139 __FILE__, __func__, __LINE__,
140 "Cartesian MPI communicator required, passed communicator is not "
141 "cartesian");
142 }
143
144 // get the local coordinates, periodicity and global size from the MPI
145 // communicator
146 MPI_Cart_get(this->globalComm, dimension, this->global_dims.data(),
147 this->periodicity.data(), this->local_coords.data());
148
149 // get the local rank from the MPI communicator
150 MPI_Cart_rank(this->globalComm, this->local_coords.data(), &this->localRank);
151
152 // create sub-communicators
153
154 if (communicators.at(1) != MPI_COMM_NULL)
155 MPI_Comm_free(&communicators.at(1));
156 if (communicators.at(2) != MPI_COMM_NULL)
157 MPI_Comm_free(&communicators.at(2));
158
159 // z-plane
160 MPI_Comm_split(this->globalComm, this->local_coords.at(2),
161 this->local_coords.at(0) +
162 this->local_coords.at(1) * this->global_dims.at(0),
163 &communicators.at(2));
164
165 // y-column
166 MPI_Comm_split(this->globalComm,
167 this->local_coords.at(2) * this->global_dims.at(1) +
168 this->local_coords.at(1),
169 this->local_coords.at(0), &communicators.at(1));
170
171 // only cell itself
172 communicators.at(0) = MPI_COMM_SELF;
173
174 // determine correct MPI data type for template T
175 if (std::is_same<T, double>::value)
176 MPIDataTypeT = MPI_DOUBLE;
177 else if (std::is_same<T, float>::value)
178 MPIDataTypeT = MPI_FLOAT;
179 else if (std::is_same<T, int>::value)
180 MPIDataTypeT = MPI_INT;
181 else if (std::is_same<T, long>::value)
182 MPIDataTypeT = MPI_LONG;
183 else {
185 __FILE__, __func__, __LINE__,
186 "Invalid data type for boundaries given (T)");
187 }
188
189 // determine correct MPI data type for template W
190 if (std::is_same<W, double>::value)
191 MPIDataTypeW = MPI_DOUBLE;
192 else if (std::is_same<W, float>::value)
193 MPIDataTypeW = MPI_FLOAT;
194 else if (std::is_same<W, int>::value)
195 MPIDataTypeW = MPI_INT;
196 else if (std::is_same<W, long>::value)
197 MPIDataTypeW = MPI_LONG;
198 else {
199 throw InvalidCommTypeException(__FILE__, __func__, __LINE__,
200 "Invalid data type for work given (W)");
201 }
202
203 // calculate neighbors
204 int rank_left, rank_right;
205
206 neighbors.clear();
207 for (int i = 0; i < dimension; ++i) {
208 MPI_Cart_shift(this->globalComm, i, 1, &rank_left, &rank_right);
209 neighbors.push_back(rank_left);
210 neighbors.push_back(rank_right);
211 }
212}
213
214template <class T, class W> void Staggered_LB<T, W>::balance(int) {
215 std::vector<Point<T>> newVertices = this->vertices;
216 int dimension = this->getDimension();
217
218 // store original vertices
219 this->prevVertices = this->vertices;
220
221 // loop over all available dimensions
222 for (int i = 0; i < dimension; ++i) {
223 W work_local_plane;
224#ifdef ALL_DEBUG_ENABLED
225 MPI_Barrier(this->globalComm);
226 if (this->localRank == 0)
227 std::cout << "ALL::Staggered_LB::balance(): before work computation..."
228 << std::endl;
229#endif
230 // collect work from all processes in the same plane
231 MPI_Allreduce(this->work.data(), &work_local_plane, 1, MPIDataTypeW,
232 MPI_SUM, communicators[i]);
233
234#ifdef ALL_DEBUG_ENABLED
235 MPI_Barrier(this->globalComm);
236 if (this->localRank == 0)
237 std::cout << "ALL::Staggered_LB::balance(): before work distribution..."
238 << std::endl;
239#endif
240 // correct right border:
241
242 W remote_work;
243 T local_size;
244 T remote_size;
245 // determine neighbors
246 int rank_left, rank_right;
247
248 MPI_Cart_shift(this->globalComm, i, 1, &rank_left, &rank_right);
249
250 // collect work from right neighbor plane
251 MPI_Request sreq, rreq;
252 MPI_Status sstat, rstat;
253
254 MPI_Irecv(&remote_work, 1, MPIDataTypeW, rank_right, 0, this->globalComm,
255 &rreq);
256 MPI_Isend(&work_local_plane, 1, MPIDataTypeW, rank_left, 0,
257 this->globalComm, &sreq);
258 MPI_Wait(&sreq, &sstat);
259 MPI_Wait(&rreq, &rstat);
260
261 // collect size in dimension from right neighbor plane
262
263 local_size = this->prevVertices.at(1)[i] - this->prevVertices.at(0)[i];
264
265 MPI_Irecv(&remote_size, 1, MPIDataTypeT, rank_right, 0, this->globalComm,
266 &rreq);
267 MPI_Isend(&local_size, 1, MPIDataTypeT, rank_left, 0, this->globalComm,
268 &sreq);
269 MPI_Wait(&sreq, &sstat);
270 MPI_Wait(&rreq, &rstat);
271
272 // automatic selection of gamma to guarantee stability of method (choice of
273 // user gamma ignored)
274 this->gamma =
275 std::max(4.1, 2.0 * (1.0 + std::max(local_size, remote_size) /
276 std::min(local_size, remote_size)));
277
278#ifdef ALL_DEBUG_ENABLED
279 MPI_Barrier(this->globalComm);
280 if (this->localRank == 0)
281 std::cout << "ALL::Staggered_LB::balance(): before shift calculation..."
282 << std::endl;
283#endif
284
285 T shift = Functions::borderShift1d(
286 rank_right, this->local_coords.at(i), this->global_dims.at(i),
287 work_local_plane, remote_work, local_size, remote_size, this->gamma,
288 this->minSize[i]);
289
290#ifdef ALL_DEBUG_ENABLED
291 MPI_Barrier(this->globalComm);
292 if (this->localRank == 0)
293 std::cout << "ALL::Staggered_LB::balance(): before shift distibution..."
294 << std::endl;
295#endif
296 // send shift to right neighbors
297 T remote_shift = (T)0;
298
299 MPI_Irecv(&remote_shift, 1, MPIDataTypeT, rank_left, 0, this->globalComm,
300 &rreq);
301 MPI_Isend(&shift, 1, MPIDataTypeT, rank_right, 0, this->globalComm, &sreq);
302 MPI_Wait(&sreq, &sstat);
303 MPI_Wait(&rreq, &rstat);
304
305#ifdef ALL_DEBUG_ENABLED
306 MPI_Barrier(this->globalComm);
307 std::cout << this->localRank << ": shift = " << shift
308 << " remoteShift = " << remote_shift
309 << " vertices: " << this->vertices.at(0)[i] << ", "
310 << this->vertices.at(1)[i] << std::endl;
311#endif
312
313 // for now: test case for simple program
314
315 // if a left neighbor exists: shift left border
316 if (rank_left != MPI_PROC_NULL && this->local_coords[i] != 0)
317 newVertices.at(0)[i] = this->prevVertices.at(0)[i] + remote_shift;
318 else
319 newVertices.at(0)[i] = this->prevVertices.at(0)[i];
320
321 // if a right neighbor exists: shift right border
322 if (rank_right != MPI_PROC_NULL &&
323 this->local_coords[i] != this->global_dims[i] - 1)
324 newVertices.at(1)[i] = this->prevVertices.at(1)[i] + shift;
325 else
326 newVertices.at(1)[i] = this->prevVertices.at(1)[i];
327
328 // check if vertices are crossed and throw exception if something went wrong
329 if (newVertices.at(1)[i] < newVertices.at(0)[i]) {
330 std::cout << "ERROR on process: " << this->localRank << std::endl;
332 __FILE__, __func__, __LINE__,
333 "Lower border of process larger than upper border of process!");
334 }
335
336 this->setVertices(newVertices);
337
338#ifdef ALL_DEBUG_ENABLED
339 MPI_Barrier(this->globalComm);
340 if (this->localRank == 0)
341 std::cout << "ALL::Staggered_LB::balance(): before neighbor search..."
342 << std::endl;
343#endif
344 find_neighbors();
345 }
346}
347
348template <class T, class W> void Staggered_LB<T, W>::find_neighbors() {
349 auto work = this->getWork();
350 auto vertices = this->prevVertices;
351 auto shifted_vertices = this->vertices;
352
353 neighbors.clear();
354 // collect work from right neighbor plane
355 MPI_Request sreq, rreq;
356 MPI_Status sstat, rstat;
357 // array to store neighbor vertices in Y/Z direction (reused)
358 T *vertices_loc = new T[4];
359 T *vertices_rem = new T[8 * this->global_dims[0] * this->global_dims[1]];
360
361 int rem_rank;
362 int rem_coords[3];
363
364 // determine neighbors
365 int rank_left, rank_right;
366
367 // offset to get the correct rank
368 int rank_offset;
369 int offset_coords[3];
370
371 // X-neighbors are static
372 nNeighbors.at(0) = nNeighbors.at(1) = 1;
373
374 // find X-neighbors
375 MPI_Cart_shift(this->globalComm, 0, 1, &rank_left, &rank_right);
376
377 // store X-neighbors
378 neighbors.push_back(rank_left);
379 neighbors.push_back(rank_right);
380
381 // find Y-neighbors to get border information from
382 MPI_Cart_shift(this->globalComm, 1, 1, &rank_left, &rank_right);
383
384 // collect border information from local column
385 vertices_loc[0] = shifted_vertices.at(0)[0];
386 vertices_loc[1] = shifted_vertices.at(1)[0];
387 MPI_Allgather(vertices_loc, 2, MPIDataTypeT,
388 vertices_rem + 2 * this->global_dims[0], 2, MPIDataTypeT,
389 communicators[1]);
390
391 // exchange local column information with upper neighbor in Y direction (cart
392 // grid)
393 MPI_Irecv(vertices_rem, 2 * this->global_dims[0], MPIDataTypeT, rank_left, 0,
394 this->globalComm, &rreq);
395 MPI_Isend(vertices_rem + 2 * this->global_dims[0], 2 * this->global_dims[0],
396 MPIDataTypeT, rank_right, 0, this->globalComm, &sreq);
397
398 // determine the offset in ranks
399 offset_coords[0] = 0;
400 offset_coords[1] = this->local_coords[1] - 1;
401 offset_coords[2] = this->local_coords[2];
402
403 rem_coords[1] = offset_coords[1];
404 rem_coords[2] = offset_coords[2];
405
406 MPI_Cart_rank(this->globalComm, offset_coords, &rank_offset);
407
408 // wait for communication
409 MPI_Wait(&sreq, &sstat);
410 MPI_Wait(&rreq, &rstat);
411
412 // iterate about neighbor borders to determine the neighborship relation
413 nNeighbors.at(2) = 0;
414 for (int x = 0; x < this->global_dims[0]; ++x) {
415 if ((vertices_rem[2 * x] <= vertices_loc[0] &&
416 vertices_loc[0] < vertices_rem[2 * x + 1]) ||
417 (vertices_rem[2 * x] < vertices_loc[1] &&
418 vertices_loc[1] <= vertices_rem[2 * x + 1]) ||
419 (vertices_rem[2 * x] >= vertices_loc[0] &&
420 vertices_loc[0] < vertices_rem[2 * x + 1] &&
421 vertices_loc[1] >= vertices_rem[2 * x + 1])) {
422 nNeighbors.at(2)++;
423 rem_coords[0] = x;
424 MPI_Cart_rank(this->globalComm, rem_coords, &rem_rank);
425 neighbors.push_back(rem_rank);
426 }
427 }
428
429 // barrier to ensure every process concluded the calculations before
430 // overwriting remote borders!
431 MPI_Barrier(this->globalComm);
432
433 // exchange local column information with lower neighbor in Y direction (cart
434 // grid)
435 MPI_Irecv(vertices_rem, 2 * this->global_dims[0], MPIDataTypeT, rank_right, 0,
436 this->globalComm, &rreq);
437 MPI_Isend(vertices_rem + 2 * this->global_dims[0], 2 * this->global_dims[0],
438 MPIDataTypeT, rank_left, 0, this->globalComm, &sreq);
439
440 // determine the offset in ranks
441 offset_coords[0] = 0;
442 offset_coords[1] = this->local_coords[1] + 1;
443 offset_coords[2] = this->local_coords[2];
444
445 rem_coords[1] = offset_coords[1];
446 rem_coords[2] = offset_coords[2];
447
448 MPI_Cart_rank(this->globalComm, offset_coords, &rank_offset);
449
450 // wait for communication
451 MPI_Wait(&sreq, &sstat);
452 MPI_Wait(&rreq, &rstat);
453
454 // iterate about neighbor borders to determine the neighborship relation
455 nNeighbors.at(3) = 0;
456 for (int x = 0; x < this->global_dims[0]; ++x) {
457 if ((vertices_rem[2 * x] <= vertices_loc[0] &&
458 vertices_loc[0] < vertices_rem[2 * x + 1]) ||
459 (vertices_rem[2 * x] < vertices_loc[1] &&
460 vertices_loc[1] <= vertices_rem[2 * x + 1]) ||
461 (vertices_rem[2 * x] >= vertices_loc[0] &&
462 vertices_loc[0] < vertices_rem[2 * x + 1] &&
463 vertices_loc[1] >= vertices_rem[2 * x + 1])) {
464 nNeighbors.at(3)++;
465 rem_coords[0] = x;
466 MPI_Cart_rank(this->globalComm, rem_coords, &rem_rank);
467 neighbors.push_back(rem_rank);
468 }
469 }
470
471 // barrier to ensure every process concluded the calculations before
472 // overwriting remote borders!
473 MPI_Barrier(this->globalComm);
474
475 // find Z-neighbors to get border information from
476 MPI_Cart_shift(this->globalComm, 2, 1, &rank_left, &rank_right);
477
478 // collect border information from local column
479 vertices_loc[0] = shifted_vertices.at(0)[0];
480 vertices_loc[1] = shifted_vertices.at(1)[0];
481 vertices_loc[2] = shifted_vertices.at(0)[1];
482 vertices_loc[3] = shifted_vertices.at(1)[1];
483
484 MPI_Barrier(this->globalComm);
485
486 MPI_Allgather(vertices_loc, 4, MPIDataTypeT,
487 vertices_rem + 4 * this->global_dims[0] * this->global_dims[1],
488 4, MPIDataTypeT, communicators[2]);
489
490 // exchange local column information with upper neighbor in Z direction (cart
491 // grid)
492 MPI_Irecv(vertices_rem, 4 * this->global_dims[0] * this->global_dims[1],
493 MPIDataTypeT, rank_left, 0, this->globalComm, &rreq);
494 MPI_Isend(vertices_rem + 4 * this->global_dims[0] * this->global_dims[1],
495 4 * this->global_dims[0] * this->global_dims[1], MPIDataTypeT,
496 rank_right, 0, this->globalComm, &sreq);
497
498 // determine the offset in ranks
499 offset_coords[0] = 0;
500 offset_coords[1] = 0;
501 offset_coords[2] = this->local_coords[2] - 1;
502
503 rem_coords[2] = offset_coords[2];
504
505 MPI_Cart_rank(this->globalComm, offset_coords, &rank_offset);
506
507 // wait for communication
508 MPI_Wait(&sreq, &sstat);
509 MPI_Wait(&rreq, &rstat);
510
511 // iterate about neighbor borders to determine the neighborship relation
512 nNeighbors.at(4) = 0;
513 for (int y = 0; y < this->global_dims[1]; ++y) {
514 for (int x = 0; x < this->global_dims[0]; ++x) {
515 if (((vertices_rem[4 * (x + y * this->global_dims[0]) + 2] <=
516 vertices_loc[2] &&
517 vertices_loc[2] <
518 vertices_rem[4 * (x + y * this->global_dims[0]) + 3]) ||
519 (vertices_rem[4 * (x + y * this->global_dims[0]) + 2] <
520 vertices_loc[3] &&
521 vertices_loc[3] <=
522 vertices_rem[4 * (x + y * this->global_dims[0]) + 3]) ||
523 (vertices_rem[4 * (x + y * this->global_dims[0]) + 2] >=
524 vertices_loc[2] &&
525 vertices_loc[2] <
526 vertices_rem[4 * (x + y * this->global_dims[0]) + 3] &&
527 vertices_loc[3] >=
528 vertices_rem[4 * (x + y * this->global_dims[0]) + 3])))
529 if (((vertices_rem[4 * (x + y * this->global_dims[0])] <=
530 vertices_loc[0] &&
531 vertices_loc[0] <
532 vertices_rem[4 * (x + y * this->global_dims[0]) + 1]) ||
533 (vertices_rem[4 * (x + y * this->global_dims[0])] <
534 vertices_loc[1] &&
535 vertices_loc[1] <=
536 vertices_rem[4 * (x + y * this->global_dims[0]) + 1]) ||
537 (vertices_rem[4 * (x + y * this->global_dims[0])] >=
538 vertices_loc[0] &&
539 vertices_loc[0] <
540 vertices_rem[4 * (x + y * this->global_dims[0]) + 1] &&
541 vertices_loc[1] >=
542 vertices_rem[4 * (x + y * this->global_dims[0]) + 1]))) {
543 nNeighbors.at(4)++;
544 rem_coords[1] = y;
545 rem_coords[0] = x;
546 MPI_Cart_rank(this->globalComm, rem_coords, &rem_rank);
547 neighbors.push_back(rem_rank);
548 }
549 }
550 }
551
552 // barrier to ensure every process concluded the calculations before
553 // overwriting remote borders!
554 MPI_Barrier(this->globalComm);
555
556 // exchange local column information with upper neighbor in Y direction (cart
557 // grid)
558 MPI_Irecv(vertices_rem, 4 * this->global_dims[0] * this->global_dims[1],
559 MPIDataTypeT, rank_right, 0, this->globalComm, &rreq);
560 MPI_Isend(vertices_rem + 4 * this->global_dims[0] * this->global_dims[1],
561 4 * this->global_dims[0] * this->global_dims[1], MPIDataTypeT,
562 rank_left, 0, this->globalComm, &sreq);
563
564 // determine the offset in ranks
565 offset_coords[0] = 0;
566 offset_coords[1] = 0;
567 offset_coords[2] = this->local_coords[2] + 1;
568
569 rem_coords[2] = offset_coords[2];
570
571 MPI_Cart_rank(this->globalComm, offset_coords, &rank_offset);
572
573 // wait for communication
574 MPI_Wait(&sreq, &sstat);
575 MPI_Wait(&rreq, &rstat);
576
577 // iterate about neighbor borders to determine the neighborship relation
578 nNeighbors.at(5) = 0;
579 for (int y = 0; y < this->global_dims[1]; ++y) {
580 for (int x = 0; x < this->global_dims[0]; ++x) {
581 if (((vertices_rem[4 * (x + y * this->global_dims[0]) + 2] <=
582 vertices_loc[2] &&
583 vertices_loc[2] <
584 vertices_rem[4 * (x + y * this->global_dims[0]) + 3]) ||
585 (vertices_rem[4 * (x + y * this->global_dims[0]) + 2] <
586 vertices_loc[3] &&
587 vertices_loc[3] <=
588 vertices_rem[4 * (x + y * this->global_dims[0]) + 3]) ||
589 (vertices_rem[4 * (x + y * this->global_dims[0]) + 2] >=
590 vertices_loc[2] &&
591 vertices_loc[2] <
592 vertices_rem[4 * (x + y * this->global_dims[0]) + 3] &&
593 vertices_loc[3] >=
594 vertices_rem[4 * (x + y * this->global_dims[0]) + 3])))
595 if (((vertices_rem[4 * (x + y * this->global_dims[0])] <=
596 vertices_loc[0] &&
597 vertices_loc[0] <
598 vertices_rem[4 * (x + y * this->global_dims[0]) + 1]) ||
599 (vertices_rem[4 * (x + y * this->global_dims[0])] <
600 vertices_loc[1] &&
601 vertices_loc[1] <=
602 vertices_rem[4 * (x + y * this->global_dims[0]) + 1]) ||
603 (vertices_rem[4 * (x + y * this->global_dims[0])] >=
604 vertices_loc[0] &&
605 vertices_loc[0] <
606 vertices_rem[4 * (x + y * this->global_dims[0]) + 1] &&
607 vertices_loc[1] >=
608 vertices_rem[4 * (x + y * this->global_dims[0]) + 1]))) {
609 nNeighbors.at(5)++;
610 rem_coords[1] = y;
611 rem_coords[0] = x;
612 MPI_Cart_rank(this->globalComm, rem_coords, &rem_rank);
613 neighbors.push_back(rem_rank);
614 }
615 }
616 }
617
618 // barrier to ensure every process concluded the calculations before
619 // overwriting remote borders!
620 MPI_Barrier(this->globalComm);
621
622 // clear up vertices array
623 delete[] vertices_loc;
624 delete[] vertices_rem;
625}
626
627// provide list of neighbors
628template <class T, class W>
630 return neighbors;
631}
632
633}//namespace ALL
634
635#endif
636// vim: sw=2 et ts=2
#define known_unused
Definition ALL_Defines.h:49
LB(const int dim, const T g)
Definition ALL_LB.hpp:49
virtual void setVertices(const std::vector< Point< T > > &vertices_in)
Definition ALL_LB.hpp:76
int localRank
local rank in the used MPI communicator
Definition ALL_LB.hpp:212
std::vector< W > work
local work
Definition ALL_LB.hpp:219
MPI_Comm globalComm
used MPI communicator
Definition ALL_LB.hpp:210
int dimension
dimension of the used vertices
Definition ALL_LB.hpp:208
std::vector< int > global_dims
dimensions of the global process grid
Definition ALL_LB.hpp:225
virtual void setWork(const std::vector< W > &w)
Definition ALL_LB.hpp:112
T gamma
correction factor
Definition ALL_LB.hpp:206
virtual int getDimension()
Definition ALL_LB.hpp:99
std::vector< int > local_coords
cartesian coordinates of the local domain in the process grid
Definition ALL_LB.hpp:227
std::vector< Point< T > > prevVertices
original vertices before previous balancing step
Definition ALL_LB.hpp:223
std::vector< Point< T > > vertices
local vertices after previous balancing step
Definition ALL_LB.hpp:221
~Staggered_LB()
default destructor
virtual void setAdditionalData(known_unused const void *data) override
void balance(int step) override
virtual W getEstimatedEfficiency() override
Staggered_LB(int d, W w, T g)
Staggered_LB()
default constructor
void setup() override
method to setup the method specific internal parameters
virtual std::vector< int > & getNeighbors() override
T borderShift1d(const int remote_rank, const int local_coord, const int global_dim, const W local_work, const W remote_work, const T local_size, const T remote_size, const T gamma, const T minSize)
Definition ALL.hpp:75