libMesh
distributed_vector.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #include "libmesh/libmesh_common.h"
21 
22 
23 
24 #ifndef LIBMESH_DISTRIBUTED_VECTOR_H
25 #define LIBMESH_DISTRIBUTED_VECTOR_H
26 
27 // Local includes
28 #include "libmesh/numeric_vector.h"
29 #include "libmesh/parallel.h"
30 
31 // C++ includes
32 #include <vector>
33 #include <algorithm>
34 #include <limits>
35 
36 namespace libMesh
37 {
38 
51 template <typename T>
52 class DistributedVector libmesh_final : public NumericVector<T>
53 {
54 public:
55 
59  explicit
60  DistributedVector (const Parallel::Communicator & comm,
61  const ParallelType = AUTOMATIC);
62 
66  explicit
67  DistributedVector (const Parallel::Communicator & comm,
68  const numeric_index_type n,
69  const ParallelType ptype = AUTOMATIC);
70 
75  DistributedVector (const Parallel::Communicator & comm,
76  const numeric_index_type n,
77  const numeric_index_type n_local,
78  const ParallelType ptype = AUTOMATIC);
79 
85  DistributedVector (const Parallel::Communicator & comm,
86  const numeric_index_type N,
87  const numeric_index_type n_local,
88  const std::vector<numeric_index_type> & ghost,
89  const ParallelType ptype = AUTOMATIC);
90 
95  ~DistributedVector ();
96 
97  virtual void close () libmesh_override;
98 
99  virtual void clear () libmesh_override;
100 
101  virtual void zero () libmesh_override;
102 
103  virtual UniquePtr<NumericVector<T>> zero_clone () const libmesh_override;
104 
105  virtual UniquePtr<NumericVector<T>> clone () const libmesh_override;
106 
107  virtual void init (const numeric_index_type N,
108  const numeric_index_type n_local,
109  const bool fast=false,
110  const ParallelType ptype=AUTOMATIC) libmesh_override;
111 
112  virtual void init (const numeric_index_type N,
113  const bool fast=false,
114  const ParallelType ptype=AUTOMATIC) libmesh_override;
115 
116  virtual void init (const numeric_index_type N,
117  const numeric_index_type n_local,
118  const std::vector<numeric_index_type> & ghost,
119  const bool fast = false,
120  const ParallelType = AUTOMATIC) libmesh_override;
121 
122  virtual void init (const NumericVector<T> & other,
123  const bool fast = false) libmesh_override;
124 
125  virtual NumericVector<T> & operator= (const T s) libmesh_override;
126 
127  virtual NumericVector<T> & operator= (const NumericVector<T> & v) libmesh_override;
128 
134  DistributedVector<T> & operator= (const DistributedVector<T> & v);
135 
136  virtual NumericVector<T> & operator= (const std::vector<T> & v) libmesh_override;
137 
138  virtual Real min () const libmesh_override;
139 
140  virtual Real max () const libmesh_override;
141 
142  virtual T sum() const libmesh_override;
143 
144  virtual Real l1_norm () const libmesh_override;
145 
146  virtual Real l2_norm () const libmesh_override;
147 
148  virtual Real linfty_norm () const libmesh_override;
149 
150  virtual numeric_index_type size () const libmesh_override;
151 
152  virtual numeric_index_type local_size() const libmesh_override;
153 
154  virtual numeric_index_type first_local_index() const libmesh_override;
155 
156  virtual numeric_index_type last_local_index() const libmesh_override;
157 
158  virtual T operator() (const numeric_index_type i) const libmesh_override;
159 
160  virtual NumericVector<T> & operator += (const NumericVector<T> & v) libmesh_override;
161 
162  virtual NumericVector<T> & operator -= (const NumericVector<T> & v) libmesh_override;
163 
164  virtual NumericVector<T> & operator /= (NumericVector<T> & v) libmesh_override;
165 
166  virtual void reciprocal() libmesh_override;
167 
168  virtual void conjugate() libmesh_override;
169 
170  virtual void set (const numeric_index_type i, const T value) libmesh_override;
171 
172  virtual void add (const numeric_index_type i, const T value) libmesh_override;
173 
174  virtual void add (const T s) libmesh_override;
175 
176  virtual void add (const NumericVector<T> & V) libmesh_override;
177 
178  virtual void add (const T a, const NumericVector<T> & v) libmesh_override;
179 
184  using NumericVector<T>::add_vector;
185 
186  virtual void add_vector (const NumericVector<T> &,
187  const SparseMatrix<T> &) libmesh_override
188  { libmesh_not_implemented(); }
189 
190  virtual void add_vector_transpose (const NumericVector<T> &,
191  const SparseMatrix<T> &) libmesh_override
192  { libmesh_not_implemented(); }
193 
194  virtual void scale (const T factor) libmesh_override;
195 
196  virtual void abs() libmesh_override;
197 
198  virtual T dot(const NumericVector<T> & V) const libmesh_override;
199 
200  virtual void localize (std::vector<T> & v_local) const libmesh_override;
201 
202  virtual void localize (NumericVector<T> & v_local) const libmesh_override;
203 
204  virtual void localize (NumericVector<T> & v_local,
205  const std::vector<numeric_index_type> & send_list) const libmesh_override;
206 
207  virtual void localize (std::vector<T> & v_local,
208  const std::vector<numeric_index_type> & indices) const libmesh_override;
209 
210  virtual void localize (const numeric_index_type first_local_idx,
211  const numeric_index_type last_local_idx,
212  const std::vector<numeric_index_type> & send_list) libmesh_override;
213 
214  virtual void localize_to_one (std::vector<T> & v_local,
215  const processor_id_type proc_id=0) const libmesh_override;
216 
217  virtual void pointwise_mult (const NumericVector<T> & vec1,
218  const NumericVector<T> & vec2) libmesh_override;
219 
220  virtual void swap (NumericVector<T> & v) libmesh_override;
221 
222 private:
223 
227  std::vector<T> _values;
228 
233 
238 
243 
248 };
249 
250 
251 //--------------------------------------------------------------------------
252 // DistributedVector inline methods
253 template <typename T>
254 inline
255 DistributedVector<T>::DistributedVector (const Parallel::Communicator & comm_in,
256  const ParallelType ptype) :
257  NumericVector<T>(comm_in, ptype),
258  _global_size (0),
259  _local_size (0),
260  _first_local_index(0),
261  _last_local_index (0)
262 {
263  this->_type = ptype;
264 }
265 
266 
267 
268 template <typename T>
269 inline
270 DistributedVector<T>::DistributedVector (const Parallel::Communicator & comm_in,
271  const numeric_index_type n,
272  const ParallelType ptype)
273  : NumericVector<T>(comm_in, ptype)
274 {
275  this->init(n, n, false, ptype);
276 }
277 
278 
279 
280 template <typename T>
281 inline
282 DistributedVector<T>::DistributedVector (const Parallel::Communicator & comm_in,
283  const numeric_index_type n,
284  const numeric_index_type n_local,
285  const ParallelType ptype)
286  : NumericVector<T>(comm_in, ptype)
287 {
288  this->init(n, n_local, false, ptype);
289 }
290 
291 
292 
293 template <typename T>
294 inline
295 DistributedVector<T>::DistributedVector (const Parallel::Communicator & comm_in,
296  const numeric_index_type n,
297  const numeric_index_type n_local,
298  const std::vector<numeric_index_type> & ghost,
299  const ParallelType ptype)
300  : NumericVector<T>(comm_in, ptype)
301 {
302  this->init(n, n_local, ghost, false, ptype);
303 }
304 
305 
306 
307 template <typename T>
308 inline
309 DistributedVector<T>::~DistributedVector ()
310 {
311  this->clear ();
312 }
313 
314 
315 
316 template <typename T>
317 inline
319  const numeric_index_type n_local,
320  const bool fast,
321  const ParallelType ptype)
322 {
323  // This function must be run on all processors at once
324  parallel_object_only();
325 
326  libmesh_assert_less_equal (n_local, n);
327 
328  if (ptype == AUTOMATIC)
329  {
330  if (n == n_local)
331  this->_type = SERIAL;
332  else
333  this->_type = PARALLEL;
334  }
335  else
336  this->_type = ptype;
337 
338  libmesh_assert ((this->_type==SERIAL && n==n_local) ||
339  this->_type==PARALLEL);
340 
341  // Clear the data structures if already initialized
342  if (this->initialized())
343  this->clear();
344 
345  // Initialize data structures
346  _values.resize(n_local);
347  _local_size = n_local;
348  _global_size = n;
349 
350  _first_local_index = 0;
351 
352 #ifdef LIBMESH_HAVE_MPI
353 
354  std::vector<numeric_index_type> local_sizes (this->n_processors(), 0);
355 
356  local_sizes[this->processor_id()] = n_local;
357 
358  this->comm().sum(local_sizes);
359 
360  // _first_local_index is the sum of _local_size
361  // for all processor ids less than ours
362  for (processor_id_type p=0; p!=this->processor_id(); p++)
363  _first_local_index += local_sizes[p];
364 
365 
366 # ifdef DEBUG
367  // Make sure all the local sizes sum up to the global
368  // size, otherwise there is big trouble!
369  numeric_index_type dbg_sum=0;
370 
371  for (processor_id_type p=0; p!=this->n_processors(); p++)
372  dbg_sum += local_sizes[p];
373 
374  libmesh_assert_equal_to (dbg_sum, n);
375 
376 # endif
377 
378 #else
379 
380  // No other options without MPI!
381  if (n != n_local)
382  libmesh_error_msg("ERROR: MPI is required for n != n_local!");
383 
384 #endif
385 
387 
388  // Set the initialized flag
389  this->_is_initialized = true;
390 
391  // Zero the components unless directed otherwise
392  if (!fast)
393  this->zero();
394 }
395 
396 
397 template <typename T>
398 inline
400  const numeric_index_type n_local,
401  const std::vector<numeric_index_type> & /*ghost*/,
402  const bool fast,
403  const ParallelType ptype)
404 {
405  // TODO: we shouldn't ignore the ghost sparsity pattern
406  this->init(n, n_local, fast, ptype);
407 }
408 
409 
410 
411 /* Default implementation for solver packages for which ghosted
412  vectors are not yet implemented. */
413 template <class T>
414 void DistributedVector<T>::init (const NumericVector<T> & other,
415  const bool fast)
416 {
417  this->init(other.size(),other.local_size(),fast,other.type());
418 }
419 
420 
421 
422 template <typename T>
423 inline
425  const bool fast,
426  const ParallelType ptype)
427 {
428  this->init(n,n,fast,ptype);
429 }
430 
431 
432 
433 template <typename T>
434 inline
435 void DistributedVector<T>::close ()
436 {
437  libmesh_assert (this->initialized());
438 
439  this->_is_closed = true;
440 }
441 
442 
443 
444 template <typename T>
445 inline
446 void DistributedVector<T>::clear ()
447 {
448  _values.clear();
449 
450  _global_size =
451  _local_size =
453  _last_local_index = 0;
454 
455 
456  this->_is_closed = this->_is_initialized = false;
457 }
458 
459 
460 
461 template <typename T>
462 inline
464 {
465  libmesh_assert (this->initialized());
466  libmesh_assert_equal_to (_values.size(), _local_size);
467  libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
468 
469  std::fill (_values.begin(),
470  _values.end(),
471  0.);
472 }
473 
474 
475 
476 template <typename T>
477 inline
478 UniquePtr<NumericVector<T>> DistributedVector<T>::zero_clone () const
479 {
480  NumericVector<T> * cloned_vector = new DistributedVector<T>(this->comm());
481  cloned_vector->init(*this);
482  return UniquePtr<NumericVector<T>>(cloned_vector);
483 }
484 
485 
486 
487 template <typename T>
488 inline
489 UniquePtr<NumericVector<T>> DistributedVector<T>::clone () const
490 {
491  NumericVector<T> * cloned_vector = new DistributedVector<T>(this->comm());
492  cloned_vector->init(*this, true);
493  *cloned_vector = *this;
494  return UniquePtr<NumericVector<T>>(cloned_vector);
495 }
496 
497 
498 
499 template <typename T>
500 inline
501 numeric_index_type DistributedVector<T>::size () const
502 {
503  libmesh_assert (this->initialized());
504  libmesh_assert_equal_to (_values.size(), _local_size);
505  libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
506 
507  return _global_size;
508 }
509 
510 
511 
512 template <typename T>
513 inline
514 numeric_index_type DistributedVector<T>::local_size () const
515 {
516  libmesh_assert (this->initialized());
517  libmesh_assert_equal_to (_values.size(), _local_size);
518  libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
519 
520  return _local_size;
521 }
522 
523 
524 
525 template <typename T>
526 inline
527 numeric_index_type DistributedVector<T>::first_local_index () const
528 {
529  libmesh_assert (this->initialized());
530  libmesh_assert_equal_to (_values.size(), _local_size);
531  libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
532 
533  return _first_local_index;
534 }
535 
536 
537 
538 template <typename T>
539 inline
540 numeric_index_type DistributedVector<T>::last_local_index () const
541 {
542  libmesh_assert (this->initialized());
543  libmesh_assert_equal_to (_values.size(), _local_size);
544  libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
545 
546  return _last_local_index;
547 }
548 
549 
550 
551 template <typename T>
552 inline
553 T DistributedVector<T>::operator() (const numeric_index_type i) const
554 {
555  libmesh_assert (this->initialized());
556  libmesh_assert_equal_to (_values.size(), _local_size);
557  libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
558  libmesh_assert ( ((i >= first_local_index()) &&
559  (i < last_local_index())) );
560 
561  return _values[i - _first_local_index];
562 }
563 
564 
565 
566 template <typename T>
567 inline
568 void DistributedVector<T>::set (const numeric_index_type i, const T value)
569 {
570  libmesh_assert (this->initialized());
571  libmesh_assert_equal_to (_values.size(), _local_size);
572  libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
573  libmesh_assert_less (i, size());
574  libmesh_assert_less (i-first_local_index(), local_size());
575 
577 }
578 
579 
580 
581 template <typename T>
582 inline
583 void DistributedVector<T>::add (const numeric_index_type i, const T value)
584 {
585  libmesh_assert (this->initialized());
586  libmesh_assert_equal_to (_values.size(), _local_size);
587  libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
588  libmesh_assert_less (i, size());
589  libmesh_assert_less (i-first_local_index(), local_size());
590 
592 }
593 
594 
595 
596 template <typename T>
597 inline
599 {
600  // This function must be run on all processors at once
601  parallel_object_only();
602 
603  libmesh_assert (this->initialized());
604  libmesh_assert_equal_to (_values.size(), _local_size);
605  libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
606 
607  Real local_min = _values.size() ?
609  for (numeric_index_type i = 1; i < _values.size(); ++i)
610  local_min = std::min(libmesh_real(_values[i]), local_min);
611 
612  this->comm().min(local_min);
613 
614  return local_min;
615 }
616 
617 
618 
619 template <typename T>
620 inline
622 {
623  // This function must be run on all processors at once
624  parallel_object_only();
625 
626  libmesh_assert (this->initialized());
627  libmesh_assert_equal_to (_values.size(), _local_size);
628  libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
629 
630  Real local_max = _values.size() ?
632  for (numeric_index_type i = 1; i < _values.size(); ++i)
633  local_max = std::max(libmesh_real(_values[i]), local_max);
634 
635  this->comm().max(local_max);
636 
637  return local_max;
638 }
639 
640 
641 template <typename T>
642 inline
644 {
645  DistributedVector<T> & v = cast_ref<DistributedVector<T> &>(other);
646 
647  std::swap(_global_size, v._global_size);
648  std::swap(_local_size, v._local_size);
649  std::swap(_first_local_index, v._first_local_index);
650  std::swap(_last_local_index, v._last_local_index);
651 
652  // This should be O(1) with any reasonable STL implementation
653  std::swap(_values, v._values);
654 }
655 
656 } // namespace libMesh
657 
658 
659 #endif // LIBMESH_DISTRIBUTED_VECTOR_H
T libmesh_real(T a)
virtual numeric_index_type n() const libmesh_override
double abs(double a)
Encapsulates the MPI_Comm object.
Definition: parallel.h:657
virtual numeric_index_type size() const =0
numeric_index_type _first_local_index
The first component stored locally.
virtual void init(const numeric_index_type, const numeric_index_type, const bool=false, const ParallelType=AUTOMATIC)=0
Change the dimension of the vector to N.
void max(T &r) const
Take a local variable and replace it with the maximum of it&#39;s values on all processors.
void scale(MeshBase &mesh, const Real xs, const Real ys=0., const Real zs=0.)
Scales the mesh.
void min(T &r) const
Take a local variable and replace it with the minimum of it&#39;s values on all processors.
numeric_index_type _global_size
The global vector size.
processor_id_type n_processors() const
uint8_t processor_id_type
Definition: id_types.h:99
Numeric vector.
Definition: dof_map.h:66
virtual numeric_index_type local_size() const libmesh_override
bool _is_initialized
true once init() has been called.
virtual bool initialized() const
The libMesh namespace provides an interface to certain functionality in the library.
void sum(T &r, const Communicator &comm=Communicator_World)
const Number zero
.
Definition: libmesh.h:178
virtual numeric_index_type size() const libmesh_override
long double max(long double a, double b)
libmesh_assert(j)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
Generic sparse matrix.
Definition: dof_map.h:65
std::vector< T > _values
Actual vector datatype to hold vector entries.
numeric_index_type _last_local_index
The last component (+1) stored locally.
ParallelType
Defines an enum for parallel data structure types.
dof_id_type numeric_index_type
Definition: id_types.h:92
void init(triangulateio &t)
Initializes the fields of t to NULL/0 as necessary.
virtual void init() libmesh_override
Initialize this matrix using the sparsity structure computed by dof_map.
ParallelType _type
Type of vector.
virtual numeric_index_type local_size() const =0
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void clear() libmesh_override
Restores the NumericVector<T> to a pristine state.
X_input swap(X_system)
const Parallel::Communicator & comm() const
void swap(Iterator &lhs, Iterator &rhs)
swap, used to implement op=
numeric_index_type _local_size
The local vector size.
virtual void zero() libmesh_override
Set all entries to zero.
virtual numeric_index_type last_local_index() const libmesh_override
virtual void add_vector_transpose(const NumericVector< T > &, const SparseMatrix< T > &) libmesh_override
Computes , i.e.
static const bool value
Definition: xdr_io.C:108
ParallelType type() const
void sum(T &r) const
Take a local variable and replace it with the sum of it&#39;s values on all processors.
virtual numeric_index_type first_local_index() const libmesh_override
long double min(long double a, double b)
bool _is_closed
Flag which tracks whether the vector&#39;s values are consistent on all processors after insertion or add...
processor_id_type processor_id() const