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