Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
petsc_parallel_block_vector.h
1 // ---------------------------------------------------------------------
2 // @f$Id: petsc_parallel_block_vector.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 2004 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__petsc_parallel_block_vector_h
18 #define __deal2__petsc_parallel_block_vector_h
19 
20 
21 #include <deal.II/base/config.h>
22 
23 #ifdef DEAL_II_WITH_PETSC
24 
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>
29 
30 DEAL_II_NAMESPACE_OPEN
31 
32 
33 namespace PETScWrappers
34 {
35  // forward declaration
36  class BlockVector;
37 
38  namespace MPI
39  {
40 
62  class BlockVector : public BlockVectorBase<Vector>
63  {
64  public:
70 
76 
81  typedef BaseClass::value_type value_type;
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;
89 
94  BlockVector ();
95 
105  explicit BlockVector (const unsigned int n_blocks,
106  const MPI_Comm &communicator,
107  const size_type block_size,
108  const size_type local_size);
109 
116  BlockVector (const BlockVector &V);
117 
130  BlockVector (const std::vector<size_type> &block_sizes,
131  const MPI_Comm &communicator,
132  const std::vector<size_type> &local_elements);
133 
138  explicit BlockVector (const std::vector<IndexSet> &parallel_partitioning,
139  const MPI_Comm &communicator = MPI_COMM_WORLD);
140 
144  BlockVector (const std::vector<IndexSet> &parallel_partitioning,
145  const std::vector<IndexSet> &ghost_indices,
146  const MPI_Comm &communicator);
147 
148 
149 
153  ~BlockVector ();
154 
161 
166  BlockVector &
167  operator= (const BlockVector &V);
168 
198  BlockVector &
200 
214  void reinit (const unsigned int n_blocks,
215  const MPI_Comm &communicator,
216  const size_type block_size,
217  const size_type local_size,
218  const bool fast = false);
219 
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);
259 
286  void reinit (const BlockVector &V,
287  const bool fast=false);
288 
293  void reinit (const std::vector<IndexSet> &parallel_partitioning,
294  const MPI_Comm &communicator);
295 
299  void reinit (const std::vector<IndexSet> &parallel_partitioning,
300  const std::vector<IndexSet> &ghost_entries,
301  const MPI_Comm &communicator);
302 
314  void reinit (const unsigned int num_blocks);
315 
319  bool has_ghost_elements() const;
320 
326  const MPI_Comm &get_mpi_communicator () const;
327 
358  void swap (BlockVector &v);
359 
363  void print (std::ostream &out,
364  const unsigned int precision = 3,
365  const bool scientific = true,
366  const bool across = true) const;
367 
371  DeclException0 (ExcIteratorRangeDoesNotMatchVectorSize);
375  DeclException0 (ExcNonMatchingBlockVectors);
376  };
377 
380  /*----------------------- Inline functions ----------------------------------*/
381 
382 
383  inline
385  {}
386 
387 
388 
389  inline
390  BlockVector::BlockVector (const unsigned int n_blocks,
391  const MPI_Comm &communicator,
392  const size_type block_size,
393  const size_type local_size)
394  {
395  reinit (n_blocks, communicator, block_size, local_size);
396  }
397 
398 
399 
400  inline
401  BlockVector::BlockVector (const std::vector<size_type> &block_sizes,
402  const MPI_Comm &communicator,
403  const std::vector<size_type> &local_elements)
404  {
405  reinit (block_sizes, communicator, local_elements, false);
406  }
407 
408 
409  inline
411  :
413  {
414  this->components.resize (v.n_blocks());
415  this->block_indices = v.block_indices;
416 
417  for (unsigned int i=0; i<this->n_blocks(); ++i)
418  this->components[i] = v.components[i];
419  }
420 
421  inline
422  BlockVector::BlockVector (const std::vector<IndexSet> &parallel_partitioning,
423  const MPI_Comm &communicator)
424  {
425  reinit(parallel_partitioning, communicator);
426  }
427 
428  inline
429  BlockVector::BlockVector (const std::vector<IndexSet> &parallel_partitioning,
430  const std::vector<IndexSet> &ghost_indices,
431  const MPI_Comm &communicator)
432  {
433  reinit(parallel_partitioning, ghost_indices, communicator);
434  }
435 
436  inline
437  BlockVector &
439  {
441  return *this;
442  }
443 
444  inline
445  BlockVector &
447  {
448  // we only allow assignment to vectors with the same number of blocks
449  // or to an empty BlockVector
450  Assert (n_blocks() == 0 || n_blocks() == v.n_blocks(),
452 
453  if (this->n_blocks() != v.n_blocks())
454  reinit(v.n_blocks());
455 
456  for (size_type i=0; i<this->n_blocks(); ++i)
457  this->components[i] = v.block(i);
458 
459  collect_sizes();
460 
461  return *this;
462  }
463 
464  inline
466  {}
467 
468 
469  inline
470  void
471  BlockVector::reinit (const unsigned int n_blocks,
472  const MPI_Comm &communicator,
473  const size_type block_size,
474  const size_type local_size,
475  const bool fast)
476  {
477  reinit(std::vector<size_type>(n_blocks, block_size),
478  communicator,
479  std::vector<size_type>(n_blocks, local_size),
480  fast);
481  }
482 
483 
484 
485  inline
486  void
487  BlockVector::reinit (const std::vector<size_type> &block_sizes,
488  const MPI_Comm &communicator,
489  const std::vector<size_type> &local_sizes,
490  const bool fast)
491  {
492  this->block_indices.reinit (block_sizes);
493  if (this->components.size() != this->n_blocks())
494  this->components.resize(this->n_blocks());
495 
496  for (unsigned int i=0; i<this->n_blocks(); ++i)
497  this->components[i].reinit(communicator, block_sizes[i],
498  local_sizes[i], fast);
499  }
500 
501 
502  inline
503  void
505  const bool fast)
506  {
507  this->block_indices = v.get_block_indices();
508  if (this->components.size() != this->n_blocks())
509  this->components.resize(this->n_blocks());
510 
511  for (unsigned int i=0; i<this->n_blocks(); ++i)
512  block(i).reinit(v.block(i), fast);
513  }
514 
515  inline
516  void
517  BlockVector::reinit (const std::vector<IndexSet> &parallel_partitioning,
518  const MPI_Comm &communicator)
519  {
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();
523 
524  this->block_indices.reinit(sizes);
525  if (this->components.size() != this->n_blocks())
526  this->components.resize(this->n_blocks());
527 
528  for (unsigned int i=0; i<this->n_blocks(); ++i)
529  block(i).reinit(parallel_partitioning[i], communicator);
530  }
531 
532  inline
533  void
534  BlockVector::reinit (const std::vector<IndexSet> &parallel_partitioning,
535  const std::vector<IndexSet> &ghost_entries,
536  const MPI_Comm &communicator)
537  {
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();
541 
542  this->block_indices.reinit(sizes);
543  if (this->components.size() != this->n_blocks())
544  this->components.resize(this->n_blocks());
545 
546  for (unsigned int i=0; i<this->n_blocks(); ++i)
547  block(i).reinit(parallel_partitioning[i], ghost_entries[i], communicator);
548  }
549 
550 
551 
552  inline
553  const MPI_Comm &
555  {
556  return block(0).get_mpi_communicator();
557  }
558 
559  inline
560  bool
562  {
563  bool ghosted=block(0).has_ghost_elements();
564 #ifdef DEBUG
565  for (unsigned int i=0; i<this->n_blocks(); ++i)
566  Assert(block(i).has_ghost_elements()==ghosted, ExcInternalError());
567 #endif
568  return ghosted;
569  }
570 
571 
572  inline
573  void
575  {
576  Assert (this->n_blocks() == v.n_blocks(),
577  ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
578 
579  for (unsigned int i=0; i<this->n_blocks(); ++i)
580  this->components[i].swap (v.components[i]);
582  }
583 
584 
585 
586  inline
587  void
588  BlockVector::print (std::ostream &out,
589  const unsigned int precision,
590  const bool scientific,
591  const bool across) const
592  {
593  for (unsigned int i=0; i<this->n_blocks(); ++i)
594  {
595  if (across)
596  out << 'C' << i << ':';
597  else
598  out << "Component " << i << std::endl;
599  this->components[i].print(out, precision, scientific, across);
600  }
601  }
602 
603 
604 
613  inline
614  void swap (BlockVector &u,
615  BlockVector &v)
616  {
617  u.swap (v);
618  }
619 
620  }
621 
622 }
623 
624 DEAL_II_NAMESPACE_CLOSE
625 
626 #endif // DEAL_II_WITH_PETSC
627 
628 #endif
void swap(BlockVector &u, BlockVector &v)
const BlockIndices & get_block_indices() const
DeclException0(ExcIteratorRangeDoesNotMatchVectorSize)
#define Assert(cond, exc)
Definition: exceptions.h:299
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
std::size_t size() 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)
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