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
petsc_block_vector.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2004 - 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_petsc_block_vector_h
17#define dealii_petsc_block_vector_h
18
19
20#include <deal.II/base/config.h>
21
22#ifdef DEAL_II_WITH_PETSC
23
29
31
32
33namespace PETScWrappers
34{
35 // forward declaration
36 class BlockVector;
37
38 namespace MPI
39 {
61 class BlockVector : public BlockVectorBase<Vector>
62 {
63 public:
68
73
85
90
97 explicit BlockVector(const unsigned int n_blocks,
98 const MPI_Comm communicator,
99 const size_type block_size,
101
106 BlockVector(const BlockVector &V);
107
115 BlockVector(const std::vector<size_type> &block_sizes,
116 const MPI_Comm communicator,
117 const std::vector<size_type> &local_elements);
118
123 explicit BlockVector(const std::vector<IndexSet> &parallel_partitioning,
124 const MPI_Comm communicator = MPI_COMM_WORLD);
125
129 BlockVector(const std::vector<IndexSet> &parallel_partitioning,
130 const std::vector<IndexSet> &ghost_indices,
131 const MPI_Comm communicator);
132
139 explicit BlockVector(Vec v);
140
144 template <size_t num_blocks>
145 explicit BlockVector(const std::array<Vec, num_blocks> &);
146
150 ~BlockVector() override;
151
157 operator=(const value_type s);
158
163 operator=(const BlockVector &V);
164
170 void
171 reinit(Vec v);
172
182 void
183 reinit(const unsigned int n_blocks,
184 const MPI_Comm communicator,
185 const size_type block_size,
187 const bool omit_zeroing_entries = false);
188
209 void
210 reinit(const std::vector<size_type> &block_sizes,
211 const MPI_Comm communicator,
212 const std::vector<size_type> &locally_owned_sizes,
213 const bool omit_zeroing_entries = false);
214
229 void
230 reinit(const BlockVector &V, const bool omit_zeroing_entries = false);
231
236 void
237 reinit(const std::vector<IndexSet> &parallel_partitioning,
238 const MPI_Comm communicator);
239
243 void
244 reinit(const std::vector<IndexSet> &parallel_partitioning,
245 const std::vector<IndexSet> &ghost_entries,
246 const MPI_Comm communicator);
247
255 void
256 reinit(
257 const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
258 & partitioners,
259 const bool make_ghosted = true);
260
268 void
270
279 void
281
288 void
289 reinit(const unsigned int num_blocks);
290
294 bool
295 has_ghost_elements() const;
296
301 get_mpi_communicator() const;
302
310 operator const Vec &() const;
311
317 Vec &
318 petsc_vector();
319
337 void
338 swap(BlockVector &v);
339
343 void
344 print(std::ostream & out,
345 const unsigned int precision = 3,
346 const bool scientific = true,
347 const bool across = true) const;
348
357
358 private:
366
370 void
372 };
373
376 /*--------------------- Inline functions --------------------------------*/
377
380 , petsc_nest_vector(nullptr)
381 {}
382
383
384
385 inline BlockVector::BlockVector(const unsigned int n_blocks,
386 const MPI_Comm communicator,
387 const size_type block_size,
388 const size_type locally_owned_size)
389 : BlockVector()
390 {
391 reinit(n_blocks, communicator, block_size, locally_owned_size);
392 }
393
394
395
397 const std::vector<size_type> &block_sizes,
398 const MPI_Comm communicator,
399 const std::vector<size_type> &local_elements)
400 : BlockVector()
401 {
402 reinit(block_sizes, communicator, local_elements, false);
403 }
404
405
406
408 : BlockVector()
409 {
411
412 this->components.resize(this->n_blocks());
413 for (unsigned int i = 0; i < this->n_blocks(); ++i)
414 this->components[i] = v.components[i];
415
416 this->collect_sizes();
417 }
418
419
420
422 const std::vector<IndexSet> &parallel_partitioning,
423 const MPI_Comm communicator)
424 : BlockVector()
425 {
426 reinit(parallel_partitioning, communicator);
427 }
428
429
430
432 const std::vector<IndexSet> &parallel_partitioning,
433 const std::vector<IndexSet> &ghost_indices,
434 const MPI_Comm communicator)
435 : BlockVector()
436 {
437 reinit(parallel_partitioning, ghost_indices, communicator);
438 }
439
440
441
443 : BlockVector()
444 {
445 this->reinit(v);
446 }
447
448
449
450 template <size_t num_blocks>
451 inline BlockVector::BlockVector(const std::array<Vec, num_blocks> &arrayV)
452 : BlockVector()
453 {
454 this->block_indices.reinit(num_blocks, 0);
455
456 this->components.resize(num_blocks);
457 for (auto i = 0; i < num_blocks; ++i)
458 this->components[i].reinit(arrayV[i]);
459 this->collect_sizes();
460 }
461
462
463
464 inline BlockVector &
466 {
468 return *this;
469 }
470
471
472
473 inline BlockVector &
475 {
476 // we only allow assignment to vectors with the same number of blocks
477 // or to an empty BlockVector
478 Assert(this->n_blocks() == 0 || this->n_blocks() == v.n_blocks(),
480
481 if (this->n_blocks() != v.n_blocks())
483
484 this->components.resize(this->n_blocks());
485 for (unsigned int i = 0; i < this->n_blocks(); ++i)
486 this->components[i] = v.components[i];
487
488 this->collect_sizes();
489
490 return *this;
491 }
492
493
494
495 inline void
496 BlockVector::reinit(const unsigned int n_blocks,
497 const MPI_Comm communicator,
498 const size_type block_size,
499 const size_type locally_owned_size,
500 const bool omit_zeroing_entries)
501 {
502 reinit(std::vector<size_type>(n_blocks, block_size),
503 communicator,
504 std::vector<size_type>(n_blocks, locally_owned_size),
505 omit_zeroing_entries);
506 }
507
508
509
510 inline void
511 BlockVector::reinit(const std::vector<size_type> &block_sizes,
512 const MPI_Comm communicator,
513 const std::vector<size_type> &locally_owned_sizes,
514 const bool omit_zeroing_entries)
515 {
516 this->block_indices.reinit(block_sizes);
517
518 this->components.resize(this->n_blocks());
519 for (unsigned int i = 0; i < this->n_blocks(); ++i)
520 this->components[i].reinit(communicator,
521 block_sizes[i],
522 locally_owned_sizes[i],
523 omit_zeroing_entries);
524
525 this->collect_sizes();
526 }
527
528
529 inline void
530 BlockVector::reinit(const BlockVector &v, const bool omit_zeroing_entries)
531 {
532 if (this->n_blocks() != v.n_blocks())
534
535 this->components.resize(this->n_blocks());
536 for (unsigned int i = 0; i < this->n_blocks(); ++i)
537 this->components[i].reinit(v.components[i], omit_zeroing_entries);
538
539 this->collect_sizes();
540 }
541
542
543
544 inline void
545 BlockVector::reinit(const std::vector<IndexSet> &parallel_partitioning,
546 const MPI_Comm communicator)
547 {
548 // update the number of blocks
549 this->block_indices.reinit(parallel_partitioning.size(), 0);
550
551 // initialize each block
552 this->components.resize(this->n_blocks());
553 for (unsigned int i = 0; i < this->n_blocks(); ++i)
554 this->components[i].reinit(parallel_partitioning[i], communicator);
555
556 // update block_indices content
557 this->collect_sizes();
558 }
559
560
561
562 inline void
563 BlockVector::reinit(const std::vector<IndexSet> &parallel_partitioning,
564 const std::vector<IndexSet> &ghost_entries,
565 const MPI_Comm communicator)
566 {
567 AssertDimension(parallel_partitioning.size(), ghost_entries.size());
568
569 // update the number of blocks
570 this->block_indices.reinit(parallel_partitioning.size(), 0);
571
572 // initialize each block
573 this->components.resize(this->n_blocks());
574 for (unsigned int i = 0; i < this->n_blocks(); ++i)
575 this->components[i].reinit(parallel_partitioning[i],
576 ghost_entries[i],
577 communicator);
578
579 // update block_indices content
580 this->collect_sizes();
581 }
582
583
584
585 inline void
587 const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
588 & partitioners,
589 const bool make_ghosted)
590 {
591 // update the number of blocks
592 this->block_indices.reinit(partitioners.size(), 0);
593
594 // initialize each block
595 this->components.resize(this->n_blocks());
596 for (unsigned int i = 0; i < this->n_blocks(); ++i)
597 this->components[i].reinit(partitioners[i], make_ghosted);
598
599 // update block_indices content
600 this->collect_sizes();
601 }
602
603
604
605 inline MPI_Comm
607 {
608 return PetscObjectComm(reinterpret_cast<PetscObject>(petsc_nest_vector));
609 }
610
611
612
613 inline bool
615 {
616 bool ghosted = block(0).has_ghost_elements();
617# ifdef DEBUG
618 for (unsigned int i = 0; i < this->n_blocks(); ++i)
620# endif
621 return ghosted;
622 }
623
624
625
626 inline void
628 {
629 std::swap(this->components, v.components);
630 std::swap(this->petsc_nest_vector, v.petsc_nest_vector);
631
633 }
634
635
636
637 inline void
638 BlockVector::print(std::ostream & out,
639 const unsigned int precision,
640 const bool scientific,
641 const bool across) const
642 {
643 for (unsigned int i = 0; i < this->n_blocks(); ++i)
644 {
645 if (across)
646 out << 'C' << i << ':';
647 else
648 out << "Component " << i << std::endl;
649 this->components[i].print(out, precision, scientific, across);
650 }
651 }
652
653
654
662 inline void
664 {
665 u.swap(v);
666 }
667 } // namespace MPI
668
669} // namespace PETScWrappers
670
671namespace internal
672{
673 namespace LinearOperatorImplementation
674 {
675 template <typename>
676 class ReinitHelper;
677
682 template <>
683 class ReinitHelper<PETScWrappers::MPI::BlockVector>
684 {
685 public:
686 template <typename Matrix>
687 static void
688 reinit_range_vector(const Matrix & matrix,
690 bool /*omit_zeroing_entries*/)
691 {
692 v.reinit(matrix.locally_owned_range_indices(),
693 matrix.get_mpi_communicator());
694 }
695
696 template <typename Matrix>
697 static void
698 reinit_domain_vector(const Matrix & matrix,
700 bool /*omit_zeroing_entries*/)
701 {
702 v.reinit(matrix.locally_owned_domain_indices(),
703 matrix.get_mpi_communicator());
704 }
705 };
706
707 } // namespace LinearOperatorImplementation
708} /* namespace internal */
709
710
714template <>
715struct is_serial_vector<PETScWrappers::MPI::BlockVector> : std::false_type
716{};
717
718
720
721#endif // DEAL_II_WITH_PETSC
722
723#endif
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
::internal::BlockVectorIterators::Iterator< BlockVectorBase, false > iterator
unsigned int n_blocks() const
const value_type * const_pointer
typename BlockType::const_reference const_reference
types::global_dof_index size_type
typename BlockType::value_type value_type
BlockVectorBase & operator=(const value_type s)
std::vector< Vector > components
::internal::BlockVectorIterators::Iterator< BlockVectorBase, true > const_iterator
typename BlockType::reference reference
BlockType & block(const unsigned int i)
std::size_t locally_owned_size() const
const BlockIndices & get_block_indices() const
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
BaseClass::const_reference const_reference
void swap(BlockVector &u, BlockVector &v)
BlockVector & operator=(const value_type s)
void compress(VectorOperation::values operation)
BaseClass::const_pointer const_pointer
bool has_ghost_elements() const
static void reinit_range_vector(const Matrix &matrix, PETScWrappers::MPI::BlockVector &v, bool)
static void reinit_domain_vector(const Matrix &matrix, PETScWrappers::MPI::BlockVector &v, bool)
#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
static ::ExceptionBase & ExcNonMatchingBlockVectors()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcIteratorRangeDoesNotMatchVectorSize()