17 #ifndef __deal2__petsc_parallel_block_vector_h
18 #define __deal2__petsc_parallel_block_vector_h
21 #include <deal.II/base/config.h>
23 #ifdef DEAL_II_WITH_PETSC
25 # include <deal.II/lac/petsc_parallel_vector.h>
26 # include <deal.II/lac/block_indices.h>
27 # include <deal.II/lac/block_vector_base.h>
28 # include <deal.II/lac/exceptions.h>
30 DEAL_II_NAMESPACE_OPEN
33 namespace PETScWrappers
82 typedef BaseClass::pointer pointer;
83 typedef BaseClass::const_pointer const_pointer;
84 typedef BaseClass::reference reference;
85 typedef BaseClass::const_reference const_reference;
86 typedef BaseClass::size_type size_type;
106 const MPI_Comm &communicator,
107 const size_type block_size,
108 const size_type local_size);
130 BlockVector (
const std::vector<size_type> &block_sizes,
131 const MPI_Comm &communicator,
132 const std::vector<size_type> &local_elements);
138 explicit BlockVector (
const std::vector<IndexSet> ¶llel_partitioning,
139 const MPI_Comm &communicator = MPI_COMM_WORLD);
144 BlockVector (
const std::vector<IndexSet> ¶llel_partitioning,
145 const std::vector<IndexSet> &ghost_indices,
146 const MPI_Comm &communicator);
215 const MPI_Comm &communicator,
216 const size_type block_size,
217 const size_type local_size,
218 const bool fast =
false);
255 void reinit (
const std::vector<size_type> &block_sizes,
256 const MPI_Comm &communicator,
257 const std::vector<size_type> &local_sizes,
258 const bool fast=
false);
287 const bool fast=
false);
293 void reinit (
const std::vector<IndexSet> ¶llel_partitioning,
294 const MPI_Comm &communicator);
299 void reinit (
const std::vector<IndexSet> ¶llel_partitioning,
300 const std::vector<IndexSet> &ghost_entries,
301 const MPI_Comm &communicator);
314 void reinit (
const unsigned int num_blocks);
363 void print (std::ostream &out,
364 const unsigned int precision = 3,
365 const bool scientific =
true,
366 const bool across =
true)
const;
391 const MPI_Comm &communicator,
392 const size_type block_size,
393 const size_type local_size)
395 reinit (n_blocks, communicator, block_size, local_size);
402 const MPI_Comm &communicator,
403 const std::vector<size_type> &local_elements)
405 reinit (block_sizes, communicator, local_elements,
false);
417 for (
unsigned int i=0; i<this->
n_blocks(); ++i)
423 const MPI_Comm &communicator)
425 reinit(parallel_partitioning, communicator);
430 const std::vector<IndexSet> &ghost_indices,
431 const MPI_Comm &communicator)
433 reinit(parallel_partitioning, ghost_indices, communicator);
456 for (size_type i=0; i<this->
n_blocks(); ++i)
472 const MPI_Comm &communicator,
473 const size_type block_size,
474 const size_type local_size,
477 reinit(std::vector<size_type>(n_blocks, block_size),
479 std::vector<size_type>(n_blocks, local_size),
488 const MPI_Comm &communicator,
489 const std::vector<size_type> &local_sizes,
496 for (
unsigned int i=0; i<this->
n_blocks(); ++i)
498 local_sizes[i], fast);
511 for (
unsigned int i=0; i<this->
n_blocks(); ++i)
518 const MPI_Comm &communicator)
520 std::vector<size_type> sizes(parallel_partitioning.size());
521 for (
unsigned int i=0; i<parallel_partitioning.size(); ++i)
522 sizes[i] = parallel_partitioning[i].
size();
528 for (
unsigned int i=0; i<this->
n_blocks(); ++i)
529 block(i).reinit(parallel_partitioning[i], communicator);
535 const std::vector<IndexSet> &ghost_entries,
536 const MPI_Comm &communicator)
538 std::vector<types::global_dof_index> sizes(parallel_partitioning.size());
539 for (
unsigned int i=0; i<parallel_partitioning.size(); ++i)
540 sizes[i] = parallel_partitioning[i].
size();
546 for (
unsigned int i=0; i<this->
n_blocks(); ++i)
547 block(i).reinit(parallel_partitioning[i], ghost_entries[i], communicator);
556 return block(0).get_mpi_communicator();
563 bool ghosted=
block(0).has_ghost_elements();
565 for (
unsigned int i=0; i<this->
n_blocks(); ++i)
579 for (
unsigned int i=0; i<this->
n_blocks(); ++i)
589 const unsigned int precision,
590 const bool scientific,
591 const bool across)
const
593 for (
unsigned int i=0; i<this->
n_blocks(); ++i)
596 out <<
'C' << i <<
':';
598 out <<
"Component " << i << std::endl;
599 this->
components[i].print(out, precision, scientific, across);
624 DEAL_II_NAMESPACE_CLOSE
626 #endif // DEAL_II_WITH_PETSC
BlockIndices block_indices
void swap(BlockVector &v)
void swap(BlockVector &u, BlockVector &v)
const BlockIndices & get_block_indices() const
DeclException0(ExcIteratorRangeDoesNotMatchVectorSize)
bool has_ghost_elements() const
BaseClass::value_type value_type
#define Assert(cond, exc)
BlockType & block(const unsigned int i)
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
unsigned int n_blocks() const
void reinit(const unsigned int n_blocks, const MPI_Comm &communicator, const size_type block_size, const size_type local_size, const bool fast=false)
BaseClass::BlockType BlockType
BlockVectorBase< Vector > BaseClass
const MPI_Comm & get_mpi_communicator() const
BlockVectorBase & operator=(const value_type s)
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
BlockVector & operator=(const value_type s)
::ExceptionBase & ExcInternalError()
std::vector< VectorType > components