31#ifndef ALL_LB_HEADER_INCLUDED
32#define ALL_LB_HEADER_INCLUDED
43template <
class T,
class W>
class LB {
52 minSize = std::vector<T>(dim, 0.0);
60 virtual ~LB() =
default;
91 __FILE__, __func__, __LINE__,
92 "Currently only three dimensional vertices are supported.");
184 double localSum = 0.0;
185 for (
auto i = this->
work.begin(); i < this->
work.end(); ++i)
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));
virtual void setMinDomainSize(const std::vector< T > &minSize)
LB(const int dim, const T g)
virtual void setVertices(const std::vector< Point< T > > &vertices_in)
void resizeVertices(const int newSize)
virtual void setCommunicator(const MPI_Comm comm)
virtual void setDimension(const int d)
virtual void setWork(const W w)
virtual void setAdditionalData(const void *data)=0
int localRank
local rank in the used MPI communicator
std::vector< W > work
local work
MPI_Comm globalComm
used MPI communicator
int dimension
dimension of the used vertices
std::vector< T > sysSize
(orthogonal) system size
virtual std::vector< Point< T > > & getPrevVertices()
std::vector< T > & getNeighborVertices()
std::vector< int > periodicity
periodicity of the MPI communicator / system
std::vector< int > global_dims
dimensions of the global process grid
virtual void setup()=0
abstract definition of the setup method
virtual const T getGamma()
virtual void setWork(const std::vector< W > &w)
virtual void setGamma(const T g)
virtual void setSysSize(const std::vector< T > &newSysSize)
virtual W getEstimatedEfficiency()=0
virtual std::vector< W > & getWork()
virtual int getDimension()
std::vector< T > neighborVertices
vertices describing neighboring domains
virtual const std::vector< T > & getSysSize()
virtual ~LB()=default
destructor
std::vector< int > local_coords
cartesian coordinates of the local domain in the process grid
virtual std::vector< Point< T > > & getVertices()
std::vector< Point< T > > prevVertices
original vertices before previous balancing step
virtual const std::vector< T > & getMinDomainSize()
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