ALL 0.9.3
A Loadbalacing Library
Loading...
Searching...
No Matches
ALL_Tensor.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_TENSOR_HEADER_INCLUDED
32#define ALL_TENSOR_HEADER_INCLUDED
33
35#include "ALL_Defines.h"
36#include "ALL_Functions.hpp"
37#include "ALL_LB.hpp"
38#include <exception>
39#include <mpi.h>
40
41namespace ALL {
42
53template <class T, class W> class Tensor_LB : public LB<T, W> {
54public:
61 Tensor_LB(int d, W w, T g) : LB<T, W>(d, g) {
62 this->setWork(w);
63 // need the lower and upper bounds
64 // of the domain for the tensor-based
65 // load-balancing scheme
66
67 // array of MPI communicators for each direction (to collect work
68 // on each plane)
69 communicators.resize(d);
70 nNeighbors.resize(2 * d);
71 }
72
74 ~Tensor_LB() = default;
75
77 void setup() override;
78
81 virtual void balance(int step) override {balance(step, MPI_SUM);}
82
83 // getter for variables (passed by reference to avoid
84 // direct access to private members of the object)
85
89 virtual std::vector<int> &getNeighbors() override;
90
93 virtual void setAdditionalData(known_unused const void *data) override {}
94
97 virtual void balance(int step, MPI_Op reductionMode);
98
102 virtual W getEstimatedEfficiency() override {return (W)-1;};
103
104private:
105 // type for MPI communication
106 MPI_Datatype MPIDataTypeT;
107 MPI_Datatype MPIDataTypeW;
108
109 // array of MPI communicators for each direction (to collect work
110 // on each plane)
111 std::vector<MPI_Comm> communicators;
112
113 // list of neighbors
114 std::vector<int> neighbors;
115 std::vector<int> nNeighbors;
116
117protected:
118
119
120};
121
122// setup routine for the tensor-based load-balancing scheme
123// requires:
124// this->globalComm (int): cartesian MPI communicator, from
125// which separate sub communicators
126// are derived in order to represent
127// each plane of domains in the system
128template <class T, class W> void Tensor_LB<T, W>::setup() {
129 int status;
130
131 // check if Communicator is cartesian
132 MPI_Topo_test(this->globalComm, &status);
133 if (status != MPI_CART) {
135 __FILE__, __func__, __LINE__,
136 "Cartesian MPI communicator required, passed communicator not "
137 "cartesian");
138 }
139
140 int dim = this->getDimension();
141
142 // get the local coordinates, periodicity and global size from the MPI
143 // communicator
144 MPI_Cart_get(this->globalComm, dim, this->global_dims.data(),
145 this->periodicity.data(), this->local_coords.data());
146
147 // get the local rank from the MPI communicator
148 MPI_Cart_rank(this->globalComm, this->local_coords.data(), &this->localRank);
149
150 // create sub-communicators
151 for (int i = 0; i < dim; ++i) {
152 MPI_Comm_split(this->globalComm, this->local_coords.at(i), 0,
153 &communicators[i]);
154 }
155
156 // determine correct MPI data type for template T
157 if (std::is_same<T, double>::value)
158 MPIDataTypeT = MPI_DOUBLE;
159 else if (std::is_same<T, float>::value)
160 MPIDataTypeT = MPI_FLOAT;
161 else if (std::is_same<T, int>::value)
162 MPIDataTypeT = MPI_INT;
163 else if (std::is_same<T, long>::value)
164 MPIDataTypeT = MPI_LONG;
165
166 // determine correct MPI data type for template W
167 if (std::is_same<W, double>::value)
168 MPIDataTypeW = MPI_DOUBLE;
169 else if (std::is_same<W, float>::value)
170 MPIDataTypeW = MPI_FLOAT;
171 else if (std::is_same<W, int>::value)
172 MPIDataTypeW = MPI_INT;
173 else if (std::is_same<W, long>::value)
174 MPIDataTypeW = MPI_LONG;
175
176 // calculate neighbors
177 int rank_left, rank_right;
178
179 neighbors.clear();
180 for (int i = 0; i < dim; ++i) {
181 MPI_Cart_shift(this->globalComm, i, 1, &rank_left, &rank_right);
182 neighbors.push_back(rank_left);
183 neighbors.push_back(rank_right);
184 nNeighbors[2 * i] = 1;
185 nNeighbors[2 * i + 1] = 1;
186 }
187}
188
189template <class T, class W> void Tensor_LB<T, W>::balance(int step, MPI_Op reductionMode) {
190 int dim = this->getDimension();
191 std::vector<Point<T>> newVertices = this->vertices;
192 this->prevVertices = this->vertices;
193
194 bool localIsSum = reductionMode == MPI_SUM;
195 bool localIsMax = reductionMode == MPI_MAX;
196
197 // loop over all available dimensions
198 for (int i = 0; i < dim; ++i) {
199 W work_local_plane;
200 // collect work from all processes in the same plane
201 MPI_Allreduce(this->getWork().data(), &work_local_plane, 1, MPIDataTypeW,
202 reductionMode, communicators.at(i));
203
204 // correct right border:
205
206 W remote_work;
207 T local_size;
208 T remote_size;
209 // determine neighbors
210 int rank_left, rank_right;
211
212 MPI_Cart_shift(this->globalComm, i, 1, &rank_left, &rank_right);
213
214 // collect work from right neighbor plane
215 MPI_Request sreq, rreq;
216 MPI_Status sstat, rstat;
217
218 MPI_Irecv(&remote_work, 1, MPIDataTypeW, rank_right, 0, this->globalComm,
219 &rreq);
220 MPI_Isend(&work_local_plane, 1, MPIDataTypeW, rank_left, 0,
221 this->globalComm, &sreq);
222 MPI_Wait(&sreq, &sstat);
223 MPI_Wait(&rreq, &rstat);
224
225 // collect size in dimension from right neighbor plane
226
227 local_size = this->prevVertices.at(1)[i] - this->prevVertices.at(0)[i];
228
229 MPI_Irecv(&remote_size, 1, MPIDataTypeT, rank_right, 0, this->globalComm,
230 &rreq);
231 MPI_Isend(&local_size, 1, MPIDataTypeT, rank_left, 0, this->globalComm,
232 &sreq);
233 MPI_Wait(&sreq, &sstat);
234 MPI_Wait(&rreq, &rstat);
235
236 // automatic selection of gamma to guarantee stability of method (choice of
237 // user gamma ignored)
238 this->gamma =
239 std::max(4.1, 2.0 * (1.0 + std::max(local_size, remote_size) /
240 std::min(local_size, remote_size)));
241
242 T shift = Functions::borderShift1d(
243 rank_right, this->local_coords.at(i), this->global_dims.at(i),
244 work_local_plane, remote_work, local_size, remote_size, this->gamma,
245 this->minSize[i]);
246
247 // send shift to right neighbors
248 T remote_shift = (T)0;
249
250 MPI_Irecv(&remote_shift, 1, MPIDataTypeT, rank_left, 0, this->globalComm,
251 &rreq);
252 MPI_Isend(&shift, 1, MPIDataTypeT, rank_right, 0, this->globalComm, &sreq);
253 MPI_Wait(&sreq, &sstat);
254 MPI_Wait(&rreq, &rstat);
255
256 // if a left neighbor exists: shift left border
257 if (rank_left != MPI_PROC_NULL && this->local_coords[i] != 0)
258 newVertices.at(0)[i] = this->prevVertices.at(0)[i] + remote_shift;
259 else
260 newVertices.at(0)[i] = this->prevVertices.at(0)[i];
261
262 // if a right neighbor exists: shift right border
263 if (rank_right != MPI_PROC_NULL &&
264 this->local_coords[i] != this->global_dims[i] - 1)
265 newVertices.at(1)[i] = this->prevVertices.at(1)[i] + shift;
266 else
267 newVertices.at(1)[i] = this->prevVertices.at(1)[i];
268 }
269
270 this->setVertices(newVertices);
271}
272
273// provide list of neighbors
274template <class T, class W> std::vector<int> &Tensor_LB<T, W>::getNeighbors() {
275 return neighbors;
276}
277
289template <class T, class W> class TensorMax_LB : public Tensor_LB<T, W> {
290public:
292 TensorMax_LB(int d, W w, T g) : Tensor_LB<T, W>(d, w, g) {}
293 virtual void balance(int step) override {Tensor_LB<T,W>::balance(step, MPI_MAX);}
294};
295
296} // namespace ALL
297
298#endif
#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
MPI_Comm globalComm
used MPI communicator
Definition ALL_LB.hpp:210
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
virtual std::vector< W > & getWork()
Definition ALL_LB.hpp:144
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
virtual void balance(int step) override
TensorMax_LB(int d, W w, T g)
virtual W getEstimatedEfficiency() override
~Tensor_LB()=default
default destructor
void setup() override
setup internal data structures and parameters
Tensor_LB()
default constuctor
virtual void setAdditionalData(known_unused const void *data) override
virtual std::vector< int > & getNeighbors() override
Tensor_LB(int d, W w, T g)
virtual void balance(int step) 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