Reference documentation for deal.II version Git 4abc4a1666 2020-07-04 19:58:34 +0200
\(\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 - 2019 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 
53 namespace TrilinosWrappers
54 {
55  namespace MPI
56  {
74  class BlockVector : public ::BlockVectorBase<MPI::Vector>
75  {
76  public:
81 
86 
90  using value_type = BaseClass::value_type;
98 
102  BlockVector() = default;
103 
110  explicit BlockVector(const std::vector<IndexSet> &parallel_partitioning,
111  const MPI_Comm &communicator = MPI_COMM_WORLD);
112 
118  BlockVector(const std::vector<IndexSet> &parallel_partitioning,
119  const std::vector<IndexSet> &ghost_values,
120  const MPI_Comm & communicator,
121  const bool vector_writable = false);
122 
127  BlockVector(const BlockVector &v);
128 
133  BlockVector(BlockVector &&v) noexcept;
134 
140  explicit BlockVector(const size_type num_blocks);
141 
145  ~BlockVector() override = default;
146 
151  BlockVector &
152  operator=(const value_type s);
153 
157  BlockVector &
158  operator=(const BlockVector &v);
159 
164  BlockVector &
165  operator=(BlockVector &&v) noexcept;
166 
177  template <typename Number>
178  BlockVector &
179  operator=(const ::BlockVector<Number> &v);
180 
189  void
190  reinit(const std::vector<IndexSet> &parallel_partitioning,
191  const MPI_Comm & communicator = MPI_COMM_WORLD,
192  const bool omit_zeroing_entries = false);
193 
211  void
212  reinit(const std::vector<IndexSet> &partitioning,
213  const std::vector<IndexSet> &ghost_values,
214  const MPI_Comm & communicator = MPI_COMM_WORLD,
215  const bool vector_writable = false);
216 
217 
232  void
233  reinit(const BlockVector &V, const bool omit_zeroing_entries = false);
234 
241  void
242  reinit(const size_type num_blocks);
243 
261  void
262  import_nonlocal_data_for_fe(const TrilinosWrappers::BlockSparseMatrix &m,
263  const BlockVector & v);
264 
271  bool
272  has_ghost_elements() const;
273 
291  void
292  swap(BlockVector &v);
293 
297  void
298  print(std::ostream & out,
299  const unsigned int precision = 3,
300  const bool scientific = true,
301  const bool across = true) const;
302 
306  DeclException0(ExcIteratorRangeDoesNotMatchVectorSize);
307 
311  DeclException0(ExcNonMatchingBlockVectors);
312  };
313 
314 
315 
316  /*-------------------------- Inline functions ---------------------------*/
318  const std::vector<IndexSet> &parallel_partitioning,
319  const MPI_Comm & communicator)
320  {
321  reinit(parallel_partitioning, communicator, false);
322  }
323 
324 
325 
327  const std::vector<IndexSet> &parallel_partitioning,
328  const std::vector<IndexSet> &ghost_values,
329  const MPI_Comm & communicator,
330  const bool vector_writable)
331  {
332  reinit(parallel_partitioning,
333  ghost_values,
334  communicator,
335  vector_writable);
336  }
337 
338 
339 
340  inline BlockVector::BlockVector(const size_type num_blocks)
341  {
342  reinit(num_blocks);
343  }
344 
345 
346 
348  : ::BlockVectorBase<MPI::Vector>()
349  {
350  this->components.resize(v.n_blocks());
351  this->block_indices = v.block_indices;
352 
353  for (size_type i = 0; i < this->n_blocks(); ++i)
354  this->components[i] = v.components[i];
355  }
356 
357 
358 
360  {
361  // initialize a minimal, valid object and swap
362  reinit(0);
363  swap(v);
364  }
365 
366 
367 
368  template <typename Number>
369  BlockVector &
370  BlockVector::operator=(const ::BlockVector<Number> &v)
371  {
372  if (n_blocks() != v.n_blocks())
373  {
374  std::vector<size_type> block_sizes(v.n_blocks(), 0);
375  block_indices.reinit(block_sizes);
376  if (components.size() != n_blocks())
377  components.resize(n_blocks());
378  }
379 
380  for (size_type i = 0; i < this->n_blocks(); ++i)
381  this->components[i] = v.block(i);
382 
383  collect_sizes();
384 
385  return *this;
386  }
387 
388 
389 
390  inline bool
392  {
393  bool ghosted = block(0).has_ghost_elements();
394 # ifdef DEBUG
395  for (unsigned int i = 0; i < this->n_blocks(); ++i)
396  Assert(block(i).has_ghost_elements() == ghosted, ExcInternalError());
397 # endif
398  return ghosted;
399  }
400 
401 
402 
403  inline void
405  {
406  std::swap(this->components, v.components);
407 
409  }
410 
411 
412 
420  inline void
422  {
423  u.swap(v);
424  }
425 
426  } /* namespace MPI */
427 
428 } /* namespace TrilinosWrappers */
429 
433 namespace internal
434 {
435  namespace LinearOperatorImplementation
436  {
437  template <typename>
438  class ReinitHelper;
439 
444  template <>
446  {
447  public:
448  template <typename Matrix>
449  static void
450  reinit_range_vector(const Matrix & matrix,
452  bool omit_zeroing_entries)
453  {
454  v.reinit(matrix.locally_owned_range_indices(),
455  matrix.get_mpi_communicator(),
456  omit_zeroing_entries);
457  }
458 
459  template <typename Matrix>
460  static void
461  reinit_domain_vector(const Matrix & matrix,
463  bool omit_zeroing_entries)
464  {
465  v.reinit(matrix.locally_owned_domain_indices(),
466  matrix.get_mpi_communicator(),
467  omit_zeroing_entries);
468  }
469  };
470 
471  } // namespace LinearOperatorImplementation
472 } /* namespace internal */
473 
474 
478 template <>
480 {};
481 
483 
484 #endif // DEAL_II_WITH_TRILINOS
485 
486 #endif
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
static void reinit_range_vector(const Matrix &matrix, TrilinosWrappers::MPI::BlockVector &v, bool omit_zeroing_entries)
Contents is actually a matrix.
types::global_dof_index size_type
Definition: cuda_kernels.h:45
static const char V
void reinit(const std::vector< IndexSet > &parallel_partitioning, const MPI_Comm &communicator=MPI_COMM_WORLD, const bool omit_zeroing_entries=false)
void swap(BlockVector &u, BlockVector &v)
BlockVector< double > BlockVector
#define Assert(cond, exc)
Definition: exceptions.h:1403
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
std::vector< MPI::Vector > components
BlockVector & operator=(const value_type s)
unsigned int n_blocks() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
static void reinit_domain_vector(const Matrix &matrix, TrilinosWrappers::MPI::BlockVector &v, bool omit_zeroing_entries)
typename BlockType::const_reference const_reference
BlockType & block(const unsigned int i)
typename BlockType::reference reference
static ::ExceptionBase & ExcInternalError()