20 #include "libmesh/libmesh_common.h" 24 #ifndef LIBMESH_DISTRIBUTED_VECTOR_H 25 #define LIBMESH_DISTRIBUTED_VECTOR_H 28 #include "libmesh/int_range.h" 29 #include "libmesh/numeric_vector.h" 30 #include "libmesh/parallel.h" 90 const std::vector<numeric_index_type> & ghost,
111 virtual void close ()
override;
113 virtual void clear ()
override;
115 virtual void zero ()
override;
117 virtual std::unique_ptr<NumericVector<T>>
zero_clone ()
const override;
119 virtual std::unique_ptr<NumericVector<T>>
clone ()
const override;
123 const bool fast=
false,
127 const bool fast=
false,
132 const std::vector<numeric_index_type> & ghost,
133 const bool fast =
false,
137 const bool fast =
false)
override;
145 virtual Real min ()
const override;
147 virtual Real max ()
const override;
149 virtual T
sum()
const override;
183 virtual void add (
const T s)
override;
197 { libmesh_not_implemented(); }
201 { libmesh_not_implemented(); }
203 virtual void scale (
const T factor)
override;
205 virtual void abs()
override;
209 virtual void localize (std::vector<T> & v_local)
const override;
214 const std::vector<numeric_index_type> & send_list)
const override;
216 virtual void localize (std::vector<T> & v_local,
217 const std::vector<numeric_index_type> & indices)
const override;
221 const std::vector<numeric_index_type> & send_list)
override;
267 template <
typename T>
274 _first_local_index(0),
275 _last_local_index (0)
282 template <
typename T>
289 this->
init(n, n,
false, ptype);
294 template <
typename T>
302 this->
init(n, n_local,
false, ptype);
307 template <
typename T>
312 const std::vector<numeric_index_type> & ghost,
316 this->
init(n, n_local, ghost,
false, ptype);
321 template <
typename T>
329 parallel_object_only();
331 libmesh_assert_less_equal (n_local, n);
354 _values.resize(n_local);
355 _local_size = n_local;
358 _first_local_index = 0;
360 #ifdef LIBMESH_HAVE_MPI 362 std::vector<numeric_index_type> local_sizes (this->n_processors(), 0);
364 local_sizes[this->processor_id()] = n_local;
366 this->comm().sum(local_sizes);
370 for (
auto p :
make_range(this->processor_id()))
371 _first_local_index += local_sizes[p];
379 for (
auto p :
make_range(this->n_processors()))
380 dbg_sum += local_sizes[p];
382 libmesh_assert_equal_to (dbg_sum, n);
389 libmesh_error_msg_if(n != n_local,
"ERROR: MPI is required for n != n_local!");
393 _last_local_index = _first_local_index + n_local;
404 template <
typename T>
408 const std::vector<numeric_index_type> & ,
413 this->
init(n, n_local, fast, ptype);
429 template <
typename T>
435 this->
init(n,n,fast,ptype);
440 template <
typename T>
446 this->_is_closed =
true;
451 template <
typename T>
460 _last_local_index = 0;
468 template <
typename T>
473 libmesh_assert_equal_to (_values.size(), _local_size);
474 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
476 std::fill (_values.begin(),
483 template <
typename T>
488 cloned_vector->
init(*
this);
489 return std::unique_ptr<NumericVector<T>>(cloned_vector);
494 template <
typename T>
499 cloned_vector->
init(*
this,
true);
500 *cloned_vector = *
this;
501 return std::unique_ptr<NumericVector<T>>(cloned_vector);
506 template <
typename T>
511 libmesh_assert_equal_to (_values.size(), _local_size);
512 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
519 template <
typename T>
524 libmesh_assert_equal_to (_values.size(), _local_size);
525 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
532 template <
typename T>
537 libmesh_assert_equal_to (_values.size(), _local_size);
538 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
540 return _first_local_index;
545 template <
typename T>
550 libmesh_assert_equal_to (_values.size(), _local_size);
551 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
553 return _last_local_index;
558 template <
typename T>
563 libmesh_assert_equal_to (_values.size(), _local_size);
564 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
566 (i < last_local_index())) );
568 return _values[i - _first_local_index];
573 template <
typename T>
578 libmesh_assert_equal_to (_values.size(), _local_size);
579 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
580 libmesh_assert_less (i, size());
581 libmesh_assert_less (i-first_local_index(), local_size());
583 std::scoped_lock lock(this->_numeric_vector_mutex);
584 _values[i - _first_local_index] =
value;
587 this->_is_closed =
false;
593 template <
typename T>
598 libmesh_assert_equal_to (_values.size(), _local_size);
599 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
600 libmesh_assert_less (i, size());
601 libmesh_assert_less (i-first_local_index(), local_size());
603 std::scoped_lock lock(this->_numeric_vector_mutex);
604 _values[i - _first_local_index] +=
value;
607 this->_is_closed =
false;
612 template <
typename T>
617 parallel_object_only();
620 libmesh_assert_equal_to (_values.size(), _local_size);
621 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
623 Real local_min = std::numeric_limits<Real>::max();
624 for (
auto v : _values)
627 this->comm().min(local_min);
634 template <
typename T>
639 parallel_object_only();
642 libmesh_assert_equal_to (_values.size(), _local_size);
643 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
645 Real local_max = -std::numeric_limits<Real>::max();
646 for (
auto v : _values)
649 this->comm().max(local_max);
655 template <
typename T>
672 template <
typename T>
677 return std::numeric_limits<typename std::vector<T>::size_type>::max();
683 #endif // LIBMESH_DISTRIBUTED_VECTOR_H
virtual Real max() const override
virtual T dot(const NumericVector< T > &V) const override
numeric_index_type _last_local_index
The last component (+1) stored locally.
virtual void scale(const T factor) override
Scale each element of the vector by the given factor.
virtual void add_vector(const NumericVector< T > &, const SparseMatrix< T > &) override
Computes , i.e.
virtual numeric_index_type last_local_index() const override
virtual void close() override
Calls the NumericVector's internal assembly routines, ensuring that the values are consistent across ...
virtual numeric_index_type size() const =0
virtual void set(const numeric_index_type i, const T value) override
Sets v(i) = value.
virtual NumericVector< T > & operator-=(const NumericVector< T > &v) override
Subtracts v from *this, .
DistributedVector & operator=(const DistributedVector &)
Copy assignment operator.
Provides a uniform interface to vector storage schemes for different linear algebra libraries...
const Parallel::Communicator & comm() const
virtual NumericVector< T > & operator/=(const NumericVector< T > &v) override
Computes the component-wise division of this vector's entries by another's, .
virtual void add(const numeric_index_type i, const T value) override
Adds value to the vector entry specified by i.
bool _is_initialized
true once init() has been called.
virtual Real linfty_norm() const override
The libMesh namespace provides an interface to certain functionality in the library.
numeric_index_type _global_size
The global vector size.
virtual std::unique_ptr< NumericVector< T > > zero_clone() const override
virtual numeric_index_type local_size() const override
virtual void init(const numeric_index_type n, const numeric_index_type n_local, const bool fast=false, const ParallelType ptype=AUTOMATIC)=0
Change the dimension of the vector to n.
virtual NumericVector< T > & operator+=(const NumericVector< T > &v) override
Adds v to *this, .
uint8_t processor_id_type
virtual void reciprocal() override
Computes the component-wise reciprocal, .
virtual void zero() override
Set all entries to zero.
virtual Real l1_norm() const override
This class provides a simple parallel, distributed vector datatype which is specific to libmesh...
dof_id_type numeric_index_type
bool _is_initialized
Flag that tells if init() has been called.
virtual ~DistributedVector()=default
virtual numeric_index_type size() const override
void init(triangulateio &t)
Initializes the fields of t to nullptr/0 as necessary.
virtual std::size_t max_allowed_id() const override
virtual void pointwise_divide(const NumericVector< T > &vec1, const NumericVector< T > &vec2) override
Computes (summation not implied) i.e.
ParallelType _type
Type of vector.
virtual std::unique_ptr< NumericVector< T > > clone() const override
virtual void pointwise_mult(const NumericVector< T > &vec1, const NumericVector< T > &vec2) override
Computes (summation not implied) i.e.
virtual void conjugate() override
Negates the imaginary component of each entry in the vector.
ParallelType type() const
virtual void localize_to_one(std::vector< T > &v_local, const processor_id_type proc_id=0) const override
Creates a local copy of the global vector in v_local only on processor proc_id.
numeric_index_type _local_size
The local vector size.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual numeric_index_type first_local_index() const override
virtual NumericVector< T > & operator*=(const NumericVector< T > &v) override
Computes the component-wise multiplication of this vector's entries by another's, ...
virtual void localize(std::vector< T > &v_local) const override
Creates a copy of the global vector in the local vector v_local.
virtual Real l2_norm() const override
virtual numeric_index_type local_size() const =0
virtual void abs() override
Sets for each entry in the vector.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
bool initialized()
Checks that library initialization has been done.
virtual void init(const numeric_index_type N, const numeric_index_type n_local, const bool fast=false, const ParallelType ptype=AUTOMATIC) override
Change the dimension of the vector to n.
virtual void clear() override
Restores the NumericVector<T> to a pristine state.
std::vector< T > _values
Actual vector datatype to hold vector entries.
virtual Real min() const override
virtual T sum() const override
virtual void swap(NumericVector< T > &v) override
Swaps the contents of this with v.
virtual void add_vector_transpose(const NumericVector< T > &, const SparseMatrix< T > &) override
Computes , i.e.
DistributedVector(const Parallel::Communicator &comm, const ParallelType=AUTOMATIC)
Dummy-Constructor.
numeric_index_type _first_local_index
The first component stored locally.
ParallelType
Defines an enum for parallel data structure types.
virtual T operator()(const numeric_index_type i) const override