Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
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
28
29# include <functional>
30
32
33// forward declaration
34# ifndef DOXYGEN
35template <typename Number>
36class BlockVectorBase;
37
38namespace TrilinosWrappers
39{
40 // forward declaration
41 namespace MPI
42 {
43 class BlockVector;
44 }
46} // namespace TrilinosWrappers
47# endif
48
54namespace 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
153 operator=(const value_type s);
154
159 operator=(const BlockVector &v);
160
166 operator=(BlockVector &&v) noexcept;
167
178 template <typename Number>
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 {
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>
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
453namespace internal
454{
455 namespace LinearOperatorImplementation
456 {
457 template <typename>
458 class ReinitHelper;
459
464 template <>
465 class ReinitHelper<TrilinosWrappers::MPI::BlockVector>
466 {
467 public:
468 template <typename Matrix>
469 static void
470 reinit_range_vector(const Matrix & matrix,
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
481 reinit_domain_vector(const Matrix & matrix,
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
498template <>
499struct is_serial_vector<TrilinosWrappers::MPI::BlockVector> : std::false_type
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:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define DeclException0(Exception0)
Definition exceptions.h:465
#define Assert(cond, exc)
static ::ExceptionBase & ExcNonMatchingBlockVectors()
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcIteratorRangeDoesNotMatchVectorSize()
void swap(BlockVector &u, BlockVector &v)