Reference documentation for deal.II version 8.4.1
parallel_block_vector.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 2016 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 at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii__parallel_block_vector_h
17 #define dealii__parallel_block_vector_h
18 
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/base/exceptions.h>
22 #include <deal.II/lac/block_indices.h>
23 #include <deal.II/lac/block_vector_base.h>
24 #include <deal.II/lac/parallel_vector.h>
25 
26 #include <deal.II/lac/petsc_parallel_block_vector.h>
27 #include <deal.II/lac/trilinos_parallel_block_vector.h>
28 
29 
30 #include <cstdio>
31 #include <vector>
32 
33 DEAL_II_NAMESPACE_OPEN
34 
35 
36 namespace parallel
37 {
38  namespace distributed
39  {
40 
62  template <typename Number>
63  class BlockVector : public BlockVectorBase<Vector<Number> >
64  {
65  public:
70 
74  typedef typename BaseClass::BlockType BlockType;
75 
79  typedef typename BaseClass::value_type value_type;
80  typedef typename BaseClass::real_type real_type;
81  typedef typename BaseClass::pointer pointer;
82  typedef typename BaseClass::const_pointer const_pointer;
83  typedef typename BaseClass::reference reference;
84  typedef typename BaseClass::const_reference const_reference;
85  typedef typename BaseClass::size_type size_type;
86  typedef typename BaseClass::iterator iterator;
88 
99  explicit BlockVector (const size_type num_blocks = 0,
100  const size_type block_size = 0);
101 
106  BlockVector (const BlockVector<Number> &V);
107 
108 
109 #ifndef DEAL_II_EXPLICIT_CONSTRUCTOR_BUG
110 
122  template <typename OtherNumber>
123  explicit
125 #endif
126 
131  BlockVector (const std::vector<size_type> &block_sizes);
132 
137  BlockVector (const std::vector<IndexSet> &local_ranges,
138  const std::vector<IndexSet> &ghost_indices,
139  const MPI_Comm communicator);
140 
144  BlockVector (const std::vector<IndexSet> &local_ranges,
145  const MPI_Comm communicator);
146 
150  ~BlockVector ();
151 
156  BlockVector &operator = (const value_type s);
157 
162  BlockVector &
163  operator= (const BlockVector &V);
164 
169  template <class Number2>
170  BlockVector &
172 
176  BlockVector &
177  operator= (const Vector<Number> &V);
178 
179 #ifdef DEAL_II_WITH_PETSC
180 
189 #endif
190 
191 #ifdef DEAL_II_WITH_TRILINOS
192 
201  operator = (const TrilinosWrappers::MPI::BlockVector &trilinos_vec);
202 #endif
203 
217  void reinit (const size_type num_blocks,
218  const size_type block_size = 0,
219  const bool omit_zeroing_entries = false);
220 
240  void reinit (const std::vector<size_type> &N,
241  const bool omit_zeroing_entries=false);
242 
257  template <typename Number2>
258  void reinit (const BlockVector<Number2> &V,
259  const bool omit_zeroing_entries=false);
260 
285  void compress (::VectorOperation::values operation);
286 
295  void update_ghost_values () const;
296 
305  void zero_out_ghosts ();
306 
310  bool has_ghost_elements() const;
311 
317  bool all_zero () const;
318 
329  bool is_non_negative () const;
330 
334  template <typename Number2>
335  bool operator == (const BlockVector<Number2> &v) const;
336 
340  template <typename Number2>
341  bool operator != (const BlockVector<Number2> &v) const;
342 
346  template <typename Number2>
347  Number operator * (const BlockVector<Number2> &V) const;
348 
353  real_type norm_sqr () const;
354 
358  Number mean_value () const;
359 
364  real_type l1_norm () const;
365 
370  real_type l2_norm () const;
371 
377  real_type lp_norm (const real_type p) const;
378 
383  real_type linfty_norm () const;
384 
401  Number add_and_dot (const Number a,
402  const BlockVector<Number> &V,
403  const BlockVector<Number> &W);
404 
409  template <class BlockVector2>
410  void scale (const BlockVector2 &v);
411 
429  void swap (BlockVector<Number> &v);
430 
439  DeclException0 (ExcIteratorRangeDoesNotMatchVectorSize);
441  };
442 
445 #ifndef DOXYGEN
446  /*----------------------- Inline functions ----------------------------------*/
447 
448 
449  template <typename Number>
450  inline
451  BlockVector<Number>::BlockVector (const size_type n_blocks,
452  const size_type block_size)
453  {
454  reinit (n_blocks, block_size);
455  }
456 
457 
458 
459  template <typename Number>
460  inline
461  BlockVector<Number>::BlockVector (const std::vector<size_type> &n)
462  {
463  reinit (n, false);
464  }
465 
466 
467  template <typename Number>
468  inline
469  BlockVector<Number>::BlockVector (const std::vector<IndexSet> &local_ranges,
470  const std::vector<IndexSet> &ghost_indices,
471  const MPI_Comm communicator)
472  {
473  std::vector<size_type> sizes(local_ranges.size());
474  for (unsigned int i=0; i<local_ranges.size(); ++i)
475  sizes[i] = local_ranges[i].size();
476 
477  this->block_indices.reinit(sizes);
478  this->components.resize(this->n_blocks());
479 
480  for (unsigned int i=0; i<this->n_blocks(); ++i)
481  this->block(i).reinit(local_ranges[i], ghost_indices[i], communicator);
482  }
483 
484 
485  template <typename Number>
486  inline
487  BlockVector<Number>::BlockVector (const std::vector<IndexSet> &local_ranges,
488  const MPI_Comm communicator)
489  {
490  std::vector<size_type> sizes(local_ranges.size());
491  for (unsigned int i=0; i<local_ranges.size(); ++i)
492  sizes[i] = local_ranges[i].size();
493 
494  this->block_indices.reinit(sizes);
495  this->components.resize(this->n_blocks());
496 
497  for (unsigned int i=0; i<this->n_blocks(); ++i)
498  this->block(i).reinit(local_ranges[i], communicator);
499  }
500 
501 
502 
503  template <typename Number>
504  inline
506  :
507  BlockVectorBase<Vector<Number> > ()
508  {
509  this->components.resize (v.n_blocks());
510  this->block_indices = v.block_indices;
511 
512  for (size_type i=0; i<this->n_blocks(); ++i)
513  this->components[i] = v.components[i];
514  }
515 
516 
517 #ifndef DEAL_II_EXPLICIT_CONSTRUCTOR_BUG
518 
519  template <typename Number>
520  template <typename OtherNumber>
521  inline
523  {
524  reinit (v, true);
525  *this = v;
526  }
527 
528 #endif
529 
530 
531 
532  template <typename Number>
533  inline
534  void BlockVector<Number>::reinit (const size_type n_bl,
535  const size_type bl_sz,
536  const bool omit_zeroing_entries)
537  {
538  std::vector<size_type> n(n_bl, bl_sz);
539  reinit(n, omit_zeroing_entries);
540  }
541 
542 
543  template <typename Number>
544  inline
545  void BlockVector<Number>::reinit (const std::vector<size_type> &n,
546  const bool omit_zeroing_entries)
547  {
548  this->block_indices.reinit (n);
549  if (this->components.size() != this->n_blocks())
550  this->components.resize(this->n_blocks());
551 
552  for (size_type i=0; i<this->n_blocks(); ++i)
553  this->components[i].reinit(n[i], omit_zeroing_entries);
554  }
555 
556 
557 
558  template <typename Number>
559  template <typename Number2>
560  inline
562  const bool omit_zeroing_entries)
563  {
564  this->block_indices = v.get_block_indices();
565  if (this->components.size() != this->n_blocks())
566  this->components.resize(this->n_blocks());
567 
568  for (unsigned int i=0; i<this->n_blocks(); ++i)
569  this->block(i).reinit(v.block(i), omit_zeroing_entries);
570  }
571 
572 
573 
574  template <typename Number>
575  inline
577  {}
578 
579 
580 
581  template <typename Number>
582  inline
584  BlockVector<Number>::operator = (const value_type s)
585  {
586 
587  AssertIsFinite(s);
588 
589  BaseClass::operator = (s);
590  return *this;
591  }
592 
593 
594 
595  template <typename Number>
596  inline
599  {
600  // we only allow assignment to vectors with the same number of blocks
601  // or to an empty BlockVector
602  Assert (this->n_blocks() == 0 || this->n_blocks() == v.n_blocks(),
603  ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
604 
605  if (this->n_blocks() != v.n_blocks())
606  reinit(v.n_blocks(), true);
607 
608  for (size_type i=0; i<this->n_blocks(); ++i)
609  this->components[i] = v.block(i);
610 
611  this->collect_sizes();
612  return *this;
613  }
614 
615 
616 
617  template <typename Number>
618  inline
621  {
622  BaseClass::operator = (v);
623  return *this;
624  }
625 
626 
627 
628  template <typename Number>
629  template <typename Number2>
630  inline
633  {
634  reinit (v, true);
635  BaseClass::operator = (v);
636  return *this;
637  }
638 
639 
640 
641 #ifdef DEAL_II_WITH_PETSC
642 
643  template <typename Number>
644  inline
647  {
648  AssertDimension(this->n_blocks(), petsc_vec.n_blocks());
649  for (unsigned int i=0; i<this->n_blocks(); ++i)
650  this->block(i) = petsc_vec.block(i);
651 
652  return *this;
653  }
654 
655 #endif
656 
657 
658 
659 #ifdef DEAL_II_WITH_TRILINOS
660 
661  template <typename Number>
662  inline
665  {
666  AssertDimension(this->n_blocks(), trilinos_vec.n_blocks());
667  for (unsigned int i=0; i<this->n_blocks(); ++i)
668  this->block(i) = trilinos_vec.block(i);
669 
670  return *this;
671  }
672 
673 #endif
674 
675 
676 
677  template <typename Number>
678  inline
679  void
680  BlockVector<Number>::compress (::VectorOperation::values operation)
681  {
682  // start all requests for all blocks before finishing the transfers as
683  // this saves repeated synchronizations
684  for (unsigned int block=0; block<this->n_blocks(); ++block)
685  this->block(block).compress_start(block*10 + 8273, operation);
686  for (unsigned int block=0; block<this->n_blocks(); ++block)
687  this->block(block).compress_finish(operation);
688  }
689 
690 
691 
692  template <typename Number>
693  inline
694  void
696  {
697  for (unsigned int block=0; block<this->n_blocks(); ++block)
698  this->block(block).update_ghost_values_start(block*10 + 9923);
699  for (unsigned int block=0; block<this->n_blocks(); ++block)
700  this->block(block).update_ghost_values_finish();
701  }
702 
703 
704 
705  template <typename Number>
706  inline
707  void
709  {
710  for (unsigned int block=0; block<this->n_blocks(); ++block)
711  this->block(block).zero_out_ghosts();
712  }
713 
714 
715 
716  template <typename Number>
717  inline
718  bool
720  {
721  bool has_ghost_elements = false;
722  for (unsigned int block=0; block<this->n_blocks(); ++block)
723  if (this->block(block).has_ghost_elements() == true)
724  has_ghost_elements = true;
725  return has_ghost_elements;
726  }
727 
728 
729 
730  template <typename Number>
731  inline
732  bool
734  {
735  Assert (this->n_blocks() > 0, ExcEmptyObject());
736 
737  // use int instead of bool. in order to make global reduction operations
738  // work also when MPI_Init was not called, only call MPI_Allreduce
739  // commands when there is more than one processor (note that reinit()
740  // functions handle this case correctly through the job_supports_mpi()
741  // query). this is the same in all the functions below
742  int local_result = -1;
743  for (unsigned int i=0; i<this->n_blocks(); ++i)
744  local_result = std::max(local_result,
745  -static_cast<int>(this->block(i).all_zero_local()));
746 
747  if (this->block(0).partitioner->n_mpi_processes() > 1)
748  return -Utilities::MPI::max(local_result,
749  this->block(0).partitioner->get_communicator());
750  else
751  return local_result;
752  }
753 
754 
755 
756  template <typename Number>
757  inline
758  bool
760  {
761  Assert (this->n_blocks() > 0, ExcEmptyObject());
762  int local_result = -1;
763  for (unsigned int i=0; i<this->n_blocks(); ++i)
764  local_result = std::max(local_result,
765  -static_cast<int>(this->block(i).is_non_negative_local()));
766  if (this->block(0).partitioner->n_mpi_processes() > 1)
767  return Utilities::MPI::max(local_result,
768  this->block(0).partitioner->get_communicator());
769  else
770  return local_result;
771  }
772 
773 
774 
775  template <typename Number>
776  template <typename Number2>
777  inline
778  bool
780  {
781  Assert (this->n_blocks() > 0, ExcEmptyObject());
782  AssertDimension (this->n_blocks(), v.n_blocks());
783 
784  // MPI does not support bools, so use unsigned int instead. Two vectors
785  // are equal if the check for non-equal fails on all processors
786  unsigned int local_result = 0;
787  for (unsigned int i=0; i<this->n_blocks(); ++i)
788  local_result = std::max(local_result,
789  static_cast<unsigned int>(!this->block(i).vectors_equal_local(v.block(i))));
790  unsigned int result =
791  this->block(0).partitioner->n_mpi_processes() > 1
792  ?
793  Utilities::MPI::max(local_result, this->block(0).partitioner->get_communicator())
794  :
795  local_result;
796  return result==0;
797  }
798 
799 
800 
801  template <typename Number>
802  template <typename Number2>
803  inline
804  bool
806  {
807  return !(operator == (v));
808  }
809 
810 
811 
812  template <typename Number>
813  template <typename Number2>
814  inline
815  Number
817  {
818  Assert (this->n_blocks() > 0, ExcEmptyObject());
819  AssertDimension (this->n_blocks(), v.n_blocks());
820 
821  Number local_result = Number();
822  for (unsigned int i=0; i<this->n_blocks(); ++i)
823  local_result += this->block(i).inner_product_local(v.block(i));
824 
825  if (this->block(0).partitioner->n_mpi_processes() > 1)
826  return Utilities::MPI::sum (local_result,
827  this->block(0).partitioner->get_communicator());
828  else
829  return local_result;
830  }
831 
832 
833 
834  template <typename Number>
835  inline
836  typename BlockVector<Number>::real_type
838  {
839  Assert (this->n_blocks() > 0, ExcEmptyObject());
840 
841  real_type local_result = real_type();
842  for (unsigned int i=0; i<this->n_blocks(); ++i)
843  local_result += this->block(i).norm_sqr_local();
844 
845  if (this->block(0).partitioner->n_mpi_processes() > 1)
846  return Utilities::MPI::sum (local_result,
847  this->block(0).partitioner->get_communicator());
848  else
849  return local_result;
850  }
851 
852 
853 
854  template <typename Number>
855  inline
856  Number
858  {
859  Assert (this->n_blocks() > 0, ExcEmptyObject());
860 
861  Number local_result = Number();
862  for (unsigned int i=0; i<this->n_blocks(); ++i)
863  local_result += this->block(i).mean_value_local()*(real_type)this->block(i).partitioner->local_size();
864 
865  if (this->block(0).partitioner->n_mpi_processes() > 1)
866  return Utilities::MPI::sum (local_result,
867  this->block(0).partitioner->get_communicator())/
868  (real_type)this->size();
869  else
870  return local_result/(real_type)this->size();
871  }
872 
873 
874 
875  template <typename Number>
876  inline
877  typename BlockVector<Number>::real_type
879  {
880  Assert (this->n_blocks() > 0, ExcEmptyObject());
881 
882  real_type local_result = real_type();
883  for (unsigned int i=0; i<this->n_blocks(); ++i)
884  local_result += this->block(i).l1_norm_local();
885 
886  if (this->block(0).partitioner->n_mpi_processes() > 1)
887  return Utilities::MPI::sum (local_result,
888  this->block(0).partitioner->get_communicator());
889  else
890  return local_result;
891  }
892 
893 
894 
895  template <typename Number>
896  inline
897  typename BlockVector<Number>::real_type
899  {
900  return std::sqrt(norm_sqr());
901  }
902 
903 
904 
905  template <typename Number>
906  inline
907  typename BlockVector<Number>::real_type
908  BlockVector<Number>::lp_norm (const real_type p) const
909  {
910  Assert (this->n_blocks() > 0, ExcEmptyObject());
911 
912  real_type local_result = real_type();
913  for (unsigned int i=0; i<this->n_blocks(); ++i)
914  local_result += std::pow(this->block(i).lp_norm_local(p), p);
915 
916  if (this->block(0).partitioner->n_mpi_processes() > 1)
917  return std::pow (Utilities::MPI::sum(local_result,
918  this->block(0).partitioner->get_communicator()),
919  static_cast<real_type>(1.0/p));
920  else
921  return std::pow (local_result, static_cast<real_type>(1.0/p));
922  }
923 
924 
925 
926  template <typename Number>
927  inline
928  typename BlockVector<Number>::real_type
930  {
931  Assert (this->n_blocks() > 0, ExcEmptyObject());
932 
933  real_type local_result = real_type();
934  for (unsigned int i=0; i<this->n_blocks(); ++i)
935  local_result = std::max(local_result, this->block(i).linfty_norm_local());
936 
937  if (this->block(0).partitioner->n_mpi_processes() > 1)
938  return Utilities::MPI::max (local_result,
939  this->block(0).partitioner->get_communicator());
940  else
941  return local_result;
942  }
943 
944 
945 
946  template <typename Number>
947  inline
948  Number
949  BlockVector<Number>::add_and_dot (const Number a,
950  const BlockVector<Number> &V,
951  const BlockVector<Number> &W)
952  {
953  Number local_result = Number();
954  for (unsigned int i=0; i<this->n_blocks(); ++i)
955  local_result += this->block(i).add_and_dot_local(a, V.block(i), W.block(i));
956 
957  if (this->block(0).partitioner->n_mpi_processes() > 1)
958  return Utilities::MPI::sum (local_result,
959  this->block(0).partitioner->get_communicator());
960  else
961  return local_result;
962  }
963 
964 
965 
966  template <typename Number>
967  inline
968  void
970  {
971  Assert (this->n_blocks() == v.n_blocks(),
972  ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
973 
974  for (size_type i=0; i<this->n_blocks(); ++i)
975  ::swap (this->components[i], v.components[i]);
976  ::swap (this->block_indices, v.block_indices);
977  }
978 
979 
980 
981  template <typename Number>
982  template <class BlockVector2>
983  void BlockVector<Number>::scale (const BlockVector2 &v)
984  {
985  BaseClass::scale (v);
986  }
987 
988 #endif // DOXYGEN
989 
990  } // end of namespace distributed
991 
992 } // end of namespace parallel
993 
1002 template <typename Number>
1003 inline
1006 {
1007  u.swap (v);
1008 }
1009 
1010 DEAL_II_NAMESPACE_CLOSE
1011 
1012 #endif
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1052
void scale(const BlockVector2 &v)
bool operator==(const BlockVector< Number2 > &v) const
void swap(parallel::distributed::BlockVector< Number > &u, parallel::distributed::BlockVector< Number > &v)
BlockVector & operator=(const value_type s)
void compress(::VectorOperation::values operation)
DeclException0(ExcIteratorRangeDoesNotMatchVectorSize)
BlockVectorBase< Vector< Number > > BaseClass
T sum(const T &t, const MPI_Comm &mpi_communicator)
Definition: mpi.h:607
#define Assert(cond, exc)
Definition: exceptions.h:294
const BlockIndices & get_block_indices() const
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
Number add_and_dot(const Number a, const BlockVector< Number > &V, const BlockVector< Number > &W)
Number operator*(const BlockVector< Number2 > &V) const
void reinit(const size_type num_blocks, const size_type block_size=0, const bool omit_zeroing_entries=false)
void swap(BlockVector< Number > &v)
std::vector< Vector< Number > > components
real_type lp_norm(const real_type p) const
bool operator!=(const BlockVector< Number2 > &v) const
BlockType & block(const unsigned int i)
#define AssertIsFinite(number)
Definition: exceptions.h:1096
T max(const T &t, const MPI_Comm &mpi_communicator)
Definition: mpi.h:689
BlockVector(const size_type num_blocks=0, const size_type block_size=0)