Reference documentation for deal.II version GIT 8e09676776 2023-03-27 21:15:01+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 - 2020 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 
218 
233  void
234  reinit(const BlockVector &V, const bool omit_zeroing_entries = false);
235 
242  void
243  reinit(const size_type num_blocks);
244 
262  void
264  const BlockVector & v);
265 
272  bool
273  has_ghost_elements() const;
274 
292  void
293  swap(BlockVector &v);
294 
298  void
299  print(std::ostream & out,
300  const unsigned int precision = 3,
301  const bool scientific = true,
302  const bool across = true) const;
303 
308 
313  };
314 
315 
316 
317  /*-------------------------- Inline functions ---------------------------*/
319  const std::vector<IndexSet> &parallel_partitioning,
320  const MPI_Comm & communicator)
321  {
322  reinit(parallel_partitioning, communicator, false);
323  }
324 
325 
326 
328  const std::vector<IndexSet> &parallel_partitioning,
329  const std::vector<IndexSet> &ghost_values,
330  const MPI_Comm & communicator,
331  const bool vector_writable)
332  {
333  reinit(parallel_partitioning,
334  ghost_values,
335  communicator,
336  vector_writable);
337  }
338 
339 
340 
341  inline BlockVector::BlockVector(const size_type num_blocks)
342  {
343  reinit(num_blocks);
344  }
345 
346 
347 
349  : ::BlockVectorBase<MPI::Vector>()
350  {
351  this->block_indices = v.block_indices;
352 
353  this->components.resize(this->n_blocks());
354  for (unsigned int i = 0; i < this->n_blocks(); ++i)
355  this->components[i] = v.components[i];
356 
357  this->collect_sizes();
358  }
359 
360 
361 
363  {
364  // initialize a minimal, valid object and swap
365  reinit(0);
366  swap(v);
367  }
368 
369 
370 
371  template <typename Number>
372  BlockVector &
373  BlockVector::operator=(const ::BlockVector<Number> &v)
374  {
375  // we only allow assignment to vectors with the same number of blocks
376  // or to an empty BlockVector
377  Assert(this->n_blocks() == 0 || this->n_blocks() == v.n_blocks(),
378  ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
379 
380  if (this->n_blocks() != v.n_blocks())
381  this->block_indices = v.get_block_indices();
382 
383  this->components.resize(this->n_blocks());
384  for (unsigned int i = 0; i < this->n_blocks(); ++i)
385  this->components[i] = v.block(i);
386 
387  this->collect_sizes();
388 
389  return *this;
390  }
391 
392 
393 
394  inline bool
396  {
397  bool ghosted = block(0).has_ghost_elements();
398 # ifdef DEBUG
399  for (unsigned int i = 0; i < this->n_blocks(); ++i)
400  Assert(block(i).has_ghost_elements() == ghosted, ExcInternalError());
401 # endif
402  return ghosted;
403  }
404 
405 
406 
407  inline void
409  {
410  std::swap(this->components, v.components);
411 
413  }
414 
415 
416 
424  inline void
426  {
427  u.swap(v);
428  }
429 
430  } /* namespace MPI */
431 
432 } /* namespace TrilinosWrappers */
433 
437 namespace internal
438 {
439  namespace LinearOperatorImplementation
440  {
441  template <typename>
442  class ReinitHelper;
443 
448  template <>
450  {
451  public:
452  template <typename Matrix>
453  static void
454  reinit_range_vector(const Matrix & matrix,
456  bool omit_zeroing_entries)
457  {
458  v.reinit(matrix.locally_owned_range_indices(),
459  matrix.get_mpi_communicator(),
460  omit_zeroing_entries);
461  }
462 
463  template <typename Matrix>
464  static void
465  reinit_domain_vector(const Matrix & matrix,
467  bool omit_zeroing_entries)
468  {
469  v.reinit(matrix.locally_owned_domain_indices(),
470  matrix.get_mpi_communicator(),
471  omit_zeroing_entries);
472  }
473  };
474 
475  } // namespace LinearOperatorImplementation
476 } /* namespace internal */
477 
478 
482 template <>
484 {};
485 
487 
488 #endif // DEAL_II_WITH_TRILINOS
489 
490 #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)
void reinit(const std::vector< IndexSet > &parallel_partitioning, const MPI_Comm &communicator=MPI_COMM_WORLD, const bool omit_zeroing_entries=false)
BlockVector & operator=(const value_type s)
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:474
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:475
#define DeclException0(Exception0)
Definition: exceptions.h:465
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
Definition: exceptions.h:1586
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