Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
parallel_block_vector.h
1 // ---------------------------------------------------------------------
2 // @f$Id: parallel_block_vector.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 1999 - 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__parallel_block_vector_h
18 #define __deal2__parallel_block_vector_h
19 
20 
21 #include <deal.II/base/config.h>
23 #include <deal.II/lac/block_indices.h>
24 #include <deal.II/lac/block_vector_base.h>
25 #include <deal.II/lac/parallel_vector.h>
26 
27 #include <cstdio>
28 #include <vector>
29 
30 DEAL_II_NAMESPACE_OPEN
31 
32 
33 
34 namespace parallel
35 {
36  namespace distributed
37  {
38 
57  template <typename Number>
58  class BlockVector : public BlockVectorBase<Vector<Number> >
59  {
60  public:
65 
69  typedef typename BaseClass::BlockType BlockType;
70 
74  typedef typename BaseClass::value_type value_type;
75  typedef typename BaseClass::real_type real_type;
76  typedef typename BaseClass::pointer pointer;
77  typedef typename BaseClass::const_pointer const_pointer;
78  typedef typename BaseClass::reference reference;
79  typedef typename BaseClass::const_reference const_reference;
80  typedef typename BaseClass::size_type size_type;
81  typedef typename BaseClass::iterator iterator;
83 
94  explicit BlockVector (const size_type num_blocks = 0,
95  const size_type block_size = 0);
96 
101  BlockVector (const BlockVector<Number> &V);
102 
103 
104 #ifndef DEAL_II_EXPLICIT_CONSTRUCTOR_BUG
105 
117  template <typename OtherNumber>
118  explicit
120 #endif
121 
126  BlockVector (const std::vector<size_type> &block_sizes);
127 
132  BlockVector (const std::vector<IndexSet> &local_ranges,
133  const std::vector<IndexSet> &ghost_indices,
134  const MPI_Comm communicator);
135 
139  BlockVector (const std::vector<IndexSet> &local_ranges,
140  const MPI_Comm communicator);
141 
145  ~BlockVector ();
146 
152 
157  BlockVector &
158  operator= (const BlockVector &V);
159 
164  template <class Number2>
165  BlockVector &
167 
171  BlockVector &
172  operator= (const Vector<Number> &V);
173 
186  void reinit (const size_type num_blocks,
187  const size_type block_size = 0,
188  const bool fast = false);
189 
208  void reinit (const std::vector<size_type> &N,
209  const bool fast=false);
210 
225  template <typename Number2>
226  void reinit (const BlockVector<Number2> &V,
227  const bool fast=false);
228 
252  void compress (::VectorOperation::values operation);
253 
262  void update_ghost_values () const;
263 
272  void zero_out_ghosts ();
273 
277  bool has_ghost_elements() const;
278 
285  bool all_zero () const;
286 
297  bool is_non_negative () const;
298 
302  template <typename Number2>
303  bool operator == (const BlockVector<Number2> &v) const;
304 
308  template <typename Number2>
309  bool operator != (const BlockVector<Number2> &v) const;
310 
314  template <typename Number2>
315  Number operator * (const BlockVector<Number2> &V) const;
316 
321  real_type norm_sqr () const;
322 
326  Number mean_value () const;
327 
332  real_type l1_norm () const;
333 
338  real_type l2_norm () const;
339 
345  real_type lp_norm (const real_type p) const;
346 
351  real_type linfty_norm () const;
352 
361  void scale (const value_type factor) DEAL_II_DEPRECATED;
362 
367  template <class BlockVector2>
368  void scale (const BlockVector2 &v);
369 
387  void swap (BlockVector<Number> &v);
388 
395  DeclException0 (ExcIteratorRangeDoesNotMatchVectorSize);
397  };
398 
401 #ifndef DOXYGEN
402  /*----------------------- Inline functions ----------------------------------*/
403 
404 
405  template <typename Number>
406  inline
408  const size_type block_size)
409  {
410  reinit (n_blocks, block_size);
411  }
412 
413 
414 
415  template <typename Number>
416  inline
417  BlockVector<Number>::BlockVector (const std::vector<size_type> &n)
418  {
419  reinit (n, false);
420  }
421 
422 
423  template <typename Number>
424  inline
425  BlockVector<Number>::BlockVector (const std::vector<IndexSet> &local_ranges,
426  const std::vector<IndexSet> &ghost_indices,
427  const MPI_Comm communicator)
428  {
429  std::vector<size_type> sizes(local_ranges.size());
430  for (unsigned int i=0; i<local_ranges.size(); ++i)
431  sizes[i] = local_ranges[i].size();
432 
433  this->block_indices.reinit(sizes);
434  this->components.resize(this->n_blocks());
435 
436  for (unsigned int i=0; i<this->n_blocks(); ++i)
437  this->block(i).reinit(local_ranges[i], ghost_indices[i], communicator);
438  }
439 
440 
441  template <typename Number>
442  inline
443  BlockVector<Number>::BlockVector (const std::vector<IndexSet> &local_ranges,
444  const MPI_Comm communicator)
445  {
446  std::vector<size_type> sizes(local_ranges.size());
447  for (unsigned int i=0; i<local_ranges.size(); ++i)
448  sizes[i] = local_ranges[i].size();
449 
450  this->block_indices.reinit(sizes);
451  this->components.resize(this->n_blocks());
452 
453  for (unsigned int i=0; i<this->n_blocks(); ++i)
454  this->block(i).reinit(local_ranges[i], communicator);
455  }
456 
457 
458 
459  template <typename Number>
460  inline
462  :
463  BlockVectorBase<Vector<Number> > ()
464  {
465  this->components.resize (v.n_blocks());
466  this->block_indices = v.block_indices;
467 
468  for (size_type i=0; i<this->n_blocks(); ++i)
469  this->components[i] = v.components[i];
470  }
471 
472 
473 #ifndef DEAL_II_EXPLICIT_CONSTRUCTOR_BUG
474 
475  template <typename Number>
476  template <typename OtherNumber>
477  inline
479  {
480  reinit (v, true);
481  *this = v;
482  }
483 
484 #endif
485 
486 
487 
488  template <typename Number>
489  inline
490  void BlockVector<Number>::reinit (const size_type n_bl,
491  const size_type bl_sz,
492  const bool fast)
493  {
494  std::vector<size_type> n(n_bl, bl_sz);
495  reinit(n, fast);
496  }
497 
498 
499  template <typename Number>
500  inline
501  void BlockVector<Number>::reinit (const std::vector<size_type> &n,
502  const bool fast)
503  {
504  this->block_indices.reinit (n);
505  if (this->components.size() != this->n_blocks())
506  this->components.resize(this->n_blocks());
507 
508  for (size_type i=0; i<this->n_blocks(); ++i)
509  this->components[i].reinit(n[i], fast);
510  }
511 
512 
513 
514  template <typename Number>
515  template <typename Number2>
516  inline
518  const bool fast)
519  {
520  this->block_indices = v.get_block_indices();
521  if (this->components.size() != this->n_blocks())
522  this->components.resize(this->n_blocks());
523 
524  for (unsigned int i=0; i<this->n_blocks(); ++i)
525  this->block(i).reinit(v.block(i), fast);
526  }
527 
528 
529 
530  template <typename Number>
531  inline
533  {}
534 
535 
536 
537  template <typename Number>
538  inline
540  BlockVector<Number>::operator = (const value_type s)
541  {
542 
544 
545  BaseClass::operator = (s);
546  return *this;
547  }
548 
549 
550 
551  template <typename Number>
552  inline
555  {
556  // we only allow assignment to vectors with the same number of blocks
557  // or to an empty BlockVector
558  Assert (this->n_blocks() == 0 || this->n_blocks() == v.n_blocks(),
559  ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
560 
561  if (this->n_blocks() != v.n_blocks())
562  reinit(v.n_blocks(), true);
563 
564  for (size_type i=0; i<this->n_blocks(); ++i)
565  this->components[i] = v.block(i);
566 
567  this->collect_sizes();
568  return *this;
569  }
570 
571 
572 
573  template <typename Number>
574  inline
577  {
578  BaseClass::operator = (v);
579  return *this;
580  }
581 
582 
583 
584  template <typename Number>
585  template <typename Number2>
586  inline
589  {
590  reinit (v, true);
591  BaseClass::operator = (v);
592  return *this;
593  }
594 
595 
596 
597  template <typename Number>
598  inline
599  void
600  BlockVector<Number>::compress (::VectorOperation::values operation)
601  {
602  // start all requests for all blocks before finishing the transfers as
603  // this saves repeated synchronizations
604  for (unsigned int block=0; block<this->n_blocks(); ++block)
605  this->block(block).compress_start(block*10 + 8273, operation);
606  for (unsigned int block=0; block<this->n_blocks(); ++block)
607  this->block(block).compress_finish(operation);
608  }
609 
610 
611 
612  template <typename Number>
613  inline
614  void
616  {
617  for (unsigned int block=0; block<this->n_blocks(); ++block)
618  this->block(block).update_ghost_values_start(block*10 + 9923);
619  for (unsigned int block=0; block<this->n_blocks(); ++block)
620  this->block(block).update_ghost_values_finish();
621  }
622 
623 
624 
625  template <typename Number>
626  inline
627  void
629  {
630  for (unsigned int block=0; block<this->n_blocks(); ++block)
631  this->block(block).zero_out_ghosts();
632  }
633 
634 
635 
636  template <typename Number>
637  inline
638  bool
640  {
641  bool has_ghost_elements = false;
642  for (unsigned int block=0; block<this->n_blocks(); ++block)
643  if (this->block(block).has_ghost_elements() == true)
644  has_ghost_elements = true;
645  return has_ghost_elements;
646  }
647 
648 
649 
650  template <typename Number>
651  inline
652  bool
654  {
655  Assert (this->n_blocks() > 0, ExcEmptyObject());
656 
657  // use int instead of bool. in order to make global reduction operations
658  // work also when MPI_Init was not called, only call MPI_Allreduce
659  // commands when there is more than one processor (note that reinit()
660  // functions handle this case correctly through the job_supports_mpi()
661  // query). this is the same in all the functions below
662  int local_result = -1;
663  for (unsigned int i=0; i<this->n_blocks(); ++i)
664  local_result = std::max(local_result,
665  -static_cast<int>(this->block(i).all_zero_local()));
666 
667  if (this->block(0).partitioner->n_mpi_processes() > 1)
668  return -Utilities::MPI::max(local_result,
669  this->block(0).partitioner->get_communicator());
670  else
671  return local_result;
672  }
673 
674 
675 
676  template <typename Number>
677  inline
678  bool
680  {
681  Assert (this->n_blocks() > 0, ExcEmptyObject());
682  int local_result = -1;
683  for (unsigned int i=0; i<this->n_blocks(); ++i)
684  local_result = std::max(local_result,
685  -static_cast<int>(this->block(i).is_non_negative_local()));
686  if (this->block(0).partitioner->n_mpi_processes() > 1)
687  return Utilities::MPI::max(local_result,
688  this->block(0).partitioner->get_communicator());
689  else
690  return local_result;
691  }
692 
693 
694 
695  template <typename Number>
696  template <typename Number2>
697  inline
698  bool
700  {
701  Assert (this->n_blocks() > 0, ExcEmptyObject());
702  AssertDimension (this->n_blocks(), v.n_blocks());
703 
704  // MPI does not support bools, so use unsigned int instead. Two vectors
705  // are equal if the check for non-equal fails on all processors
706  unsigned int local_result = 0;
707  for (unsigned int i=0; i<this->n_blocks(); ++i)
708  local_result = std::max(local_result,
709  static_cast<unsigned int>(!this->block(i).vectors_equal_local(v.block(i))));
710  unsigned int result =
711  this->block(0).partitioner->n_mpi_processes() > 1
712  ?
713  Utilities::MPI::max(local_result, this->block(0).partitioner->get_communicator())
714  :
715  local_result;
716  return result==0;
717  }
718 
719 
720 
721  template <typename Number>
722  template <typename Number2>
723  inline
724  bool
726  {
727  return !(operator == (v));
728  }
729 
730 
731 
732  template <typename Number>
733  template <typename Number2>
734  inline
735  Number
737  {
738  Assert (this->n_blocks() > 0, ExcEmptyObject());
739  AssertDimension (this->n_blocks(), v.n_blocks());
740 
741  Number local_result = Number();
742  for (unsigned int i=0; i<this->n_blocks(); ++i)
743  local_result += this->block(i).inner_product_local(v.block(i));
744 
745  if (this->block(0).partitioner->n_mpi_processes() > 1)
746  return Utilities::MPI::sum (local_result,
747  this->block(0).partitioner->get_communicator());
748  else
749  return local_result;
750  }
751 
752 
753 
754  template <typename Number>
755  inline
756  typename BlockVector<Number>::real_type
758  {
759  Assert (this->n_blocks() > 0, ExcEmptyObject());
760 
761  real_type local_result = real_type();
762  for (unsigned int i=0; i<this->n_blocks(); ++i)
763  local_result += this->block(i).norm_sqr_local();
764 
765  if (this->block(0).partitioner->n_mpi_processes() > 1)
766  return Utilities::MPI::sum (local_result,
767  this->block(0).partitioner->get_communicator());
768  else
769  return local_result;
770  }
771 
772 
773 
774  template <typename Number>
775  inline
776  Number
778  {
779  Assert (this->n_blocks() > 0, ExcEmptyObject());
780 
781  Number local_result = Number();
782  for (unsigned int i=0; i<this->n_blocks(); ++i)
783  local_result += this->block(i).mean_value_local()*(real_type)this->block(i).partitioner->local_size();
784 
785  if (this->block(0).partitioner->n_mpi_processes() > 1)
786  return Utilities::MPI::sum (local_result,
787  this->block(0).partitioner->get_communicator())/
788  (real_type)this->size();
789  else
790  return local_result/(real_type)this->size();
791  }
792 
793 
794 
795  template <typename Number>
796  inline
797  typename BlockVector<Number>::real_type
799  {
800  Assert (this->n_blocks() > 0, ExcEmptyObject());
801 
802  real_type local_result = real_type();
803  for (unsigned int i=0; i<this->n_blocks(); ++i)
804  local_result += this->block(i).l1_norm_local();
805 
806  if (this->block(0).partitioner->n_mpi_processes() > 1)
807  return Utilities::MPI::sum (local_result,
808  this->block(0).partitioner->get_communicator());
809  else
810  return local_result;
811  }
812 
813 
814 
815  template <typename Number>
816  inline
817  typename BlockVector<Number>::real_type
819  {
820  return std::sqrt(norm_sqr());
821  }
822 
823 
824 
825  template <typename Number>
826  inline
827  typename BlockVector<Number>::real_type
828  BlockVector<Number>::lp_norm (const real_type p) const
829  {
830  Assert (this->n_blocks() > 0, ExcEmptyObject());
831 
832  real_type local_result = real_type();
833  for (unsigned int i=0; i<this->n_blocks(); ++i)
834  local_result += std::pow(this->block(i).lp_norm_local(p), p);
835 
836  if (this->block(0).partitioner->n_mpi_processes() > 1)
837  return std::pow (Utilities::MPI::sum(local_result,
838  this->block(0).partitioner->get_communicator()),
839  static_cast<real_type>(1.0/p));
840  else
841  return std::pow (local_result, static_cast<real_type>(1.0/p));
842  }
843 
844 
845 
846  template <typename Number>
847  inline
848  typename BlockVector<Number>::real_type
850  {
851  Assert (this->n_blocks() > 0, ExcEmptyObject());
852 
853  real_type local_result = real_type();
854  for (unsigned int i=0; i<this->n_blocks(); ++i)
855  local_result = std::max(local_result, this->block(i).linfty_norm_local());
856 
857  if (this->block(0).partitioner->n_mpi_processes() > 1)
858  return Utilities::MPI::max (local_result,
859  this->block(0).partitioner->get_communicator());
860  else
861  return local_result;
862  }
863 
864 
865 
866  template <typename Number>
867  inline
869  {
870  Assert (this->n_blocks() == v.n_blocks(),
871  ExcDimensionMismatch(this->n_blocks(), v.n_blocks()));
872 
873  for (size_type i=0; i<this->n_blocks(); ++i)
874  ::swap (this->components[i], v.components[i]);
875  ::swap (this->block_indices, v.block_indices);
876  }
877 
878 
879 
880  template <typename Number>
881  void BlockVector<Number>::scale (const value_type factor)
882  {
883 
885 
886  for (size_type i=0; i<this->n_blocks(); ++i)
887  this->components[i].scale(factor);
888  }
889 
890 
891 
892  template <typename Number>
893  template <class BlockVector2>
894  void BlockVector<Number>::scale (const BlockVector2 &v)
895  {
896  BaseClass::scale (v);
897  }
898 
899 #endif // DOXYGEN
900 
901  } // end of namespace distributed
902 
903 } // end of namespace parallel
904 
913 template <typename Number>
914 inline
917 {
918  u.swap (v);
919 }
920 
921 DEAL_II_NAMESPACE_CLOSE
922 
923 #endif
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:858
void reinit(const size_type num_blocks, const size_type block_size=0, const bool fast=false)
bool operator==(const BlockVector< Number2 > &v) const
void scale(const value_type factor) DEAL_II_DEPRECATED
void swap(parallel::distributed::BlockVector< Number > &u, parallel::distributed::BlockVector< Number > &v)
void scale(const double scaling_factor, Triangulation< dim, spacedim > &triangulation)
BlockVector & operator=(const value_type s)
const BlockIndices & get_block_indices() const
DeclException0(ExcIteratorRangeDoesNotMatchVectorSize)
BlockVectorBase< Vector< Number > > BaseClass
bool is_finite(const double x)
void compress() DEAL_II_DEPRECATED
T sum(const T &t, const MPI_Comm &mpi_communicator)
Definition: mpi.h:447
#define Assert(cond, exc)
Definition: exceptions.h:299
::ExceptionBase & ExcEmptyObject()
BlockType & block(const unsigned int i)
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
Number operator*(const BlockVector< Number2 > &V) const
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
void swap(BlockVector< Number > &v)
real_type lp_norm(const real_type p) const
::ExceptionBase & ExcNumberNotFinite()
bool operator!=(const BlockVector< Number2 > &v) const
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
std::vector< VectorType > components
T max(const T &t, const MPI_Comm &mpi_communicator)
Definition: mpi.h:501
BlockVector(const size_type num_blocks=0, const size_type block_size=0)