Reference documentation for deal.II version GIT 9042b9283b 2023-12-02 14:50: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\}}\)
trilinos_parallel_block_vector.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2008 - 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_trilinos_parallel_block_vector_h
17 #define dealii_trilinos_parallel_block_vector_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #ifdef DEAL_II_WITH_TRILINOS
23 
26 # include <deal.II/lac/exceptions.h>
28 
29 # include <functional>
30 
32 
33 // forward declaration
34 # ifndef DOXYGEN
35 template <typename Number>
36 class BlockVectorBase;
37 
38 namespace TrilinosWrappers
39 {
40  // forward declaration
41  namespace MPI
42  {
43  class BlockVector;
44  }
45  class BlockSparseMatrix;
46 } // namespace TrilinosWrappers
47 # endif
48 
54 namespace TrilinosWrappers
55 {
56  namespace MPI
57  {
75  class BlockVector : public ::BlockVectorBase<MPI::Vector>
76  {
77  public:
82 
87 
99 
103  BlockVector() = default;
104 
111  explicit BlockVector(const std::vector<IndexSet> &parallel_partitioning,
112  const MPI_Comm communicator = MPI_COMM_WORLD);
113 
119  BlockVector(const std::vector<IndexSet> &parallel_partitioning,
120  const std::vector<IndexSet> &ghost_values,
121  const MPI_Comm communicator,
122  const bool vector_writable = false);
123 
128  BlockVector(const BlockVector &v);
129 
134  BlockVector(BlockVector &&v) noexcept;
135 
141  explicit BlockVector(const size_type num_blocks);
142 
146  ~BlockVector() override = default;
147 
152  BlockVector &
153  operator=(const value_type s);
154 
158  BlockVector &
159  operator=(const BlockVector &v);
160 
165  BlockVector &
166  operator=(BlockVector &&v) noexcept;
167 
178  template <typename Number>
179  BlockVector &
180  operator=(const ::BlockVector<Number> &v);
181 
190  void
191  reinit(const std::vector<IndexSet> &parallel_partitioning,
192  const MPI_Comm communicator = MPI_COMM_WORLD,
193  const bool omit_zeroing_entries = false);
194 
212  void
213  reinit(const std::vector<IndexSet> &partitioning,
214  const std::vector<IndexSet> &ghost_values,
215  const MPI_Comm communicator = MPI_COMM_WORLD,
216  const bool vector_writable = false);
217 
228  void
229  reinit(
230  const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
231  &partitioners,
232  const bool make_ghosted = true,
233  const bool vector_writable = false);
234 
249  void
250  reinit(const BlockVector &V, const bool omit_zeroing_entries = false);
251 
258  void
259  reinit(const size_type num_blocks);
260 
278  void
280  const BlockVector &v);
281 
288  bool
289  has_ghost_elements() const;
290 
308  void
309  swap(BlockVector &v);
310 
314  void
315  print(std::ostream &out,
316  const unsigned int precision = 3,
317  const bool scientific = true,
318  const bool across = true) const;
319 
324 
329  };
330 
331 
332 
333  /*-------------------------- Inline functions ---------------------------*/
335  const std::vector<IndexSet> &parallel_partitioning,
336  const MPI_Comm communicator)
337  {
338  reinit(parallel_partitioning, communicator, false);
339  }
340 
341 
342 
344  const std::vector<IndexSet> &parallel_partitioning,
345  const std::vector<IndexSet> &ghost_values,
346  const MPI_Comm communicator,
347  const bool vector_writable)
348  {
349  reinit(parallel_partitioning,
350  ghost_values,
351  communicator,
352  vector_writable);
353  }
354 
355 
356 
357  inline BlockVector::BlockVector(const size_type num_blocks)
358  {
359  reinit(num_blocks);
360  }
361 
362 
363 
365  : ::BlockVectorBase<MPI::Vector>()
366  {
367  this->block_indices = v.block_indices;
368 
369  this->components.resize(this->n_blocks());
370  for (unsigned int i = 0; i < this->n_blocks(); ++i)
371  this->components[i] = v.components[i];
372 
373  this->collect_sizes();
374  }
375 
376 
377 
379  {
380  // initialize a minimal, valid object and swap
381  reinit(0);
382  swap(v);
383  }
384 
385 
386 
387  template <typename Number>
388  BlockVector &
389  BlockVector::operator=(const ::BlockVector<Number> &v)
390  {
391  // we only allow assignment to vectors with the same number of blocks
392  // or to an empty BlockVector
393  Assert(this->n_blocks() == 0 || this->n_blocks() == v.n_blocks(),
394  ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
395 
396  if (this->n_blocks() != v.n_blocks())
397  this->block_indices = v.get_block_indices();
398 
399  this->components.resize(this->n_blocks());
400  for (unsigned int i = 0; i < this->n_blocks(); ++i)
401  this->components[i] = v.block(i);
402 
403  this->collect_sizes();
404 
405  return *this;
406  }
407 
408 
409 
410  inline bool
412  {
413  bool ghosted = block(0).has_ghost_elements();
414 # ifdef DEBUG
415  for (unsigned int i = 0; i < this->n_blocks(); ++i)
416  Assert(block(i).has_ghost_elements() == ghosted, ExcInternalError());
417 # endif
418  return ghosted;
419  }
420 
421 
422 
423  inline void
425  {
426  std::swap(this->components, v.components);
427 
429  }
430 
431 
432 
440  inline void
442  {
443  u.swap(v);
444  }
445 
446  } /* namespace MPI */
447 
448 } /* namespace TrilinosWrappers */
449 
453 namespace internal
454 {
455  namespace LinearOperatorImplementation
456  {
457  template <typename>
458  class ReinitHelper;
459 
464  template <>
466  {
467  public:
468  template <typename Matrix>
469  static void
472  bool omit_zeroing_entries)
473  {
474  v.reinit(matrix.locally_owned_range_indices(),
475  matrix.get_mpi_communicator(),
476  omit_zeroing_entries);
477  }
478 
479  template <typename Matrix>
480  static void
483  bool omit_zeroing_entries)
484  {
485  v.reinit(matrix.locally_owned_domain_indices(),
486  matrix.get_mpi_communicator(),
487  omit_zeroing_entries);
488  }
489  };
490 
491  } // namespace LinearOperatorImplementation
492 } /* namespace internal */
493 
494 
498 template <>
500 {};
501 
503 
504 #endif // DEAL_II_WITH_TRILINOS
505 
506 #endif
::internal::BlockVectorIterators::Iterator< BlockVectorBase, false > iterator
unsigned int n_blocks() const
typename BlockType::const_reference const_reference
types::global_dof_index size_type
typename BlockType::value_type value_type
std::vector< MPI::Vector > components
::internal::BlockVectorIterators::Iterator< BlockVectorBase, true > const_iterator
typename BlockType::reference reference
BlockType & block(const unsigned int i)
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
void swap(BlockVector &u, BlockVector &v)
BlockVector & operator=(const value_type s)
void reinit(const std::vector< IndexSet > &parallel_partitioning, const MPI_Comm communicator=MPI_COMM_WORLD, const bool omit_zeroing_entries=false)
void import_nonlocal_data_for_fe(const TrilinosWrappers::BlockSparseMatrix &m, const BlockVector &v)
static void reinit_range_vector(const Matrix &matrix, TrilinosWrappers::MPI::BlockVector &v, bool omit_zeroing_entries)
static void reinit_domain_vector(const Matrix &matrix, TrilinosWrappers::MPI::BlockVector &v, bool omit_zeroing_entries)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
#define DeclException0(Exception0)
Definition: exceptions.h:472
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
Definition: exceptions.h:1631
static ::ExceptionBase & ExcIteratorRangeDoesNotMatchVectorSize()
static ::ExceptionBase & ExcNonMatchingBlockVectors()
@ matrix
Contents is actually a matrix.
static const char V
BlockVector< double > BlockVector
void swap(BlockVector &u, BlockVector &v)
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618