ALL 0.9.3
A Loadbalacing Library
Loading...
Searching...
No Matches
ALL_Voronoi.cpp
Go to the documentation of this file.
1/*
2 Copyright 2020-2020 Stephan Schulz, Forschungszentrum Juelich GmbH, Germany
3 Copyright 2018-2020 Rene Halver, Forschungszentrum Juelich GmbH, Germany
4 Copyright 2018-2020 Godehard Sutmann, Forschungszentrum Juelich GmbH, Germany
5
6 Redistribution and use in source and binary forms, with or without
7 modification, are permitted provided that the following conditions are met:
8
9 1. Redistributions of source code must retain the above copyright notice,
10 this list of conditions and the following disclaimer.
11
12 2. Redistributions in binary form must reproduce the above copyright notice,
13 this list of conditions and the following disclaimer in the documentation
14 and/or other materials provided with the distribution.
15
16 3. Neither the name of the copyright holder nor the names of its contributors
17 may be used to endorse or promote products derived from this software without
18 specific prior written permission.
19
20 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
21 AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
22 IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
23 ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
24 LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
25 CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
26 SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
27 INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
28 CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
29 ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
30 POSSIBILITY OF SUCH DAMAGE.
31 */
32
33#include <stdlib.h>
34#include <stdio.h>
35#include <mpi.h>
36
37//#define ALL_VTK_OUTPUT
38#include <ALL.hpp>
39
40// Run fun in order of ranks
41// Todo(s.schulz): This seems to only work roughly with the result width an 32 ranks, with up to 16 it seems to work correctly.
42// Adding sleep(1) also orders everything correctly. So this is probably a flushing problem.
43// It also exists for the cout stream with endl.
44#define MPI_RUN_ORDER(comm, rank, max_ranks, fun) {int MPI_RO_IT;\
45 for(MPI_RO_IT=0;MPI_RO_IT<max_ranks;MPI_RO_IT++)\
46 {\
47 if(MPI_RO_IT==rank)\
48 {\
49 fun;\
50 MPI_Barrier(comm);\
51 } else {\
52 MPI_Barrier(comm);\
53 }\
54 }\
55}
56
57// Quick and dirty helper function. Assumes comm, rank and max_ranks
58// CAVEAT: the function call must be wrapped in () if it contains a comma
59#define MPI_RUN_ORDER_DEF(fun) MPI_RUN_ORDER(MPI_COMM_WORLD, MyRank, MaximumRank, fun)
60
61void print_width(int rank, double width, double bottom, double top)
62{
63 printf("[%03d] Result Width: %10.6f (%10.6f -- %10.6f)\n", rank, width, bottom, top);
64 fflush(stdout);
65}
66
67void print_testing_output(int rank, std::vector<ALL::Point<double>>& vertices, int timestep)
68{
69// printf("[%4d,%03d] Result Width: %10.6f %10.6f %10.6f\n",
70// timestep,
71// rank,
72// vertices.at(1)[0]-vertices.at(0)[0],
73// vertices.at(1)[1]-vertices.at(0)[1],
74// vertices.at(1)[2]-vertices.at(0)[2]);
75// fflush(stdout);
76 for(int Vertex=0; Vertex<vertices.size(); Vertex++)
77 {
78 printf("[%4d,%03d,%02d] Result Vertex: %10.6f %10.6f %10.6f\n",
79 timestep,
80 rank,
81 Vertex,
82 vertices.at(Vertex)[0],
83 vertices.at(Vertex)[1],
84 vertices.at(Vertex)[2]);
85 fflush(stdout);
86 }
87}
88
89void print_loc(int rank, int* loc, int* size)
90{
91 printf("[%03d] Location: (%3d,%3d,%3d)/(%3d,%3d,%3d)\n", rank, loc[0], loc[1], loc[2], size[0], size[1], size[2]);
92 fflush(stdout);
93}
94
95void print_domain(int rank, double* verts)
96{
97 printf("[%03d] Lower: %g\t%g\t%g\n", rank, verts[0], verts[1], verts[2]);
98 printf("[%03d] Upper: %g\t%g\t%g\n", rank, verts[3], verts[4], verts[5]);
99 fflush(stdout);
100}
101
102void print_work(int rank, double work)
103{
104 printf("[%03d] Work: %g\n", rank, work);
105 fflush(stdout);
106}
107
108void convert_verts(std::vector<ALL::Point<double>>* vv, double* verts)
109{
110 verts[0] = vv->at(0)[0];
111 verts[1] = vv->at(0)[1];
112 verts[2] = vv->at(0)[2];
113 verts[3] = vv->at(1)[0];
114 verts[4] = vv->at(1)[1];
115 verts[5] = vv->at(1)[2];
116}
117
118int main(int argc, char** argv)
119{
120 MPI_Init(&argc,&argv);
121 int CurrentStep = 0;
122 const int NumberOfSteps = 50;
123
124 const int Dimensions = 3;
125 const int LoadbalancerGamma = 0; // ignored for staggered method
126 ALL::ALL<double, double> *jall = new ALL::ALL<double, double>(ALL::VORONOI, Dimensions, LoadbalancerGamma);
127 int MyLocation[3] = {0};
128 // All domains are placed along a line in z direction, even though they are three dimensional
129 MPI_Comm_rank(MPI_COMM_WORLD,&MyLocation[2]);
130 int MyRank = MyLocation[2];
131
132 int NumberOfProcesses[3] = {1,1,1};
133 MPI_Comm_size(MPI_COMM_WORLD, &NumberOfProcesses[2]);
134 int MaximumRank = NumberOfProcesses[2];
135
136 if(MyRank==0)
137 {
138 printf("Ranks: %d\nNumber of Steps: %d\n", MaximumRank, NumberOfSteps);
139 fflush(stdout);
140 }
141 MPI_Barrier(MPI_COMM_WORLD);
142 MPI_RUN_ORDER_DEF((print_loc(MyRank, MyLocation, NumberOfProcesses)));
143
144 // For a cartesian communicator this is not required, but we are using
145 // MPI_COMM_WORLD here.
146 std::vector<int> MyLocationVector(MyLocation, MyLocation+3);
147 std::vector<int> NumberOfProcessesVector(NumberOfProcesses, NumberOfProcesses+3);
148 jall->setProcGridParams(MyLocationVector, NumberOfProcessesVector);
149
150 jall->setCommunicator(MPI_COMM_WORLD);
151
152 // We also set the optional process tag for the output.
153 // This can be useful if we want to know which of 'our'
154 // ranks is which in the output produces by the library.
155 // The ranks used inside the library do not necessarily
156 // match our own numbering.
157 jall->setProcTag(MyRank);
158
159 jall->setup();
160 // Todo(s.schulz): document: what exactly must be set before setup()?
161
162 // A first domain distribution must be given to the balancer.
163 // We use the provided ALL::Point class to define the vertices,
164 // but a simple double array can also be used. We need 2 vertices
165 // which correspond to the two opposing corners.
166 std::vector<ALL::Point<double>> DomainVertices(2, ALL::Point<double>(3));
167 const double DomainSize = 1.0; // Domain size
168 // We create a cubic domain initially
169 for(int VertexIndex=0; VertexIndex<2; VertexIndex++)
170 {
171 for(int DimensionIndex=0; DimensionIndex<Dimensions; DimensionIndex++)
172 {
173 DomainVertices.at(VertexIndex)[DimensionIndex] = (MyLocation[DimensionIndex]+VertexIndex) * DomainSize;
174 }
175 }
176 double VertexArray[6];
177 convert_verts(&DomainVertices, VertexArray);
178 MPI_RUN_ORDER_DEF((print_domain(MyRank, VertexArray)));
179 jall->setVertices(DomainVertices);
180
181 // Calculate the work of our domain. Here we just use
182 double MyWork = (double) MyRank + 1.;
183 jall->setWork(MyWork);
184 MPI_RUN_ORDER_DEF((print_work(MyRank,MyWork)));
185 for(CurrentStep=0; CurrentStep<NumberOfSteps; CurrentStep++)
186 {
187 // In a real code we need to set the updated work in each
188 // iteration before calling balance()
189 if(MyRank==0)
190 {
191 printf("Starting step: %d/%d\n", CurrentStep+1, NumberOfSteps);
192 fflush(stdout);
193 }
194#ifdef ALL_VTK_OUTPUT_EXAMPLE
195 jall->printVTKoutlines(CurrentStep);
196#endif
197 jall->balance();
198
199 std::vector<ALL::Point<double>> NewVertices = jall->getVertices();
200 //MPI_RUN_ORDER(MPI_COMM_WORLD, MyRank, MaximumRank, (print_width(MyRank, NewVertices.at(1)[2]-NewVertices.at(0)[2], NewVertices.at(0)[2], NewVertices.at(1)[2])));
201 MPI_RUN_ORDER_DEF((print_testing_output(MyRank, NewVertices, CurrentStep+1)));
202 // Maybe print our new domain? Or not..
203 //convert_verts(&NewVertices, VertexArray);
204 //MPI_RUN_ORDER_DEF((print_domain(MyRank, VertexArray)));
205 //jall->getWork(MyWork);
206 //MPI_RUN_ORDER_DEF((print_work(MyRank,MyWork)));
207 MPI_Barrier(MPI_COMM_WORLD);
208 }
209#ifdef ALL_VTK_OUTPUT_EXAMPLE
210 jall->printVTKoutlines(CurrentStep);
211#endif
212
213 delete jall;
214
215 MPI_Finalize();
216 return EXIT_SUCCESS;
217}
int main(int argc, char **argv)
void print_loc(int rank, int *loc, int *size)
void convert_verts(std::vector< ALL::Point< double > > *vv, double *verts)
#define MPI_RUN_ORDER_DEF(fun)
void print_testing_output(int rank, std::vector< ALL::Point< double > > &vertices, int timestep)
void print_width(int rank, double width, double bottom, double top)
void print_domain(int rank, double *verts)
void print_work(int rank, double work)
void setProcGridParams(const std::vector< int > &loc, const std::vector< int > &size)
Definition ALL.hpp:494
void setCommunicator(const MPI_Comm comm)
Definition ALL.hpp:211
std::vector< Point< T > > & getVertices()
Definition ALL.hpp:423
void setWork(const W work)
Definition ALL.hpp:281
void setup()
Definition ALL.hpp:292
std::vector< Point< T > > & balance()
Definition ALL.hpp:299
void setVertices(const std::vector< Point< T > > &inp)
Definition ALL.hpp:195
void printVTKoutlines(const int step)
Definition ALL.hpp:519
void setProcTag(int tag)
Definition ALL.hpp:504
@ VORONOI
voronoi cell based load balancing
Definition ALL.hpp:89