Reference documentation for deal.II version GIT 4e92747666 2022-12-02 08:40: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\}}\)
petsc_block_vector.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2021 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_petsc_block_vector_h
17 #define dealii_petsc_block_vector_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #ifdef DEAL_II_WITH_PETSC
23 
26 # include <deal.II/lac/exceptions.h>
29 
31 
32 
33 namespace PETScWrappers
34 {
35  // forward declaration
36  class BlockVector;
37 
38  namespace MPI
39  {
61  class BlockVector : public BlockVectorBase<Vector>
62  {
63  public:
68 
73 
85 
89  BlockVector() = default;
90 
97  explicit BlockVector(const unsigned int n_blocks,
98  const MPI_Comm & communicator,
99  const size_type block_size,
101 
106  BlockVector(const BlockVector &V);
107 
115  BlockVector(const std::vector<size_type> &block_sizes,
116  const MPI_Comm & communicator,
117  const std::vector<size_type> &local_elements);
118 
123  explicit BlockVector(const std::vector<IndexSet> &parallel_partitioning,
124  const MPI_Comm &communicator = MPI_COMM_WORLD);
125 
129  BlockVector(const std::vector<IndexSet> &parallel_partitioning,
130  const std::vector<IndexSet> &ghost_indices,
131  const MPI_Comm & communicator);
132 
133 
134 
138  ~BlockVector() override = default;
139 
144  BlockVector &
145  operator=(const value_type s);
146 
150  BlockVector &
151  operator=(const BlockVector &V);
152 
162  void
163  reinit(const unsigned int n_blocks,
164  const MPI_Comm & communicator,
165  const size_type block_size,
167  const bool omit_zeroing_entries = false);
168 
189  void
190  reinit(const std::vector<size_type> &block_sizes,
191  const MPI_Comm & communicator,
192  const std::vector<size_type> &locally_owned_sizes,
193  const bool omit_zeroing_entries = false);
194 
209  void
210  reinit(const BlockVector &V, const bool omit_zeroing_entries = false);
211 
216  void
217  reinit(const std::vector<IndexSet> &parallel_partitioning,
218  const MPI_Comm & communicator);
219 
223  void
224  reinit(const std::vector<IndexSet> &parallel_partitioning,
225  const std::vector<IndexSet> &ghost_entries,
226  const MPI_Comm & communicator);
227 
234  void
235  reinit(const unsigned int num_blocks);
236 
240  bool
241  has_ghost_elements() const;
242 
247  const MPI_Comm &
248  get_mpi_communicator() const;
249 
267  void
268  swap(BlockVector &v);
269 
273  void
274  print(std::ostream & out,
275  const unsigned int precision = 3,
276  const bool scientific = true,
277  const bool across = true) const;
278 
287  };
288 
291  /*--------------------- Inline functions --------------------------------*/
292 
293  inline BlockVector::BlockVector(const unsigned int n_blocks,
294  const MPI_Comm & communicator,
295  const size_type block_size,
296  const size_type locally_owned_size)
297  {
298  reinit(n_blocks, communicator, block_size, locally_owned_size);
299  }
300 
301 
302 
304  const std::vector<size_type> &block_sizes,
305  const MPI_Comm & communicator,
306  const std::vector<size_type> &local_elements)
307  {
308  reinit(block_sizes, communicator, local_elements, false);
309  }
310 
311 
314  {
315  this->block_indices = v.block_indices;
316 
317  this->components.resize(this->n_blocks());
318  for (unsigned int i = 0; i < this->n_blocks(); ++i)
319  this->components[i] = v.components[i];
320 
321  this->collect_sizes();
322  }
323 
325  const std::vector<IndexSet> &parallel_partitioning,
326  const MPI_Comm & communicator)
327  {
328  reinit(parallel_partitioning, communicator);
329  }
330 
332  const std::vector<IndexSet> &parallel_partitioning,
333  const std::vector<IndexSet> &ghost_indices,
334  const MPI_Comm & communicator)
335  {
336  reinit(parallel_partitioning, ghost_indices, communicator);
337  }
338 
339  inline BlockVector &
341  {
343  return *this;
344  }
345 
346  inline BlockVector &
348  {
349  // we only allow assignment to vectors with the same number of blocks
350  // or to an empty BlockVector
351  Assert(this->n_blocks() == 0 || this->n_blocks() == v.n_blocks(),
352  ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
353 
354  if (this->n_blocks() != v.n_blocks())
355  this->block_indices = v.block_indices;
356 
357  this->components.resize(this->n_blocks());
358  for (unsigned int i = 0; i < this->n_blocks(); ++i)
359  this->components[i] = v.components[i];
360 
361  this->collect_sizes();
362 
363  return *this;
364  }
365 
366 
367 
368  inline void
369  BlockVector::reinit(const unsigned int n_blocks,
370  const MPI_Comm & communicator,
371  const size_type block_size,
372  const size_type locally_owned_size,
373  const bool omit_zeroing_entries)
374  {
375  reinit(std::vector<size_type>(n_blocks, block_size),
376  communicator,
377  std::vector<size_type>(n_blocks, locally_owned_size),
378  omit_zeroing_entries);
379  }
380 
381 
382 
383  inline void
384  BlockVector::reinit(const std::vector<size_type> &block_sizes,
385  const MPI_Comm & communicator,
386  const std::vector<size_type> &locally_owned_sizes,
387  const bool omit_zeroing_entries)
388  {
389  this->block_indices.reinit(block_sizes);
390 
391  this->components.resize(this->n_blocks());
392  for (unsigned int i = 0; i < this->n_blocks(); ++i)
393  this->components[i].reinit(communicator,
394  block_sizes[i],
395  locally_owned_sizes[i],
396  omit_zeroing_entries);
397 
398  this->collect_sizes();
399  }
400 
401 
402  inline void
403  BlockVector::reinit(const BlockVector &v, const bool omit_zeroing_entries)
404  {
405  if (this->n_blocks() != v.n_blocks())
406  this->block_indices = v.get_block_indices();
407 
408  this->components.resize(this->n_blocks());
409  for (unsigned int i = 0; i < this->n_blocks(); ++i)
410  this->components[i].reinit(v.components[i], omit_zeroing_entries);
411 
412  this->collect_sizes();
413  }
414 
415  inline void
416  BlockVector::reinit(const std::vector<IndexSet> &parallel_partitioning,
417  const MPI_Comm & communicator)
418  {
419  // update the number of blocks
420  this->block_indices.reinit(parallel_partitioning.size(), 0);
421 
422  // initialize each block
423  this->components.resize(this->n_blocks());
424  for (unsigned int i = 0; i < this->n_blocks(); ++i)
425  this->components[i].reinit(parallel_partitioning[i], communicator);
426 
427  // update block_indices content
428  this->collect_sizes();
429  }
430 
431  inline void
432  BlockVector::reinit(const std::vector<IndexSet> &parallel_partitioning,
433  const std::vector<IndexSet> &ghost_entries,
434  const MPI_Comm & communicator)
435  {
436  AssertDimension(parallel_partitioning.size(), ghost_entries.size());
437 
438  // update the number of blocks
439  this->block_indices.reinit(parallel_partitioning.size(), 0);
440 
441  // initialize each block
442  this->components.resize(this->n_blocks());
443  for (unsigned int i = 0; i < this->n_blocks(); ++i)
444  this->components[i].reinit(parallel_partitioning[i],
445  ghost_entries[i],
446  communicator);
447 
448  // update block_indices content
449  this->collect_sizes();
450  }
451 
452 
453 
454  inline const MPI_Comm &
456  {
457  return block(0).get_mpi_communicator();
458  }
459 
460  inline bool
462  {
463  bool ghosted = block(0).has_ghost_elements();
464 # ifdef DEBUG
465  for (unsigned int i = 0; i < this->n_blocks(); ++i)
466  Assert(block(i).has_ghost_elements() == ghosted, ExcInternalError());
467 # endif
468  return ghosted;
469  }
470 
471 
472  inline void
474  {
475  std::swap(this->components, v.components);
476 
478  }
479 
480 
481 
482  inline void
483  BlockVector::print(std::ostream & out,
484  const unsigned int precision,
485  const bool scientific,
486  const bool across) const
487  {
488  for (unsigned int i = 0; i < this->n_blocks(); ++i)
489  {
490  if (across)
491  out << 'C' << i << ':';
492  else
493  out << "Component " << i << std::endl;
494  this->components[i].print(out, precision, scientific, across);
495  }
496  }
497 
498 
499 
507  inline void
509  {
510  u.swap(v);
511  }
512 
513  } // namespace MPI
514 
515 } // namespace PETScWrappers
516 
517 namespace internal
518 {
519  namespace LinearOperatorImplementation
520  {
521  template <typename>
522  class ReinitHelper;
523 
528  template <>
530  {
531  public:
532  template <typename Matrix>
533  static void
534  reinit_range_vector(const Matrix & matrix,
536  bool /*omit_zeroing_entries*/)
537  {
538  v.reinit(matrix.locally_owned_range_indices(),
539  matrix.get_mpi_communicator());
540  }
541 
542  template <typename Matrix>
543  static void
544  reinit_domain_vector(const Matrix & matrix,
546  bool /*omit_zeroing_entries*/)
547  {
548  v.reinit(matrix.locally_owned_domain_indices(),
549  matrix.get_mpi_communicator());
550  }
551  };
552 
553  } // namespace LinearOperatorImplementation
554 } /* namespace internal */
555 
556 
560 template <>
562 {};
563 
564 
566 
567 #endif // DEAL_II_WITH_PETSC
568 
569 #endif
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
const BlockIndices & get_block_indices() const
::internal::BlockVectorIterators::Iterator< BlockVectorBase, false > iterator
unsigned int n_blocks() const
const value_type * const_pointer
typename BlockType::const_reference const_reference
types::global_dof_index size_type
typename BlockType::value_type value_type
BlockVectorBase & operator=(const value_type s)
std::vector< Vector > components
::internal::BlockVectorIterators::Iterator< BlockVectorBase, true > const_iterator
typename BlockType::reference reference
BlockType & block(const unsigned int i)
std::size_t locally_owned_size() const
BaseClass::value_type value_type
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
void reinit(const unsigned int n_blocks, const MPI_Comm &communicator, const size_type block_size, const size_type locally_owned_size, const bool omit_zeroing_entries=false)
BaseClass::const_reference const_reference
void swap(BlockVector &u, BlockVector &v)
~BlockVector() override=default
BlockVector & operator=(const value_type s)
BaseClass::const_pointer const_pointer
const MPI_Comm & get_mpi_communicator() const
Definition: vector.h:109
bool has_ghost_elements() const
static void reinit_range_vector(const Matrix &matrix, PETScWrappers::MPI::BlockVector &v, bool)
static void reinit_domain_vector(const Matrix &matrix, PETScWrappers::MPI::BlockVector &v, bool)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:458
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:459
#define DeclException0(Exception0)
Definition: exceptions.h:464
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
Definition: exceptions.h:1501
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1695
static ::ExceptionBase & ExcNonMatchingBlockVectors()
static ::ExceptionBase & ExcIteratorRangeDoesNotMatchVectorSize()
constexpr int block_size
Definition: cuda_size.h:29
@ matrix
Contents is actually a matrix.
static const char V
BlockVector< double > BlockVector
std::enable_if_t< IsBlockVector< VectorType >::value, unsigned int > n_blocks(const VectorType &vector)
Definition: operators.h:50
void swap(BlockVector &u, BlockVector &v)