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"
61 const ParallelType = AUTOMATIC);
69 const ParallelType ptype = AUTOMATIC);
77 const ParallelType ptype = AUTOMATIC);
86 const std::vector<numeric_index_type> & ghost,
87 const ParallelType ptype = AUTOMATIC);
104 virtual void close ()
override;
106 virtual void clear ()
override;
108 virtual void zero ()
override;
110 virtual std::unique_ptr<NumericVector<T>>
zero_clone ()
const override;
112 virtual std::unique_ptr<NumericVector<T>>
clone ()
const override;
116 const bool fast=
false,
117 const ParallelType ptype=AUTOMATIC)
override;
120 const bool fast=
false,
121 const ParallelType ptype=AUTOMATIC)
override;
125 const std::vector<numeric_index_type> & ghost,
126 const bool fast =
false,
127 const ParallelType = AUTOMATIC)
override;
130 const bool fast =
false)
override;
138 virtual Real min ()
const override;
140 virtual Real max ()
const override;
142 virtual T
sum()
const override;
176 virtual void add (
const T s)
override;
190 { libmesh_not_implemented(); }
194 { libmesh_not_implemented(); }
196 virtual void scale (
const T factor)
override;
198 virtual void abs()
override;
202 virtual void localize (std::vector<T> & v_local)
const override;
207 const std::vector<numeric_index_type> & send_list)
const override;
209 virtual void localize (std::vector<T> & v_local,
210 const std::vector<numeric_index_type> & indices)
const override;
214 const std::vector<numeric_index_type> & send_list)
override;
260 template <
typename T>
263 const ParallelType ptype) :
267 _first_local_index(0),
268 _last_local_index (0)
275 template <
typename T>
279 const ParallelType ptype)
282 this->
init(n, n,
false, ptype);
287 template <
typename T>
292 const ParallelType ptype)
295 this->
init(n, n_local,
false, ptype);
300 template <
typename T>
305 const std::vector<numeric_index_type> & ghost,
306 const ParallelType ptype)
309 this->
init(n, n_local, ghost,
false, ptype);
314 template <
typename T>
319 const ParallelType ptype)
322 parallel_object_only();
324 libmesh_assert_less_equal (n_local, n);
326 if (ptype == AUTOMATIC)
329 this->_type = SERIAL;
331 this->_type = PARALLEL;
333 else if (ptype == GHOSTED &&
335 this->_type = SERIAL;
339 libmesh_assert ((this->_type==SERIAL && n==n_local) ||
340 this->_type==PARALLEL);
347 _values.resize(n_local);
348 _local_size = n_local;
351 _first_local_index = 0;
353 #ifdef LIBMESH_HAVE_MPI
355 std::vector<numeric_index_type> local_sizes (this->n_processors(), 0);
357 local_sizes[this->processor_id()] = n_local;
359 this->comm().sum(local_sizes);
363 for (
auto p : make_range(this->processor_id()))
364 _first_local_index += local_sizes[p];
372 for (
auto p : make_range(this->n_processors()))
373 dbg_sum += local_sizes[p];
375 libmesh_assert_equal_to (dbg_sum, n);
382 libmesh_error_msg_if(n != n_local,
"ERROR: MPI is required for n != n_local!");
386 _last_local_index = _first_local_index + n_local;
397 template <
typename T>
401 const std::vector<numeric_index_type> & ,
403 const ParallelType ptype)
406 this->init(n, n_local, fast, ptype);
422 template <
typename T>
426 const ParallelType ptype)
428 this->init(n,n,fast,ptype);
433 template <
typename T>
439 this->_is_closed =
true;
444 template <
typename T>
453 _last_local_index = 0;
461 template <
typename T>
466 libmesh_assert_equal_to (_values.size(), _local_size);
467 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
469 std::fill (_values.begin(),
476 template <
typename T>
481 cloned_vector->
init(*
this);
482 return std::unique_ptr<NumericVector<T>>(cloned_vector);
487 template <
typename T>
492 cloned_vector->
init(*
this,
true);
493 *cloned_vector = *
this;
494 return std::unique_ptr<NumericVector<T>>(cloned_vector);
499 template <
typename T>
504 libmesh_assert_equal_to (_values.size(), _local_size);
505 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
512 template <
typename T>
517 libmesh_assert_equal_to (_values.size(), _local_size);
518 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
525 template <
typename T>
530 libmesh_assert_equal_to (_values.size(), _local_size);
531 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
533 return _first_local_index;
538 template <
typename T>
543 libmesh_assert_equal_to (_values.size(), _local_size);
544 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
546 return _last_local_index;
551 template <
typename T>
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())) );
561 return _values[i - _first_local_index];
566 template <
typename T>
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());
576 std::scoped_lock lock(this->_numeric_vector_mutex);
577 _values[i - _first_local_index] = value;
580 this->_is_closed =
false;
586 template <
typename T>
591 libmesh_assert_equal_to (_values.size(), _local_size);
592 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
593 libmesh_assert_less (i, size());
594 libmesh_assert_less (i-first_local_index(), local_size());
596 std::scoped_lock lock(this->_numeric_vector_mutex);
597 _values[i - _first_local_index] += value;
600 this->_is_closed =
false;
605 template <
typename T>
610 parallel_object_only();
613 libmesh_assert_equal_to (_values.size(), _local_size);
614 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
616 Real local_min = std::numeric_limits<Real>::max();
617 for (
auto v : _values)
620 this->comm().min(local_min);
627 template <
typename T>
632 parallel_object_only();
635 libmesh_assert_equal_to (_values.size(), _local_size);
636 libmesh_assert_equal_to ((_last_local_index - _first_local_index), _local_size);
638 Real local_max = -std::numeric_limits<Real>::max();
639 for (
auto v : _values)
642 this->comm().max(local_max);
648 template <
typename T>
663 template <
typename T>
668 return std::numeric_limits<typename std::vector<T>::size_type>::max();
674 #endif // LIBMESH_DISTRIBUTED_VECTOR_H
virtual Real max() const override
获取向量中的最大值,或者在复数情况下获取最大的实部。
virtual T dot(const NumericVector< T > &V) const override
计算 ,即 (*this) 与向量 v 的点积。
numeric_index_type _last_local_index
本地存储的最后一个分量(+1)。
virtual void scale(const T factor) override
缩放向量的每个元素。
virtual void add_vector(const NumericVector< T > &, const SparseMatrix< T > &) override
计算 , 即将 SparseMatrix A 和 NumericVector v 的乘积添加到 this。
virtual numeric_index_type last_local_index() const override
获取实际存储在该处理器上的最后一个向量元素的索引+1。
virtual void close() override
调用 NumericVector 的内部组装函数,确保值在处理器之间一致。
virtual numeric_index_type size() const =0
获取向量的大小。
virtual void set(const numeric_index_type i, const T value) override
设置 v(i) = value 。 请注意,此方法的库实现是线程安全的, 例如,将在写入向量之前锁定 _numeric_vector_mutex 。
virtual NumericVector< T > & operator-=(const NumericVector< T > &v) override
将 v 从 *this 减去, 。等价于 u.add(-1, v)。
DistributedVector & operator=(const DistributedVector &)
复制赋值运算符。我们不能默认实现它(尽管它本质上实现了默认行为), 因为生成的默认版本尝试自动调用基类(NumericVector)的复制赋值运算符, 而选择使其成为纯虚拟函数,出于其他设计原因。 ...
virtual NumericVector< T > & operator/=(const NumericVector< T > &v) override
计算此向量条目与另一个向量的分量除法, 。
virtual void add(const numeric_index_type i, const T value) override
将 value 添加到由 i 指定的向量条目。 请注意,此方法的库实现是线程安全的, 例如,将在写入向量之前锁定 _numeric_vector_mutex 。
virtual Real linfty_norm() const override
获取向量的 -范数,即向量的最大绝对值。
numeric_index_type _global_size
全局向量大小。
virtual std::unique_ptr< NumericVector< T > > zero_clone() const override
返回一个智能指针,指向具有相同类型、大小和分区的此向量的副本,但所有条目都为零。
virtual numeric_index_type local_size() const override
获取向量的本地大小,即 index_stop - index_start。
virtual void init(const numeric_index_type n, const numeric_index_type n_local, const bool fast=false, const ParallelType ptype=AUTOMATIC)=0
更改向量的维度为 n 。如果可能的话 ,该向量的保留内存保持不变。 如果 n==0 ,所有内存都将被释放。因此,如果要调整向量的大小并释放不需要的内存, 必须首先调用 init(0) ,然后调用 ini...
virtual NumericVector< T > & operator+=(const NumericVector< T > &v) override
将向量加上 v , 。等价于 u.add(1, v)。
uint8_t processor_id_type
这是一个通用的稀疏矩阵类。该类包含了必须在派生类中覆盖的纯虚拟成员。 使用一个公共的基类允许从不同的求解器包中以不同的格式统一访问稀疏矩阵。
virtual void reciprocal() override
计算每个向量条目的分量倒数, 。
virtual void zero() override
将所有条目设置为零。等同于 v = 0,但更明显且更快。
virtual Real l1_norm() const override
获取向量的 -范数,即条目的绝对值之和。
该类提供了一个简单的并行分布式向量数据类型, 专门用于 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
获取向量的大小。
virtual std::size_t max_allowed_id() const override
返回 NumericVector 可以包含的最大条目数(在所有处理器上)。
virtual void pointwise_divide(const NumericVector< T > &vec1, const NumericVector< T > &vec2) override
计算该向量与另一个向量的逐点除法。
virtual std::unique_ptr< NumericVector< T > > clone() const override
返回一个包装了此向量副本的智能指针。
virtual void pointwise_mult(const NumericVector< T > &vec1, const NumericVector< T > &vec2) override
比较该向量与另一个向量的全局相对差异。
virtual void conjugate() override
反转向量中每个条目的虚部。
virtual void localize_to_one(std::vector< T > &v_local, const processor_id_type proc_id=0) const override
在处理器 proc_id 上创建全局向量的本地副本。 默认情况下,数据发送到处理器 0。此方法对于从一个处理器输出数据非常有用。
numeric_index_type _local_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
计算此向量条目与另一个向量的条目之间的分量乘法, 。
virtual void localize(std::vector< T > &v_local) const override
创建全局向量的副本并存储在本地向量 v_local 中。
virtual Real l2_norm() const override
获取向量的 -范数,即条目平方和的平方根。
virtual numeric_index_type local_size() const =0
获取向量的本地大小,即 index_stop - index_start。
virtual void abs() override
设置 ,对向量中的每个条目进行绝对值操作。
ParallelType type() const
获取向量的类型。
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
更改向量的维度为 n 。如果可能的话 ,该向量的保留内存保持不变。 如果 n==0 ,所有内存都将被释放。因此,如果要调整向量的大小并释放不需要的内存, 必须首先调用 init(0) ,然后调用 ini...
virtual void clear() override
将 NumericVector<T> 恢复到原始状态。
std::vector< T > _values
实际的向量数据类型,用于保存向量条目。
virtual Real min() const override
获取向量中的最小值,或者在复数情况下获取最小的实部。
virtual T sum() const override
获取向量中所有值的总和。
virtual void swap(NumericVector< T > &v) override
交换该向量的内容与向量 v 的内容。子类应提供足够的间接性以使此操作成为 O(1) 的头部交换操作。
virtual void add_vector_transpose(const NumericVector< T > &, const SparseMatrix< T > &) override
计算 , 即将矩阵 A 的转置与 NumericVector v 的乘积添加到 this。
DistributedVector(const Parallel::Communicator &comm, const ParallelType=AUTOMATIC)
虚构造函数。维度为 0。
numeric_index_type _first_local_index
本地存储的第一个分量。
virtual T operator()(const numeric_index_type i) const override
获取向量的第 i 个条目的副本。