ALL 0.9.3
A Loadbalacing Library
Loading...
Searching...
No Matches
ALL_Staggered.cpp
Go to the documentation of this file.
1/*
2 Copyright 2020-2020 Stephan Schulz, Forschungszentrum Juelich GmbH, Germany
3 Copyright 2018-2024 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 <time.h>
36#include <mpi.h>
37#include <string>
38#include <fstream>
39#include <iostream>
40
41// #define ALL_VTK_OUTPUT
42
43#include <ALL.hpp>
44
45// Run fun in order of ranks
46// 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.
47// Adding sleep(1) also orders everything correctly. So this is probably a flushing problem.
48// It also exists for the cout stream with endl.
49#define MPI_RUN_ORDER(comm, rank, max_ranks, fun) \
50 { \
51 int MPI_RO_IT; \
52 for (MPI_RO_IT = 0; MPI_RO_IT < max_ranks; MPI_RO_IT++) \
53 { \
54 if (MPI_RO_IT == rank) \
55 { \
56 fun; \
57 MPI_Barrier(comm); \
58 } \
59 else \
60 { \
61 MPI_Barrier(comm); \
62 } \
63 } \
64 }
65
66// Quick and dirty helper function. Assumes comm, rank and max_ranks
67// CAVEAT: the function call must be wrapped in () if it contains a comma
68#define MPI_RUN_ORDER_DEF(fun) MPI_RUN_ORDER(MPI_COMM_WORLD, MyRank, MaximumRank, fun)
69
70// void millisleep(unsigned long ms)
71//{
72// struct timespec SleepTime;
73// SleepTime.tv_sec = ms/1000;
74// SleepTime.tv_nsec = (ms%1000)*1000000L;
75// struct timespec RemainingTime;
76// while(nanosleep(&SleepTime, &RemainingTime))
77// {
78// SleepTime.tv_sec = RemainingTime.tv_sec;
79// SleepTime.tv_nsec = RemainingTime.tv_nsec;
80// }
81// }
82
83void print_width(int rank, double width, double bottom, double top)
84{
85 printf("[%03d] Result Width: %10.6f (%10.6f -- %10.6f)\n", rank, width, bottom, top);
86 fflush(stdout);
87}
88
89void print_testing_output(int rank, std::vector<ALL::Point<double>> &vertices, int timestep)
90{
91 // printf("[%4d,%03d] Result Width: %10.6f %10.6f %10.6f\n",
92 // timestep,
93 // rank,
94 // vertices.at(1)[0]-vertices.at(0)[0],
95 // vertices.at(1)[1]-vertices.at(0)[1],
96 // vertices.at(1)[2]-vertices.at(0)[2]);
97 // fflush(stdout);
98 printf("[%4d,%03d,%02d] Result Vertex: %10.6f %10.6f %10.6f\n",
99 timestep,
100 rank,
101 0,
102 vertices.at(0)[0],
103 vertices.at(0)[1],
104 vertices.at(0)[2]);
105 fflush(stdout);
106 printf("[%4d,%03d,%02d] Result Vertex: %10.6f %10.6f %10.6f\n",
107 timestep,
108 rank,
109 1,
110 vertices.at(1)[0],
111 vertices.at(1)[1],
112 vertices.at(1)[2]);
113 fflush(stdout);
114}
115
116void print_loc(int rank, int *loc, int *size)
117{
118 printf("[%03d] Location: (%3d,%3d,%3d)/(%3d,%3d,%3d)\n", rank, loc[0], loc[1], loc[2], size[0], size[1], size[2]);
119 fflush(stdout);
120}
121
122void print_domain(int rank, double *verts)
123{
124 printf("[%03d] Lower: %g\t%g\t%g\n", rank, verts[0], verts[1], verts[2]);
125 printf("[%03d] Upper: %g\t%g\t%g\n", rank, verts[3], verts[4], verts[5]);
126 fflush(stdout);
127}
128
129void print_work(int rank, double work)
130{
131 printf("[%03d] Work: %g\n", rank, work);
132 fflush(stdout);
133}
134
135void print_binary(int step, int rank, double work, std::vector<ALL::Point<double>> &vertices, int *loc, int *size, MPI_File fh)
136{
137 int nranks = size[0] * size[1] * size[2];
138
139 MPI_Offset offset = (step * size[0] * size[1] * size[2] * 11 + rank * 11) * sizeof(double);
140
141 double buf[11];
142
143 buf[0] = (double)rank;
144 buf[1] = work;
145 buf[2] = vertices.at(0)[0];
146 buf[3] = vertices.at(0)[1];
147 buf[4] = vertices.at(0)[2];
148 buf[5] = vertices.at(1)[0];
149 buf[6] = vertices.at(1)[1];
150 buf[7] = vertices.at(1)[2];
151 buf[8] = (double)(loc[0]);
152 buf[9] = (double)(loc[1]);
153 buf[10] = (double)(loc[2]);
154
155 MPI_Status state;
156 MPI_File_write_at(fh, offset, buf, 11, MPI_DOUBLE, &state);
157}
158
159void convert_verts(std::vector<ALL::Point<double>> *vv, double *verts)
160{
161 verts[0] = vv->at(0)[0];
162 verts[1] = vv->at(0)[1];
163 verts[2] = vv->at(0)[2];
164 verts[3] = vv->at(1)[0];
165 verts[4] = vv->at(1)[1];
166 verts[5] = vv->at(1)[2];
167}
168
169int main(int argc, char **argv)
170{
171 MPI_Init(&argc, &argv);
172 int CurrentStep = 0;
173 const int NumberOfSteps = 50;
174
175 // set the output directory for the binary file
176 std::string outputDir;
177 if (argc > 1)
178 {
179 outputDir.assign(argv[1]);
180 }
181 else
182 outputDir.assign(".");
183
184 const int Dimensions = 3;
185 const int LoadbalancerGamma = 0; // ignored for staggered method
186 // Todo(s.schulz): is gamma also ignored for tensor method?
187 // Todo(s.schulz): maybe change to automatic destruction on scope end
188 // The ALL::TENSOR method can be used just like the staggered method.
189 // Just exchange ALL::STAGGERED with ALL::TENSOR in the constructor.
190 ALL::ALL<double, double> *jall = new ALL::ALL<double, double>(ALL::STAGGERED, Dimensions, LoadbalancerGamma);
191 int MyLocation[3] = {0};
192 // All domains are placed along a line in z direction, even though they are three dimensional
193 MPI_Comm_rank(MPI_COMM_WORLD, &MyLocation[2]);
194 int MyRank = MyLocation[2];
195
196 int NumberOfProcesses[3] = {1, 1, 1};
197 MPI_Comm_size(MPI_COMM_WORLD, &NumberOfProcesses[2]);
198 int MaximumRank = NumberOfProcesses[2];
199
200 std::stringstream ss;
201 ss << outputDir << "/" << "ALL_Staggered_" << NumberOfProcesses[2] << ".bin";
202 std::string file = ss.str();
203 MPI_File fh;
204 MPI_File_delete(file.c_str(), MPI_INFO_NULL);
205 MPI_File_open(MPI_COMM_WORLD, file.c_str(), MPI_MODE_CREATE | MPI_MODE_WRONLY, MPI_INFO_NULL, &fh);
206 MPI_File_set_size(fh, NumberOfSteps * NumberOfProcesses[2] * 11 * sizeof(MPI_DOUBLE));
207 MPI_File_preallocate(fh, NumberOfSteps * NumberOfProcesses[2] * 11 * sizeof(MPI_DOUBLE));
208
209 if (MyRank == 0)
210 {
211 printf("Ranks: %d\nNumber of Steps: %d\n", MaximumRank, NumberOfSteps);
212 fflush(stdout);
213 }
214 MPI_Barrier(MPI_COMM_WORLD);
215 MPI_RUN_ORDER_DEF((print_loc(MyRank, MyLocation, NumberOfProcesses)));
216
217 // For a cartesian communicator this is not required, but we are using
218 // MPI_COMM_WORLD here.
219 std::vector<int> MyLocationVector(MyLocation, MyLocation + 3);
220 std::vector<int> NumberOfProcessesVector(NumberOfProcesses, NumberOfProcesses + 3);
221 jall->setProcGridParams(MyLocationVector, NumberOfProcessesVector);
222
223 // If we want to set a minimum domain size:
224 std::vector<double> MinimumDomainSize{0.1, 0.1, 0.1};
225 jall->setMinDomainSize(MinimumDomainSize);
226
227 jall->setCommunicator(MPI_COMM_WORLD);
228
229 // We also set the optional process tag for the output.
230 // This can be useful if we want to know which of 'our'
231 // ranks is which in the output produces by the library.
232 // The ranks used inside the library do not necessarily
233 // match our own numbering.
234 jall->setProcTag(MyRank);
235
236 jall->setup();
237 // Todo(s.schulz): document: what exactly must be set before setup()?
238
239 // A first domain distribution must be given to the balancer.
240 // We use the provided ALL::Point class to define the vertices,
241 // but a simple double array can also be used. We need 2 vertices
242 // which correspond to the two opposing corners.
243 std::vector<ALL::Point<double>> DomainVertices(2, ALL::Point<double>(3));
244 const double DomainSize = 1.0; // Domain size
245 // We create a cubic domain initially
246 for (int VertexIndex = 0; VertexIndex < 2; VertexIndex++)
247 {
248 for (int DimensionIndex = 0; DimensionIndex < Dimensions; DimensionIndex++)
249 {
250 DomainVertices.at(VertexIndex)[DimensionIndex] = (MyLocation[DimensionIndex] + VertexIndex) * DomainSize;
251 }
252 }
253 double VertexArray[6];
254 convert_verts(&DomainVertices, VertexArray);
255 MPI_RUN_ORDER_DEF((print_domain(MyRank, VertexArray)));
256 jall->setVertices(DomainVertices);
257
258 // Calculate the work of our domain. Here we just use
259 double MyWork = (double)MyRank + 1.;
260 jall->setWork(MyWork);
261 MPI_RUN_ORDER_DEF((print_work(MyRank, MyWork)));
262 for (CurrentStep = 0; CurrentStep < NumberOfSteps; CurrentStep++)
263 {
264 // In a real code we need to set the updated work in each
265 // iteration before calling balance()
266 if (MyRank == 0)
267 {
268 printf("Starting step: %d/%d\n", CurrentStep + 1, NumberOfSteps);
269 fflush(stdout);
270 }
271#ifdef ALL_VTK_OUTPUT_EXAMPLE
272 try
273 {
274 jall->printVTKoutlines(CurrentStep);
275 }
277 {
278 std::cout << e.what() << std::endl;
279 std::exit(2);
280 // Maybe also treat this error in some way, e.g. whether the simulation
281 // should abort if no output could be written or not.
282 }
283#endif
284 jall->balance();
285
286 std::vector<ALL::Point<double>> NewVertices = jall->getVertices();
287 // 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])));
288 // MPI_RUN_ORDER_DEF((print_testing_output(MyRank, NewVertices, CurrentStep + 1)));
289 print_binary(CurrentStep, MyRank, MyWork, NewVertices, MyLocation, NumberOfProcesses, fh);
290
291 // Maybe print our new domain? Or not..
292 // convert_verts(&NewVertices, VertexArray);
293 // MPI_RUN_ORDER_DEF((print_domain(MyRank, VertexArray)));
294 // jall->getWork(MyWork);
295 // MPI_RUN_ORDER_DEF((print_work(MyRank,MyWork)));
296 // MPI_Barrier(MPI_COMM_WORLD);
297 }
298#ifdef ALL_VTK_OUTPUT_EXAMPLE
299 try
300 {
301 jall->printVTKoutlines(CurrentStep);
302 }
304 {
305 std::cout << e.what() << std::endl;
306 }
307#endif
308
309 delete jall;
310 MPI_File_close(&fh);
311
312 /* Check correctness of binary file (if in doubt)*/
313 /*
314 if (MyRank == 0)
315 {
316 std::string filename = ss.str();
317 std::ifstream infile(filename, std::ios::binary);
318
319 double d;
320 for (int nSteps = 0; nSteps < NumberOfSteps; ++nSteps)
321 for (int n = 0; n < NumberOfProcesses[2]; ++n)
322 {
323 for (int i = 0; i < 11; ++i)
324 {
325 infile.read(reinterpret_cast<char *>(&d), sizeof d);
326 std::cerr << d << " ";
327 }
328 std::cerr << '\n';
329 }
330 infile.close();
331 }
332 */
333
334 MPI_Finalize();
335 return EXIT_SUCCESS;
336}
void print_loc(int rank, int *loc, int *size)
int main(int argc, char **argv)
void print_binary(int step, int rank, double work, std::vector< ALL::Point< double > > &vertices, int *loc, int *size, MPI_File fh)
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 setMinDomainSize(const std::vector< T > &minSize)
Definition ALL.hpp:510
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
@ STAGGERED
staggered grid load balancing
Definition ALL.hpp:82
virtual const char * what() const