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