deal.II version GIT relicensing-2890-gf173123fa3 2025-03-22 01:40:00+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\}}\)
Loading...
Searching...
No Matches
trilinos_parallel_block_vector.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2012 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_trilinos_parallel_block_vector_h
16#define dealii_trilinos_parallel_block_vector_h
17
18
19#include <deal.II/base/config.h>
20
21#ifdef DEAL_II_WITH_TRILINOS
22
27
28# include <functional>
29
31
32// forward declaration
33# ifndef DOXYGEN
34template <typename Number>
35class BlockVectorBase;
36
37namespace TrilinosWrappers
38{
39 // forward declaration
40 namespace MPI
41 {
42 class BlockVector;
43 }
45} // namespace TrilinosWrappers
46# endif
47
53namespace TrilinosWrappers
54{
55 namespace MPI
56 {
74 class BlockVector : public ::BlockVectorBase<MPI::Vector>
75 {
76 public:
81
86
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
152 operator=(const value_type s);
153
165 operator=(const BlockVector &v);
166
172 operator=(BlockVector &&v) noexcept;
173
184 template <typename Number>
186 operator=(const ::BlockVector<Number> &v);
187
196 void
197 reinit(const std::vector<IndexSet> &parallel_partitioning,
198 const MPI_Comm communicator = MPI_COMM_WORLD,
199 const bool omit_zeroing_entries = false);
200
218 void
219 reinit(const std::vector<IndexSet> &partitioning,
220 const std::vector<IndexSet> &ghost_values,
221 const MPI_Comm communicator = MPI_COMM_WORLD,
222 const bool vector_writable = false);
223
234 void
235 reinit(
236 const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
237 &partitioners,
238 const bool make_ghosted = true,
239 const bool vector_writable = false);
240
255 void
256 reinit(const BlockVector &V, const bool omit_zeroing_entries = false);
257
264 void
265 reinit(const size_type num_blocks);
266
284 void
286 const BlockVector &v);
287
294 bool
295 has_ghost_elements() const;
296
314 void
315 swap(BlockVector &v) noexcept;
316
320 void
321 print(std::ostream &out,
322 const unsigned int precision = 3,
323 const bool scientific = true,
324 const bool across = true) const;
325
330
335 };
336
337
338
339 /*-------------------------- Inline functions ---------------------------*/
341 const std::vector<IndexSet> &parallel_partitioning,
342 const MPI_Comm communicator)
343 {
344 reinit(parallel_partitioning, communicator, false);
345 }
346
347
348
350 const std::vector<IndexSet> &parallel_partitioning,
351 const std::vector<IndexSet> &ghost_values,
352 const MPI_Comm communicator,
353 const bool vector_writable)
354 {
355 reinit(parallel_partitioning,
356 ghost_values,
357 communicator,
358 vector_writable);
359 }
360
361
362
363 inline BlockVector::BlockVector(const size_type num_blocks)
364 {
365 reinit(num_blocks);
366 }
367
368
369
372 {
374
375 this->components.resize(this->n_blocks());
376 for (unsigned int i = 0; i < this->n_blocks(); ++i)
377 this->components[i] = v.components[i];
378
379 this->collect_sizes();
380 }
381
382
383
385 {
386 // initialize a minimal, valid object and swap
387 reinit(0);
388 swap(v);
389 }
390
391
392
393 template <typename Number>
395 BlockVector::operator=(const ::BlockVector<Number> &v)
396 {
397 // we only allow assignment to vectors with the same number of blocks
398 // or to an empty BlockVector
399 Assert(this->n_blocks() == 0 || this->n_blocks() == v.n_blocks(),
400 ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
401
402 if (this->n_blocks() != v.n_blocks())
403 this->block_indices = v.get_block_indices();
404
405 this->components.resize(this->n_blocks());
406 for (unsigned int i = 0; i < this->n_blocks(); ++i)
407 this->components[i] = v.block(i);
408
409 this->collect_sizes();
410
411 return *this;
412 }
413
414
415
416 inline bool
418 {
419 bool ghosted = block(0).has_ghost_elements();
420 if constexpr (running_in_debug_mode())
421 {
422 for (unsigned int i = 0; i < this->n_blocks(); ++i)
423 Assert(block(i).has_ghost_elements() == ghosted,
425 }
426 return ghosted;
427 }
428
429
430
431 inline void
433 {
434 std::swap(this->components, v.components);
435
436 ::swap(this->block_indices, v.block_indices);
437 }
438
439
440
448 inline void
449 swap(BlockVector &u, BlockVector &v) noexcept
450 {
451 u.swap(v);
452 }
453
454 } /* namespace MPI */
455
456} /* namespace TrilinosWrappers */
457
461namespace internal
462{
463 namespace LinearOperatorImplementation
464 {
465 template <typename>
466 class ReinitHelper;
467
472 template <>
473 class ReinitHelper<TrilinosWrappers::MPI::BlockVector>
474 {
475 public:
476 template <typename Matrix>
477 static void
478 reinit_range_vector(const Matrix &matrix,
480 bool omit_zeroing_entries)
481 {
482 v.reinit(matrix.locally_owned_range_indices(),
483 matrix.get_mpi_communicator(),
484 omit_zeroing_entries);
485 }
486
487 template <typename Matrix>
488 static void
489 reinit_domain_vector(const Matrix &matrix,
491 bool omit_zeroing_entries)
492 {
493 v.reinit(matrix.locally_owned_domain_indices(),
494 matrix.get_mpi_communicator(),
495 omit_zeroing_entries);
496 }
497 };
498
499 } // namespace LinearOperatorImplementation
500} /* namespace internal */
501
502
506template <>
507struct is_serial_vector<TrilinosWrappers::MPI::BlockVector> : std::false_type
508{};
509
511
512#endif // DEAL_II_WITH_TRILINOS
513
514#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) noexcept
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:40
constexpr bool running_in_debug_mode()
Definition config.h:78
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:41
#define DeclException0(Exception0)
#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) noexcept