Reference documentation for deal.II version GIT c1ddcf411d 2023-11-30 16:30:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
la_parallel_block_vector.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 2023 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_la_parallel_block_vector_h
17 #define dealii_la_parallel_block_vector_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 
29 
30 #include <cstdio>
31 #include <vector>
32 
34 
35 
36 // Forward declarations
37 #ifndef DOXYGEN
38 # ifdef DEAL_II_WITH_PETSC
39 namespace PETScWrappers
40 {
41  namespace MPI
42  {
43  class BlockVector;
44  }
45 } // namespace PETScWrappers
46 # endif
47 
48 # ifdef DEAL_II_WITH_TRILINOS
49 namespace TrilinosWrappers
50 {
51  namespace MPI
52  {
53  class BlockVector;
54  }
55 } // namespace TrilinosWrappers
56 # endif
57 #endif
58 
59 namespace LinearAlgebra
60 {
61  namespace distributed
62  {
84  template <typename Number>
85  class BlockVector : public BlockVectorBase<Vector<Number>>
86  {
87  public:
100  static constexpr unsigned int communication_block_size = 20;
101 
106 
110  using BlockType = typename BaseClass::BlockType;
111 
116  using real_type = typename BaseClass::real_type;
117  using pointer = typename BaseClass::pointer;
119  using reference = typename BaseClass::reference;
121  using size_type = typename BaseClass::size_type;
122  using iterator = typename BaseClass::iterator;
124 
140  explicit BlockVector(const size_type num_blocks = 0,
141  const size_type block_size = 0);
142 
148 
155  template <typename OtherNumber>
157 
162  BlockVector(const std::vector<size_type> &block_sizes);
163 
168  BlockVector(const std::vector<IndexSet> &local_ranges,
169  const std::vector<IndexSet> &ghost_indices,
170  const MPI_Comm communicator);
171 
175  BlockVector(const std::vector<IndexSet> &local_ranges,
176  const MPI_Comm communicator);
177 
190  const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
191  &partitioners,
192  const MPI_Comm &comm_sm = MPI_COMM_SELF);
193 
202  ~BlockVector() = default;
203 
208  BlockVector &
210 
215  BlockVector &
217 
222  template <class Number2>
223  BlockVector &
225 
229  BlockVector &
231 
232 #ifdef DEAL_II_WITH_PETSC
242 #endif
243 
244 #ifdef DEAL_II_WITH_TRILINOS
255 #endif
256 
270  void
271  reinit(const size_type num_blocks,
272  const size_type block_size = 0,
273  const bool omit_zeroing_entries = false);
274 
294  void
295  reinit(const std::vector<size_type> &block_sizes,
296  const bool omit_zeroing_entries = false);
297 
312  template <typename Number2>
313  void
315  const bool omit_zeroing_entries = false);
316 
334  void
335  reinit(const std::vector<IndexSet> &local_ranges,
336  const std::vector<IndexSet> &ghost_indices,
337  const MPI_Comm communicator);
338 
342  void
343  reinit(const std::vector<IndexSet> &local_ranges,
344  const MPI_Comm communicator);
345 
359  void
361  const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
362  &partitioners,
363  const MPI_Comm &comm_sm = MPI_COMM_SELF);
364 
371  void
373  const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
374  &partitioners,
375  const bool make_ghosted,
376  const MPI_Comm &comm_sm = MPI_COMM_SELF);
377 
401  void
403 
412  void
414 
423  void
425 
429  bool
431 
435  void
436  set_ghost_state(const bool ghosted) const;
437 
452  template <typename Number2>
453  void
455 
460  template <typename OtherNumber>
461  void
462  add(const std::vector<size_type> &indices,
463  const ::Vector<OtherNumber> &values);
464 
469  void
470  sadd(const Number s, const BlockVector<Number> &V);
471 
477  bool
478  all_zero() const;
479 
483  Number
484  mean_value() const;
485 
490  real_type
491  lp_norm(const real_type p) const;
492 
510  void
523  operator*=(const Number factor);
524 
529  operator/=(const Number factor);
530 
536 
542 
551  void
554  const VectorOperation::values operation,
555  const std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
556  &communication_pattern = {});
557 
561  DEAL_II_DEPRECATED void
563  VectorOperation::values operation,
564  std::shared_ptr<const Utilities::MPI::CommunicationPatternBase>
565  communication_pattern = {})
566  {
567  import_elements(V, operation, communication_pattern);
568  }
569 
573  Number
575 
592  template <typename FullMatrixType>
593  void
595  const BlockVector<Number> &V,
596  const bool symmetric = false) const;
597 
614  template <typename FullMatrixType>
615  Number
617  const BlockVector<Number> &V,
618  const bool symmetric = false) const;
619 
629  template <typename FullMatrixType>
630  void
632  const FullMatrixType &matrix,
633  const Number s = Number(0.),
634  const Number b = Number(1.)) const;
635 
639  void
640  add(const Number a);
641 
645  void
646  add(const Number a, const BlockVector<Number> &V);
647 
651  void
652  add(const Number a,
653  const BlockVector<Number> &V,
654  const Number b,
655  const BlockVector<Number> &W);
656 
661  void
662  add(const std::vector<size_type> &indices,
663  const std::vector<Number> &values);
664 
669  void
670  sadd(const Number s, const Number a, const BlockVector<Number> &V);
671 
677  void
678  scale(const BlockVector<Number> &scaling_factors);
679 
683  void
684  equ(const Number a, const BlockVector<Number> &V);
685 
690  real_type
691  l1_norm() const;
692 
697  real_type
698  l2_norm() const;
699 
703  real_type
704  norm_sqr() const;
705 
710  real_type
711  linfty_norm() const;
712 
733  Number
734  add_and_dot(const Number a,
735  const BlockVector<Number> &V,
736  const BlockVector<Number> &W);
737 
742  virtual size_type
743  size() const override;
744 
756  ::IndexSet
758 
762  void
763  print(std::ostream &out,
764  const unsigned int precision = 3,
765  const bool scientific = true,
766  const bool across = true) const;
767 
771  std::size_t
786 
792  };
793 
796  } // end of namespace distributed
797 
798 } // end of namespace LinearAlgebra
799 
800 
808 template <typename Number>
809 inline void
812 {
813  u.swap(v);
814 }
815 
816 
821 template <typename Number>
823  : std::false_type
824 {};
825 
826 
828 
829 #ifdef DEAL_II_MSVC
830 # include <deal.II/lac/la_parallel_block_vector.templates.h>
831 #endif
832 
833 #endif
::internal::BlockVectorIterators::Iterator< BlockVectorBase, false > iterator
typename BlockType::const_reference const_reference
typename BlockType::value_type value_type
::internal::BlockVectorIterators::Iterator< BlockVectorBase, true > const_iterator
typename BlockType::reference reference
typename BlockType::real_type real_type
void swap(LinearAlgebra::distributed::BlockVector< Number > &u, LinearAlgebra::distributed::BlockVector< Number > &v)
BlockVector & operator=(const BlockVector &V)
BlockVector & operator=(const Vector< Number > &V)
BlockVector(const std::vector< std::shared_ptr< const Utilities::MPI::Partitioner >> &partitioners, const MPI_Comm &comm_sm=MPI_COMM_SELF)
void swap(BlockVector< Number > &v)
Number operator*(const BlockVector< Number > &V) const
BlockVector(const size_type num_blocks=0, const size_type block_size=0)
void add(const std::vector< size_type > &indices, const std::vector< Number > &values)
void reinit(const std::vector< std::shared_ptr< const Utilities::MPI::Partitioner >> &partitioners, const MPI_Comm &comm_sm=MPI_COMM_SELF)
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
BlockVector< Number > & operator=(const TrilinosWrappers::MPI::BlockVector &trilinos_vec)
void reinit(const std::vector< IndexSet > &local_ranges, const MPI_Comm communicator)
void sadd(const Number s, const Number a, const BlockVector< Number > &V)
static constexpr unsigned int communication_block_size
BlockVector & operator=(const value_type s)
void add(const Number a, const BlockVector< Number > &V, const Number b, const BlockVector< Number > &W)
BlockVector & operator=(const BlockVector< Number2 > &V)
BlockVector(const std::vector< size_type > &block_sizes)
BlockVector(const BlockVector< OtherNumber > &v)
BlockVector(const BlockVector< Number > &V)
void reinit(const size_type num_blocks, const size_type block_size=0, const bool omit_zeroing_entries=false)
void set_ghost_state(const bool ghosted) const
void add(const Number a, const BlockVector< Number > &V)
void add(const std::vector< size_type > &indices, const ::Vector< OtherNumber > &values)
BlockVector< Number > & operator+=(const BlockVector< Number > &V)
BlockVector(const std::vector< IndexSet > &local_ranges, const MPI_Comm communicator)
void scale(const BlockVector< Number > &scaling_factors)
void reinit(const std::vector< IndexSet > &local_ranges, const std::vector< IndexSet > &ghost_indices, const MPI_Comm communicator)
Number multivector_inner_product_with_metric(const FullMatrixType &matrix, const BlockVector< Number > &V, const bool symmetric=false) const
BlockVector< Number > & operator-=(const BlockVector< Number > &V)
BlockVector< Number > & operator=(const PETScWrappers::MPI::BlockVector &petsc_vec)
void reinit(const std::vector< std::shared_ptr< const Utilities::MPI::Partitioner >> &partitioners, const bool make_ghosted, const MPI_Comm &comm_sm=MPI_COMM_SELF)
typename BaseClass::const_reference const_reference
real_type lp_norm(const real_type p) const
void reinit(const BlockVector< Number2 > &V, const bool omit_zeroing_entries=false)
Number add_and_dot(const Number a, const BlockVector< Number > &V, const BlockVector< Number > &W)
void equ(const Number a, const BlockVector< Number > &V)
void sadd(const Number s, const BlockVector< Number > &V)
void copy_locally_owned_data_from(const BlockVector< Number2 > &src)
BlockVector(const std::vector< IndexSet > &local_ranges, const std::vector< IndexSet > &ghost_indices, const MPI_Comm communicator)
BlockVector< Number > & operator/=(const Number factor)
typename BaseClass::const_pointer const_pointer
typename BaseClass::const_iterator const_iterator
void multivector_inner_product(FullMatrixType &matrix, const BlockVector< Number > &V, const bool symmetric=false) const
void mmult(BlockVector< Number > &V, const FullMatrixType &matrix, const Number s=Number(0.), const Number b=Number(1.)) const
BlockVector< Number > & operator*=(const Number factor)
void reinit(const std::vector< size_type > &block_sizes, const bool omit_zeroing_entries=false)
void compress(VectorOperation::values operation)
void import_elements(const LinearAlgebra::ReadWriteVector< Number > &V, const VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
virtual size_type size() const override
#define DEAL_II_DEPRECATED
Definition: config.h:178
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcIteratorRangeDoesNotMatchVectorSize()
#define DeclException0(Exception0)
Definition: exceptions.h:472
static ::ExceptionBase & ExcVectorTypeNotCompatible()
constexpr int block_size
Definition: cuda_size.h:29
@ matrix
Contents is actually a matrix.
@ symmetric
Matrix is symmetric.
static const char V
BlockVector< double > BlockVector
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)