ALL 0.9.3
A Loadbalacing Library
Loading...
Searching...
No Matches
ALL_LB.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_LB_HEADER_INCLUDED
32#define ALL_LB_HEADER_INCLUDED
33
34#include "ALL_Point.hpp"
35#include <mpi.h>
36#include <vector>
37#include <iostream>
38
39namespace ALL {
40
43template <class T, class W> class LB {
44public:
49 LB(const int dim, const T g) : gamma(g) {
50 // set allowed minimum size to zero
51 setDimension(dim);
52 minSize = std::vector<T>(dim, 0.0);
53 periodicity.resize(dim);
54 local_coords.resize(dim);
55 global_dims.resize(dim);
56 neighborVertices.resize(0);
57 }
58
60 virtual ~LB() = default;
61
63 virtual void setup() = 0;
64
66 virtual void balance(const int step) = 0;
67
71 virtual std::vector<int> &getNeighbors() = 0;
72
76 virtual void setVertices(const std::vector<Point<T>> &vertices_in) {
77 // as a basic assumption the number of resulting vertices
78 // is equal to the number of input vertices:
79 // exceptions:
80 // VORONOI
81 vertices = vertices_in;
82 prevVertices = vertices_in;
83 }
84
88 virtual void setDimension(const int d) {
89 if (d != 3)
91 __FILE__, __func__, __LINE__,
92 "Currently only three dimensional vertices are supported.");
93 else
94 dimension = d;
95 }
96
99 virtual int getDimension() { return dimension; }
100
103 virtual void setGamma(const T g) { gamma = g; }
104
107 virtual const T getGamma() { return gamma; }
108
112 virtual void setWork(const std::vector<W> &w) { work = w; }
113
116 virtual void setWork(const W w) {
117 work.resize(1);
118 work.at(0) = w;
119 }
120
123 virtual void setMinDomainSize(const std::vector<T> &minSize) {
124 this->minSize = minSize;
125 }
126
129 virtual const std::vector<T> &getMinDomainSize() { return minSize; }
130
133 virtual void setSysSize(const std::vector<T> &newSysSize) {
134 sysSize = newSysSize;
135 }
136
139 virtual const std::vector<T> &getSysSize() { return sysSize; }
140
144 virtual std::vector<W> &getWork() { return work; }
145
148 virtual std::vector<Point<T>> &getVertices() { return vertices; }
149
155 virtual std::vector<Point<T>> &getPrevVertices() { return prevVertices; };
156
159 virtual void setCommunicator(const MPI_Comm comm) {
160 // store communicator
161 globalComm = comm;
162 // get local rank within communicator
163 MPI_Comm_rank(comm, &localRank);
164 }
165
168 int getNVertices() { return vertices.size(); }
169
173 virtual void setAdditionalData(const void *data) = 0;
174
178 virtual W getEstimatedEfficiency() = 0;
179
183 {
184 double localSum = 0.0;
185 for (auto i = this->work.begin(); i < this->work.end(); ++i)
186 {
187 localSum += *i;
188 }
189 double globalMin;
190 double globalMax;
191 MPI_Allreduce(&localSum, &globalMin, 1, MPI_DOUBLE, MPI_MIN, this->globalComm);
192 MPI_Allreduce(&localSum, &globalMax, 1, MPI_DOUBLE, MPI_MAX, this->globalComm);
193 return (1.0 - (double)(globalMax - globalMin) /
194 (double)(globalMax + globalMin));
195 }
196
197
201 std::vector<T> &getNeighborVertices() {
202 return neighborVertices; }
203
204protected:
210 MPI_Comm globalComm;
215 std::vector<T> minSize;
217 std::vector<T> sysSize;
219 std::vector<W> work;
221 std::vector<Point<T>> vertices;
223 std::vector<Point<T>> prevVertices;
225 std::vector<int> global_dims;
227 std::vector<int> local_coords;
229 std::vector<int> periodicity;
231 std::vector<T> neighborVertices;
234 void resizeVertices(const int newSize) { vertices.resize(newSize); }
235
236};
237
238}//namespace ALL
239
240#endif // ALL_LB_HEADER_INCLUDED
virtual void setMinDomainSize(const std::vector< T > &minSize)
Definition ALL_LB.hpp:123
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
void resizeVertices(const int newSize)
Definition ALL_LB.hpp:234
virtual void setCommunicator(const MPI_Comm comm)
Definition ALL_LB.hpp:159
virtual void setDimension(const int d)
Definition ALL_LB.hpp:88
virtual void setWork(const W w)
Definition ALL_LB.hpp:116
virtual void setAdditionalData(const void *data)=0
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
std::vector< T > minSize
Definition ALL_LB.hpp:215
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< T > sysSize
(orthogonal) system size
Definition ALL_LB.hpp:217
virtual std::vector< Point< T > > & getPrevVertices()
Definition ALL_LB.hpp:155
int getNVertices()
Definition ALL_LB.hpp:168
std::vector< T > & getNeighborVertices()
Definition ALL_LB.hpp:201
std::vector< int > periodicity
periodicity of the MPI communicator / system
Definition ALL_LB.hpp:229
std::vector< int > global_dims
dimensions of the global process grid
Definition ALL_LB.hpp:225
virtual void setup()=0
abstract definition of the setup method
virtual const T getGamma()
Definition ALL_LB.hpp:107
virtual void setWork(const std::vector< W > &w)
Definition ALL_LB.hpp:112
virtual void setGamma(const T g)
Definition ALL_LB.hpp:103
virtual void setSysSize(const std::vector< T > &newSysSize)
Definition ALL_LB.hpp:133
double getEfficiency()
Definition ALL_LB.hpp:182
virtual W getEstimatedEfficiency()=0
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< T > neighborVertices
vertices describing neighboring domains
Definition ALL_LB.hpp:231
virtual const std::vector< T > & getSysSize()
Definition ALL_LB.hpp:139
virtual ~LB()=default
destructor
std::vector< int > local_coords
cartesian coordinates of the local domain in the process grid
Definition ALL_LB.hpp:227
virtual std::vector< Point< T > > & getVertices()
Definition ALL_LB.hpp:148
std::vector< Point< T > > prevVertices
original vertices before previous balancing step
Definition ALL_LB.hpp:223
virtual const std::vector< T > & getMinDomainSize()
Definition ALL_LB.hpp:129
virtual void balance(const int step)=0
abstract definition of the balancing method
virtual std::vector< int > & getNeighbors()=0
std::vector< Point< T > > vertices
local vertices after previous balancing step
Definition ALL_LB.hpp:221
Definition ALL.hpp:75