Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
trilinos_parallel_block_vector.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2008 - 2018 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 
24 # include <deal.II/lac/block_indices.h>
25 # include <deal.II/lac/block_vector_base.h>
26 # include <deal.II/lac/exceptions.h>
27 # include <deal.II/lac/trilinos_vector.h>
28 
29 # include <functional>
30 
31 DEAL_II_NAMESPACE_OPEN
32 
33 // forward declaration
34 template <typename Number>
35 class BlockVectorBase;
36 
41 namespace TrilinosWrappers
42 {
43  // forward declaration
44  namespace MPI
45  {
46  class BlockVector;
47  }
48  class BlockSparseMatrix;
49 
50 
51  namespace MPI
52  {
71  class BlockVector : public ::BlockVectorBase<MPI::Vector>
72  {
73  public:
78 
83 
87  using value_type = BaseClass::value_type;
88  using pointer = BaseClass::pointer;
89  using const_pointer = BaseClass::const_pointer;
90  using reference = BaseClass::reference;
91  using const_reference = BaseClass::const_reference;
92  using size_type = BaseClass::size_type;
95 
99  BlockVector() = default;
100 
107  explicit BlockVector(const std::vector<IndexSet> &parallel_partitioning,
108  const MPI_Comm &communicator = MPI_COMM_WORLD);
109 
115  BlockVector(const std::vector<IndexSet> &parallel_partitioning,
116  const std::vector<IndexSet> &ghost_values,
117  const MPI_Comm & communicator,
118  const bool vector_writable = false);
119 
124  BlockVector(const BlockVector &v);
125 
130  BlockVector(BlockVector &&v) noexcept;
131 
137  explicit BlockVector(const size_type num_blocks);
138 
142  ~BlockVector() override = default;
143 
148  BlockVector &
149  operator=(const value_type s);
150 
154  BlockVector &
155  operator=(const BlockVector &v);
156 
161  BlockVector &
162  operator=(BlockVector &&v) noexcept;
163 
174  template <typename Number>
175  BlockVector &
176  operator=(const ::BlockVector<Number> &v);
177 
186  void
187  reinit(const std::vector<IndexSet> &parallel_partitioning,
188  const MPI_Comm & communicator = MPI_COMM_WORLD,
189  const bool omit_zeroing_entries = false);
190 
208  void
209  reinit(const std::vector<IndexSet> &partitioning,
210  const std::vector<IndexSet> &ghost_values,
211  const MPI_Comm & communicator = MPI_COMM_WORLD,
212  const bool vector_writable = false);
213 
214 
229  void
230  reinit(const BlockVector &V, const bool omit_zeroing_entries = false);
231 
238  void
239  reinit(const size_type num_blocks);
240 
258  void
260  const BlockVector & v);
261 
268  bool
269  has_ghost_elements() const;
270 
288  void
289  swap(BlockVector &v);
290 
294  void
295  print(std::ostream & out,
296  const unsigned int precision = 3,
297  const bool scientific = true,
298  const bool across = true) const;
299 
304 
309  };
310 
311 
312 
313  /*-------------------------- Inline functions ---------------------------*/
315  const std::vector<IndexSet> &parallel_partitioning,
316  const MPI_Comm & communicator)
317  {
318  reinit(parallel_partitioning, communicator, false);
319  }
320 
321 
322 
324  const std::vector<IndexSet> &parallel_partitioning,
325  const std::vector<IndexSet> &ghost_values,
326  const MPI_Comm & communicator,
327  const bool vector_writable)
328  {
329  reinit(parallel_partitioning,
330  ghost_values,
331  communicator,
332  vector_writable);
333  }
334 
335 
336 
337  inline BlockVector::BlockVector(const size_type num_blocks)
338  {
339  reinit(num_blocks);
340  }
341 
342 
343 
345  : ::BlockVectorBase<MPI::Vector>()
346  {
347  this->components.resize(v.n_blocks());
348  this->block_indices = v.block_indices;
349 
350  for (size_type i = 0; i < this->n_blocks(); ++i)
351  this->components[i] = v.components[i];
352  }
353 
354 
355 
357  {
358  // initialize a minimal, valid object and swap
359  reinit(0);
360  swap(v);
361  }
362 
363 
364 
365  template <typename Number>
366  BlockVector &
367  BlockVector::operator=(const ::BlockVector<Number> &v)
368  {
369  if (n_blocks() != v.n_blocks())
370  {
371  std::vector<size_type> block_sizes(v.n_blocks(), 0);
372  block_indices.reinit(block_sizes);
373  if (components.size() != n_blocks())
374  components.resize(n_blocks());
375  }
376 
377  for (size_type i = 0; i < this->n_blocks(); ++i)
378  this->components[i] = v.block(i);
379 
380  collect_sizes();
381 
382  return *this;
383  }
384 
385 
386 
387  inline bool
389  {
390  bool ghosted = block(0).has_ghost_elements();
391 # ifdef DEBUG
392  for (unsigned int i = 0; i < this->n_blocks(); ++i)
393  Assert(block(i).has_ghost_elements() == ghosted, ExcInternalError());
394 # endif
395  return ghosted;
396  }
397 
398 
399 
400  inline void
402  {
403  std::swap(this->components, v.components);
404 
406  }
407 
408 
409 
418  inline void
420  {
421  u.swap(v);
422  }
423 
424  } /* namespace MPI */
425 
426 } /* namespace TrilinosWrappers */
427 
431 namespace internal
432 {
433  namespace LinearOperatorImplementation
434  {
435  template <typename>
436  class ReinitHelper;
437 
442  template <>
444  {
445  public:
446  template <typename Matrix>
447  static void
448  reinit_range_vector(const Matrix & matrix,
450  bool omit_zeroing_entries)
451  {
452  v.reinit(matrix.locally_owned_range_indices(),
453  matrix.get_mpi_communicator(),
454  omit_zeroing_entries);
455  }
456 
457  template <typename Matrix>
458  static void
459  reinit_domain_vector(const Matrix & matrix,
461  bool omit_zeroing_entries)
462  {
463  v.reinit(matrix.locally_owned_domain_indices(),
464  matrix.get_mpi_communicator(),
465  omit_zeroing_entries);
466  }
467  };
468 
469  } // namespace LinearOperatorImplementation
470 } /* namespace internal */
471 
472 
478 template <>
480 {};
481 
482 DEAL_II_NAMESPACE_CLOSE
483 
484 #endif // DEAL_II_WITH_TRILINOS
485 
486 #endif
static ::ExceptionBase & ExcIteratorRangeDoesNotMatchVectorSize()
static ::ExceptionBase & ExcNonMatchingBlockVectors()
void reinit(const std::vector< IndexSet > &parallel_partitioning, const MPI_Comm &communicator=MPI_COMM_WORLD, const bool omit_zeroing_entries=false)
static void reinit_range_vector(const Matrix &matrix, Vector &v, bool omit_zeroing_entries)
static void reinit_domain_vector(const Matrix &matrix, Vector &v, bool omit_zeroing_entries)
void swap(BlockVector &u, BlockVector &v)
#define Assert(cond, exc)
Definition: exceptions.h:1407
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
#define DeclException0(Exception0)
Definition: exceptions.h:473
std::vector< MPI::Vector > components
void swap(Vector< Number > &u, Vector< Number > &v)
Definition: vector.h:1376
BlockVector & operator=(const value_type s)
unsigned int n_blocks() const
void import_nonlocal_data_for_fe(const TrilinosWrappers::BlockSparseMatrix &m, const BlockVector &v)
LinearAlgebra::distributed::BlockVector< Number > BlockVector
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
static ::ExceptionBase & ExcInternalError()