linear_algebra_distribution.cc
Go to the documentation of this file.
1// LIC// ====================================================================
2// LIC// This file forms part of oomph-lib, the object-oriented,
3// LIC// multi-physics finite-element library, available
4// LIC// at http://www.oomph-lib.org.
5// LIC//
6// LIC// Copyright (C) 2006-2022 Matthias Heil and Andrew Hazel
7// LIC//
8// LIC// This library is free software; you can redistribute it and/or
9// LIC// modify it under the terms of the GNU Lesser General Public
10// LIC// License as published by the Free Software Foundation; either
11// LIC// version 2.1 of the License, or (at your option) any later version.
12// LIC//
13// LIC// This library is distributed in the hope that it will be useful,
14// LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16// LIC// Lesser General Public License for more details.
17// LIC//
18// LIC// You should have received a copy of the GNU Lesser General Public
19// LIC// License along with this library; if not, write to the Free Software
20// LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
21// LIC// 02110-1301 USA.
22// LIC//
23// LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
24// LIC//
25// LIC//====================================================================
27
28namespace oomph
29{
30 //============================================================================
31 /// Sets the distribution. Takes first_row, local_nrow and
32 /// global_nrow as arguments. If global_nrow is not provided or equal to
33 /// 0 then it is computed automatically
34 //============================================================================
36 const unsigned& first_row,
37 const unsigned& local_nrow,
38 const unsigned& global_nrow)
39 {
40 // copy the communicator
41 delete Comm_pt;
42 Comm_pt = new OomphCommunicator(*comm_pt);
43
44 // get the rank and the number of processors
45 int my_rank = Comm_pt->my_rank();
46 int nproc = Comm_pt->nproc();
47
48 // resize the storage
49 First_row.clear();
50 First_row.resize(nproc);
51 Nrow_local.clear();
52 Nrow_local.resize(nproc);
53
54 // set first row and local nrow for this processor
55 First_row[my_rank] = first_row;
56 Nrow_local[my_rank] = local_nrow;
57
58#ifdef OOMPH_HAS_MPI
59 // gather the First_row vector
60 unsigned my_nrow_local = Nrow_local[my_rank];
61 MPI_Allgather(&my_nrow_local,
62 1,
63 MPI_UNSIGNED,
64 &Nrow_local[0],
65 1,
66 MPI_UNSIGNED,
67 Comm_pt->mpi_comm());
68
69 // gather the Nrow_local vector
70 unsigned my_first_row = First_row[my_rank];
71 MPI_Allgather(&my_first_row,
72 1,
73 MPI_UNSIGNED,
74 &First_row[0],
75 1,
76 MPI_UNSIGNED,
77 Comm_pt->mpi_comm());
78#endif
79
80 // if global nrow is not provided then compute by summing local_nrow over
81 // all processors
82 if (global_nrow == 0)
83 {
84 if (nproc == 1)
85 {
86 Nrow = local_nrow;
87 }
88 else
89 {
90 Nrow = 0;
91 for (int p = 0; p < nproc; p++)
92 {
93 Nrow += Nrow_local[p];
94 }
95 }
96 }
97 else
98 {
99 Nrow = global_nrow;
100 }
101
102 // the distribution is true, unless we only have 1 processor
103 if (nproc != 1)
104 {
105 Distributed = true;
106 }
107 else
108 {
109 Distributed = false;
110 }
111
112#ifdef OOMPH_HAS_MPI
113#ifdef PARANOID
114 // paranoid check that the distribution works
115
116
117 // check that none of the processors partition overlap
118 for (int p = 0; p < nproc; p++)
119 {
120 if (Nrow_local[p] > 0)
121 {
122 for (int pp = p + 1; pp < nproc; pp++)
123 {
124 if (Nrow_local[pp] > 0)
125 {
126 if ((First_row[p] >= First_row[pp] &&
127 First_row[p] < First_row[pp] + Nrow_local[pp]) ||
128 (First_row[p] + Nrow_local[p] - 1 >= First_row[pp] &&
129 First_row[p] + Nrow_local[p] - 1 <
130 First_row[pp] + Nrow_local[pp]))
131 {
132 std::ostringstream error_message;
133 error_message
134 << "The distributed rows on processor " << p
135 << " and processor " << pp << " overlap.\n"
136 << "Processor " << p << " : first_row = " << First_row[p]
137 << ", nrow = " << Nrow_local[p] << ".\n"
138 << "Processor " << pp << " : first_row = " << First_row[pp]
139 << ", nrow = " << Nrow_local[pp] << ".\n";
140 throw OomphLibWarning(
141 error_message.str(),
142 "LinearAlgebraDistribution::distribute(...)",
143 OOMPH_EXCEPTION_LOCATION);
144 }
145 }
146 }
147 }
148 }
149
150 // check that no processor has a row with a global row index greater than
151 // the number of global rows
152 for (int p = 0; p < nproc; p++)
153 {
154 if (First_row[p] + Nrow_local[p] > Nrow)
155 {
156 std::ostringstream error_message;
157 error_message << "Processor " << p << " contains rows " << First_row[p]
158 << " to " << First_row[p] + Nrow_local[p] - 1
159 << " but there are only " << Nrow << " to be distributed."
160 << std::endl;
161 throw OomphLibWarning(error_message.str(),
162 "LinearAlgebraDistribution::distribute(...)",
163 OOMPH_EXCEPTION_LOCATION);
164 }
165 }
166#endif
167#endif
168 }
169
170 //============================================================================
171 /// Uniformly distribute global_nrow over all processors where
172 /// processors 0 holds approximately the first
173 /// global_nrow/n_proc, processor 1 holds the next
174 /// global_nrow/n_proc and so on...
175 //============================================================================
177 const unsigned& global_nrow,
178 const bool& distribute)
179 {
180 // copy the communicator
181 delete Comm_pt;
182 Comm_pt = new OomphCommunicator(*comm_pt);
183
184 // delete existing storage
185 First_row.clear();
186 Nrow_local.clear();
187
188 // set global nrow
189 Nrow = global_nrow;
190
191 // store the distributed flag
192 Distributed = distribute;
193
194#ifdef OOMPH_HAS_MPI
195
196 // if distributed object then compute uniform distribution
197 if (distribute == true)
198 {
199 // get the number of processors
200 int nproc = Comm_pt->nproc();
201
202 // resize the storage
203 First_row.resize(nproc);
204 Nrow_local.resize(nproc);
205
206 // compute first row
207 for (int p = 0; p < nproc; p++)
208 {
209 First_row[p] =
210 unsigned((double(p) / double(nproc)) * double(global_nrow));
211 }
212
213 // compute local nrow
214 for (int p = 0; p < nproc - 1; p++)
215 {
216 Nrow_local[p] = First_row[p + 1] - First_row[p];
217 }
218 Nrow_local[nproc - 1] = global_nrow - First_row[nproc - 1];
219 }
220#endif
221 }
222
223 //============================================================================
224 /// helper method for the =assignment operator and copy constructor
225 //============================================================================
227 const LinearAlgebraDistribution& new_dist)
228 {
229 // delete the existing storage
230 First_row.clear();
231 Nrow_local.clear();
232
233 // if new_dist is not setup
234 if (new_dist.communicator_pt() == 0)
235 {
236 delete Comm_pt;
237 Comm_pt = 0;
238 Distributed = true;
239 Nrow = 0;
240 }
241 else
242 {
243 // copy the communicator
244 delete Comm_pt;
245 Comm_pt = new OomphCommunicator(*new_dist.communicator_pt());
246
247 // the new distribution is distributed
248 if (new_dist.distributed())
249 {
250 First_row = new_dist.first_row_vector();
251 Nrow_local = new_dist.nrow_local_vector();
252
253 Distributed = true;
254 }
255 // else if the new ditribution is not distributed
256 else
257 {
258 Distributed = false;
259 }
260 Nrow = new_dist.nrow();
261 }
262 }
263
264
265 //============================================================================
266 /// operator==
267 //============================================================================
269 const LinearAlgebraDistribution& other_dist) const
270 {
271#ifdef OOMPH_HAS_MPI
272 // compare the communicators
273 if (!((*Comm_pt) == (*other_dist.communicator_pt())))
274 {
275 return false;
276 }
277
278 // compare Distributed
279 if (Distributed != other_dist.distributed())
280 {
281 return false;
282 }
283
284 // if not distributed compare nrow
285 if (!Distributed)
286 {
287 if (other_dist.nrow() == Nrow)
288 {
289 return true;
290 }
291 return false;
292 }
293
294 // compare
295 bool flag = true;
296 int nproc = Comm_pt->nproc();
297 for (int i = 0; i < nproc && flag == true; i++)
298 {
299 if (other_dist.first_row(i) != First_row[i] ||
300 other_dist.nrow_local(i) != Nrow_local[i])
301 {
302 flag = false;
303 }
304 }
305 return flag;
306#else
307 if (other_dist.nrow() == Nrow)
308 {
309 return true;
310 }
311 return false;
312#endif
313 }
314
315 //=============================================================================
316 /// output operator
317 //=============================================================================
318 std::ostream& operator<<(std::ostream& stream,
320 {
321 stream << "nrow()=" << dist.nrow() << ", first_row()=" << dist.first_row()
322 << ", nrow_local()=" << dist.nrow_local()
323 << ", distributed()=" << dist.distributed() << std::endl;
324 return stream;
325 }
326
327 //=============================================================================
328 /// Namespace for helper functions for LinearAlgebraDistributions
329 //=============================================================================
330 namespace LinearAlgebraDistributionHelpers
331 {
332 //===========================================================================
333 /// Takes a vector of LinearAlgebraDistribution objects and
334 /// concatenates them such that the nrow_local of the out_distribution
335 /// is the sum of the nrow_local of all the in_distributions and the number
336 /// of global rows of the out_distribution is the sum of the number of
337 /// global rows of all the in_distributions. This results in a permutation
338 /// of the rows in the out_distribution. Think of this in terms of
339 /// DoubleVectors, if we have DoubleVectors with distributions A and B,
340 /// distributed across two processors (p0 and p1), A: [a0] (on p0) B:
341 /// [b0] (on p0)
342 /// [a1] (on p1) [b1] (on P1),
343 ///
344 /// then the out_distribution is
345 /// [a0 (on p0)
346 /// b0] (on p0)
347 /// [a1 (on p1)
348 /// b1] (on p1),
349 ///
350 /// as opposed to
351 /// [a0 (on p0)
352 /// a1] (on p0)
353 /// [b0 (on p1)
354 /// b1] (on p1).
355 ///
356 /// Note (1): The out_distribution may not be uniformly distributed even
357 /// if the in_distributions are uniform distributions.
358 /// Try this out with two distributions of global rows 3 and 5, uniformly
359 /// distributed across two processors. Compare this against a distribution
360 /// of global row 8 distributed across two processors.
361 ///
362 /// Note (2): There is no equivalent function which takes a Vector of
363 /// LinearAlgebraDistribution objects (as opposed to pointers), there should
364 /// not be one since we do not want to invoke the assignment operator when
365 /// creating the Vector of LinearAlgebraDistribution objects.
366 //===========================================================================
368 const Vector<LinearAlgebraDistribution*>& in_distribution_pt,
369 LinearAlgebraDistribution& out_distribution)
370 {
371 // How many distributions are in in_distribution?
372 unsigned ndistributions = in_distribution_pt.size();
373
374#ifdef PARANOID
375
376 // Are any of the in_distribution pointers null?
377 // We do not want null pointers.
378 for (unsigned dist_i = 0; dist_i < ndistributions; dist_i++)
379 {
380 if (in_distribution_pt[dist_i] == 0)
381 {
382 std::ostringstream error_message;
383 error_message << "The pointer in_distribution_pt[" << dist_i
384 << "] is null.\n";
385 throw OomphLibError(error_message.str(),
386 OOMPH_CURRENT_FUNCTION,
387 OOMPH_EXCEPTION_LOCATION);
388 }
389 }
390
391 /// /////////////////////////////
392
393 // Check that all in distributions are built.
394 for (unsigned dist_i = 0; dist_i < ndistributions; dist_i++)
395 {
396 if (!in_distribution_pt[dist_i]->built())
397 {
398 std::ostringstream error_message;
399 error_message << "The in_distribution_pt[" << dist_i
400 << "] is not built.\n";
401 throw OomphLibError(error_message.str(),
402 OOMPH_CURRENT_FUNCTION,
403 OOMPH_EXCEPTION_LOCATION);
404 }
405 }
406
407 /// /////////////////////////////
408
409 // Check that all communicators to concatenate are the same
410 // by comparing all communicators against the first one.
411 const OomphCommunicator first_comm =
412 *(in_distribution_pt[0]->communicator_pt());
413
414 for (unsigned dist_i = 0; dist_i < ndistributions; dist_i++)
415 {
416 // Get the Communicator for the current vector.
417 const OomphCommunicator another_comm =
418 *(in_distribution_pt[dist_i]->communicator_pt());
419
420 if (!(first_comm == another_comm))
421 {
422 std::ostringstream error_message;
423 error_message << "The communicator in position " << dist_i << " is \n"
424 << "not the same as the first one.\n";
425 throw OomphLibError(error_message.str(),
426 OOMPH_CURRENT_FUNCTION,
427 OOMPH_EXCEPTION_LOCATION);
428 }
429 }
430
431 /// /////////////////////////////
432
433 // Ensure that all distributions are either distributed or not.
434 // This is because we use the distributed() function from the first
435 // distribution to indicate if the result distribution should be
436 // distributed or not.
437 bool first_distributed = in_distribution_pt[0]->distributed();
438 for (unsigned dist_i = 0; dist_i < ndistributions; dist_i++)
439 {
440 // Is the current distribution distributed?
441 bool another_distributed = in_distribution_pt[dist_i]->distributed();
442 if (first_distributed != another_distributed)
443 {
444 std::ostringstream error_message;
445 error_message
446 << "The distribution in position " << dist_i << " has a different\n"
447 << "different distributed boolean than the distribution.\n";
448 throw OomphLibError(error_message.str(),
449 OOMPH_CURRENT_FUNCTION,
450 OOMPH_EXCEPTION_LOCATION);
451 }
452 }
453
454 /// /////////////////////////////
455
456 // Check that the out distribution is not built.
457 if (out_distribution.built())
458 {
459 std::ostringstream error_message;
460 error_message << "The out distribution is built.\n"
461 << "Please clear it.\n";
462 throw OomphLibError(error_message.str(),
463 OOMPH_CURRENT_FUNCTION,
464 OOMPH_EXCEPTION_LOCATION);
465 }
466
467#endif
468
469 // Get the communicator pointer
470 const OomphCommunicator* const comm_pt =
471 in_distribution_pt[0]->communicator_pt();
472
473 // Number of processors
474 unsigned nproc = comm_pt->nproc();
475
476 // First determine the out_nrow_local and the out_nrow (the global nrow)
477 unsigned out_nrow_local = 0;
478 unsigned out_nrow = 0;
479 for (unsigned in_dist_i = 0; in_dist_i < ndistributions; in_dist_i++)
480 {
481 out_nrow_local += in_distribution_pt[in_dist_i]->nrow_local();
482 out_nrow += in_distribution_pt[in_dist_i]->nrow();
483 }
484
485 // Now calculate the first row for this processor. We need to know the
486 // out_nrow_local for all the other processors, MPI_Allgather must be
487 // used.
488 unsigned out_first_row = 0;
489
490 // Distributed case: We need to gather all the out_nrow_local from all
491 // processors, the out_first_row for this processors is the partial sum up
492 // of the out_nrow_local to my_rank.
493 bool distributed = in_distribution_pt[0]->distributed();
494 if (distributed)
495 {
496#ifdef OOMPH_HAS_MPI
497 // My rank
498 unsigned my_rank = comm_pt->my_rank();
499
500 unsigned* out_nrow_local_all = new unsigned[nproc];
501 MPI_Allgather(&out_nrow_local,
502 1,
503 MPI_UNSIGNED,
504 out_nrow_local_all,
505 1,
506 MPI_UNSIGNED,
507 comm_pt->mpi_comm());
508 for (unsigned proc_i = 0; proc_i < my_rank; proc_i++)
509 {
510 out_first_row += out_nrow_local_all[proc_i];
511 }
512 delete[] out_nrow_local_all;
513#endif
514 }
515 // Non-distributed case
516 else
517 {
518 out_first_row = 0;
519 }
520
521 // Build the distribution.
522 if (nproc == 1 || !distributed)
523 {
524 // Some ambiguity here -- on a single processor it doesn't really
525 // matter if we think of ourselves as distributed or not but this
526 // follows the pattern used elsewhere.
527 out_distribution.build(comm_pt, out_nrow, false);
528 }
529 else
530 {
531 out_distribution.build(
532 comm_pt, out_first_row, out_nrow_local, out_nrow);
533 }
534 } // function concatenate
535
536 } // end of namespace LinearAlgebraDistributionHelpers
537
538} // namespace oomph
cstr elem_len * i
Definition: cfortran.h:603
Describes the distribution of a distributable linear algebra type object. Typically this is a contain...
bool distributed() const
access function to the distributed - indicates whether the distribution is serial or distributed
unsigned first_row() const
access function for the first row on this processor. If not distributed then this is just zero.
OomphCommunicator * communicator_pt() const
const access to the communicator pointer
bool operator==(const LinearAlgebraDistribution &other_dist) const
== Operator
Vector< unsigned > Nrow_local
the number of local rows on the processor
void build(const OomphCommunicator *const comm_pt, const unsigned &first_row, const unsigned &nrow_local, const unsigned &nrow=0)
Sets the distribution. Takes first_row, nrow_local and nrow as arguments. If nrow is not provided or ...
OomphCommunicator * Comm_pt
the pointer to the MPI communicator object in this distribution
unsigned Nrow
the number of global rows
Vector< unsigned > First_row
the first row on this processor
bool built() const
if the communicator_pt is null then the distribution is not setup then false is returned,...
bool Distributed
flag to indicate whether this distribution describes an object that is distributed over the processor...
Vector< unsigned > first_row_vector() const
return the first_row Vector
unsigned nrow() const
access function to the number of global rows.
unsigned nrow_local() const
access function for the num of local rows on this processor. If no MPI then Nrow is returned.
Vector< unsigned > nrow_local_vector() const
return the nrow_local Vector
An oomph-lib wrapper to the MPI_Comm communicator object. Just contains an MPI_Comm object (which is ...
Definition: communicator.h:54
An OomphLibError object which should be thrown when an run-time error is encountered....
An OomphLibWarning object which should be created as a temporary object to issue a warning....
A slight extension to the standard template vector class so that we can include "graceful" array rang...
Definition: Vector.h:58
void concatenate(const Vector< LinearAlgebraDistribution * > &in_distribution_pt, LinearAlgebraDistribution &out_distribution)
Takes a vector of LinearAlgebraDistribution objects and concatenates them such that the nrow_local of...
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
std::ostream & operator<<(std::ostream &out, const DoubleVector &v)
output operator