Reference documentation for deal.II version GIT 6d02bd1105 2022-05-22 03:35:01+00:00
\(\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\}}\)
multigrid.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 2020 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_multigrid_h
17 #define dealii_multigrid_h
18 
19 
20 #include <deal.II/base/config.h>
21 
25 
27 
29 
31 #include <deal.II/lac/vector.h>
32 
34 
35 #include <vector>
36 
38 
39 #ifdef signals
40 # error \
41  "The name 'signals' is already defined. You are most likely using the QT library \
42 and using the 'signals' keyword. You can either #include the Qt headers (or any conflicting headers) \
43 *after* the deal.II headers or you can define the 'QT_NO_KEYWORDS' macro and use the 'Q_SIGNALS' macro."
44 #endif
45 
48 
49 namespace mg
50 {
60  struct Signals
61  {
67  boost::signals2::signal<void(const bool before)> transfer_to_mg;
68 
74  boost::signals2::signal<void(const bool before)> transfer_to_global;
75 
85  boost::signals2::signal<void(const bool before, const unsigned int level)>
87 
96  boost::signals2::signal<void(const bool before, const unsigned int level)>
98 
104  boost::signals2::signal<void(const bool before, const unsigned int level)>
106 
115  boost::signals2::signal<void(const bool before, const unsigned int level)>
117 
123  boost::signals2::signal<void(const bool before, const unsigned int level)>
125 
131  boost::signals2::signal<void(const bool before, const unsigned int level)>
133 
138  boost::signals2::signal<void(const bool before, const unsigned int level)>
140  };
141 } // namespace mg
142 
160 template <typename VectorType>
161 class Multigrid : public Subscriptor
162 {
163 public:
167  enum Cycle
168  {
174  f_cycle
175  };
176 
177  using vector_type = VectorType;
178  using const_vector_type = const VectorType;
179 
194  const unsigned int minlevel = 0,
195  const unsigned int maxlevel = numbers::invalid_unsigned_int,
196  Cycle cycle = v_cycle);
197 
201  void
202  reinit(const unsigned int minlevel, const unsigned int maxlevel);
203 
208  void
209  cycle();
210 
221  void
223 
236  void
239 
246  void
248 
261  void
264 
268  unsigned int
269  get_maxlevel() const;
270 
274  unsigned int
275  get_minlevel() const;
276 
282  void
283  set_maxlevel(const unsigned int);
284 
300  void
301  set_minlevel(const unsigned int level, bool relative = false);
302 
307 
311  boost::signals2::connection
313  const std::function<void(const bool, const unsigned int)> &slot);
314 
318  boost::signals2::connection
320  const std::function<void(const bool, const unsigned int)> &slot);
321 
325  boost::signals2::connection
327  const std::function<void(const bool, const unsigned int)> &slot);
328 
332  boost::signals2::connection
334  const std::function<void(const bool, const unsigned int)> &slot);
335 
339  boost::signals2::connection
341  const std::function<void(const bool, const unsigned int)> &slot);
342 
346  boost::signals2::connection
348  const std::function<void(const bool, const unsigned int)> &slot);
349 
353  boost::signals2::connection
355  const std::function<void(const bool, const unsigned int)> &slot);
356 
357 private:
362 
369  void
370  level_v_step(const unsigned int level);
371 
379  void
380  level_step(const unsigned int level, Cycle cycle);
381 
386 
390  unsigned int minlevel;
391 
395  unsigned int maxlevel;
396 
397 public:
403 
408 
409 private:
414 
419 
420 
425 
431 
437 
443 
449 
456 
464 
471 
478 
479  template <int dim, class OtherVectorType, class TRANSFER>
480  friend class PreconditionMG;
481 };
482 
483 
498 template <int dim, typename VectorType, class TRANSFER>
500 {
501 public:
506  PreconditionMG(const DoFHandler<dim> &dof_handler,
508  const TRANSFER & transfer);
509 
514  PreconditionMG(const std::vector<const DoFHandler<dim> *> &dof_handler,
516  const TRANSFER & transfer);
517 
521  bool
522  empty() const;
523 
530  template <class OtherVectorType>
531  void
532  vmult(OtherVectorType &dst, const OtherVectorType &src) const;
533 
538  template <class OtherVectorType>
539  void
540  vmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
541 
547  template <class OtherVectorType>
548  void
549  Tvmult(OtherVectorType &dst, const OtherVectorType &src) const;
550 
556  template <class OtherVectorType>
557  void
558  Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
559 
566  IndexSet
567  locally_owned_range_indices(const unsigned int block = 0) const;
568 
575  IndexSet
576  locally_owned_domain_indices(const unsigned int block = 0) const;
577 
581  MPI_Comm
583 
587  boost::signals2::connection
588  connect_transfer_to_mg(const std::function<void(bool)> &slot);
589 
593  boost::signals2::connection
594  connect_transfer_to_global(const std::function<void(bool)> &slot);
595 
601 
605  const Multigrid<VectorType> &
606  get_multigrid() const;
607 
608 private:
612  std::vector<SmartPointer<const DoFHandler<dim>,
615 
620  std::vector<const DoFHandler<dim> *> dof_handler_vector_raw;
621 
627 
633 
639 
644 };
645 
648 #ifndef DOXYGEN
649 /* --------------------------- inline functions --------------------- */
650 
651 
652 template <typename VectorType>
654  const MGCoarseGridBase<VectorType> &coarse,
655  const MGTransferBase<VectorType> & transfer,
656  const MGSmootherBase<VectorType> & pre_smooth,
657  const MGSmootherBase<VectorType> &post_smooth,
658  const unsigned int min_level,
659  const unsigned int max_level,
660  Cycle cycle)
661  : cycle_type(cycle)
662  , matrix(&matrix, typeid(*this).name())
663  , coarse(&coarse, typeid(*this).name())
664  , transfer(&transfer, typeid(*this).name())
665  , pre_smooth(&pre_smooth, typeid(*this).name())
666  , post_smooth(&post_smooth, typeid(*this).name())
667  , edge_out(nullptr, typeid(*this).name())
668  , edge_in(nullptr, typeid(*this).name())
669  , edge_down(nullptr, typeid(*this).name())
670  , edge_up(nullptr, typeid(*this).name())
671 {
672  if (max_level == numbers::invalid_unsigned_int)
673  maxlevel = matrix.get_maxlevel();
674  else
675  maxlevel = max_level;
676  reinit(min_level, maxlevel);
677 }
678 
679 
680 
681 template <typename VectorType>
682 inline unsigned int
684 {
685  return maxlevel;
686 }
687 
688 
689 
690 template <typename VectorType>
691 inline unsigned int
693 {
694  return minlevel;
695 }
696 
697 
698 /* --------------------------- inline functions --------------------- */
699 
700 
701 namespace internal
702 {
703  namespace PreconditionMGImplementation
704  {
705  template <int dim,
706  typename VectorType,
707  class TRANSFER,
708  typename OtherVectorType>
709  typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
710  vmult(
711  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
712  ::Multigrid<VectorType> & multigrid,
713  const TRANSFER & transfer,
714  OtherVectorType & dst,
715  const OtherVectorType & src,
716  const bool uses_dof_handler_vector,
717  const typename ::mg::Signals &signals,
718  int)
719  {
720  signals.transfer_to_mg(true);
721  if (uses_dof_handler_vector)
722  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
723  else
724  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
725  signals.transfer_to_mg(false);
726 
727  multigrid.cycle();
728 
729  signals.transfer_to_global(true);
730  if (uses_dof_handler_vector)
731  transfer.copy_from_mg(dof_handler_vector, dst, multigrid.solution);
732  else
733  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
734  signals.transfer_to_global(false);
735  }
736 
737  template <int dim,
738  typename VectorType,
739  class TRANSFER,
740  typename OtherVectorType>
741  void
742  vmult(
743  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
744  ::Multigrid<VectorType> & multigrid,
745  const TRANSFER & transfer,
746  OtherVectorType & dst,
747  const OtherVectorType & src,
748  const bool uses_dof_handler_vector,
749  const typename ::mg::Signals &signals,
750  ...)
751  {
752  (void)uses_dof_handler_vector;
753  Assert(!uses_dof_handler_vector, ExcInternalError());
754 
755  signals.transfer_to_mg(true);
756  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
757  signals.transfer_to_mg(false);
758 
759  multigrid.cycle();
760 
761  signals.transfer_to_global(true);
762  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
763  signals.transfer_to_global(false);
764  }
765 
766  template <int dim,
767  typename VectorType,
768  class TRANSFER,
769  typename OtherVectorType>
770  typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
771  vmult_add(
772  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
773  ::Multigrid<VectorType> & multigrid,
774  const TRANSFER & transfer,
775  OtherVectorType & dst,
776  const OtherVectorType & src,
777  const bool uses_dof_handler_vector,
778  const typename ::mg::Signals &signals,
779  int)
780  {
781  signals.transfer_to_mg(true);
782  if (uses_dof_handler_vector)
783  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
784  else
785  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
786  signals.transfer_to_mg(false);
787 
788  multigrid.cycle();
789 
790  signals.transfer_to_global(true);
791  if (uses_dof_handler_vector)
792  transfer.copy_from_mg_add(dof_handler_vector, dst, multigrid.solution);
793  else
794  transfer.copy_from_mg_add(*dof_handler_vector[0],
795  dst,
796  multigrid.solution);
797  signals.transfer_to_global(false);
798  }
799 
800  template <int dim,
801  typename VectorType,
802  class TRANSFER,
803  typename OtherVectorType>
804  void
805  vmult_add(
806  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
807  ::Multigrid<VectorType> & multigrid,
808  const TRANSFER & transfer,
809  OtherVectorType & dst,
810  const OtherVectorType & src,
811  const bool uses_dof_handler_vector,
812  const typename ::mg::Signals &signals,
813  ...)
814  {
815  (void)uses_dof_handler_vector;
816  Assert(!uses_dof_handler_vector, ExcInternalError());
817 
818  signals.transfer_to_mg(true);
819  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
820  signals.transfer_to_mg(false);
821 
822  multigrid.cycle();
823 
824  signals.transfer_to_global(true);
825  transfer.copy_from_mg_add(*dof_handler_vector[0],
826  dst,
827  multigrid.solution);
828  signals.transfer_to_global(false);
829  }
830  } // namespace PreconditionMGImplementation
831 } // namespace internal
832 
833 template <int dim, typename VectorType, class TRANSFER>
835  const DoFHandler<dim> &dof_handler,
837  const TRANSFER & transfer)
838  : dof_handler_vector(1, &dof_handler)
839  , dof_handler_vector_raw(1, &dof_handler)
840  , multigrid(&mg)
841  , transfer(&transfer)
842  , uses_dof_handler_vector(false)
843 {}
844 
845 template <int dim, typename VectorType, class TRANSFER>
847  const std::vector<const DoFHandler<dim> *> &dof_handler,
849  const TRANSFER & transfer)
850  : dof_handler_vector(dof_handler.size())
851  , dof_handler_vector_raw(dof_handler.size())
852  , multigrid(&mg)
853  , transfer(&transfer)
854  , uses_dof_handler_vector(true)
855 {
856  for (unsigned int i = 0; i < dof_handler.size(); ++i)
857  {
858  dof_handler_vector[i] = dof_handler[i];
859  dof_handler_vector_raw[i] = dof_handler[i];
860  }
861 }
862 
863 template <int dim, typename VectorType, class TRANSFER>
864 inline bool
866 {
867  return false;
868 }
869 
870 template <int dim, typename VectorType, class TRANSFER>
871 template <class OtherVectorType>
872 void
874  OtherVectorType & dst,
875  const OtherVectorType &src) const
876 {
877  internal::PreconditionMGImplementation::vmult(dof_handler_vector_raw,
878  *multigrid,
879  *transfer,
880  dst,
881  src,
882  uses_dof_handler_vector,
883  this->signals,
884  0);
885 }
886 
887 
888 template <int dim, typename VectorType, class TRANSFER>
889 IndexSet
891  const unsigned int block) const
892 {
893  AssertIndexRange(block, dof_handler_vector.size());
894  return dof_handler_vector[block]->locally_owned_dofs();
895 }
896 
897 
898 template <int dim, typename VectorType, class TRANSFER>
899 IndexSet
901  const unsigned int block) const
902 {
903  AssertIndexRange(block, dof_handler_vector.size());
904  return dof_handler_vector[block]->locally_owned_dofs();
905 }
906 
907 
908 
909 template <int dim, typename VectorType, class TRANSFER>
910 MPI_Comm
912 {
913  // currently parallel GMG works with parallel triangulations only,
914  // so it should be a safe bet to use it to query MPI communicator:
915  const Triangulation<dim> &tria = dof_handler_vector[0]->get_triangulation();
916  const parallel::TriangulationBase<dim> *ptria =
917  dynamic_cast<const parallel::TriangulationBase<dim> *>(&tria);
918  Assert(ptria != nullptr, ExcInternalError());
919  return ptria->get_communicator();
920 }
921 
922 
923 
924 template <int dim, typename VectorType, class TRANSFER>
925 boost::signals2::connection
927  const std::function<void(bool)> &slot)
928 {
929  return this->signals.transfer_to_mg.connect(slot);
930 }
931 
932 
933 
934 template <int dim, typename VectorType, class TRANSFER>
935 boost::signals2::connection
937  const std::function<void(bool)> &slot)
938 {
939  return this->signals.transfer_to_global.connect(slot);
940 }
941 
942 
943 
944 template <int dim, typename VectorType, class TRANSFER>
945 template <class OtherVectorType>
946 void
948  OtherVectorType & dst,
949  const OtherVectorType &src) const
950 {
951  internal::PreconditionMGImplementation::vmult_add(dof_handler_vector_raw,
952  *multigrid,
953  *transfer,
954  dst,
955  src,
956  uses_dof_handler_vector,
957  this->signals,
958  0);
959 }
960 
961 
962 template <int dim, typename VectorType, class TRANSFER>
963 template <class OtherVectorType>
964 void
966  const OtherVectorType &) const
967 {
968  Assert(false, ExcNotImplemented());
969 }
970 
971 
972 template <int dim, typename VectorType, class TRANSFER>
973 template <class OtherVectorType>
974 void
976  OtherVectorType &,
977  const OtherVectorType &) const
978 {
979  Assert(false, ExcNotImplemented());
980 }
981 
982 
983 template <int dim, typename VectorType, class TRANSFER>
986 {
987  return *this->multigrid;
988 }
989 
990 
991 template <int dim, typename VectorType, class TRANSFER>
992 const Multigrid<VectorType> &
994 {
995  return *this->multigrid;
996 }
997 
998 #endif // DOXYGEN
999 
1001 
1002 #endif
@ v_cycle
The V-cycle.
Definition: multigrid.h:170
@ w_cycle
The W-cycle.
Definition: multigrid.h:172
@ f_cycle
The F-cycle.
Definition: multigrid.h:174
MGLevelObject< VectorType > defect2
Definition: multigrid.h:418
void cycle()
void set_edge_flux_matrices(const MGMatrixBase< VectorType > &edge_down, const MGMatrixBase< VectorType > &edge_up)
SmartPointer< const MGTransferBase< VectorType >, Multigrid< VectorType > > transfer
Definition: multigrid.h:436
unsigned int minlevel
Definition: multigrid.h:390
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_up
Definition: multigrid.h:477
mg::Signals signals
Definition: multigrid.h:361
SmartPointer< const MGMatrixBase< VectorType > > edge_in
Definition: multigrid.h:463
SmartPointer< const MGMatrixBase< VectorType > > edge_out
Definition: multigrid.h:455
MGLevelObject< VectorType > solution
Definition: multigrid.h:407
boost::signals2::connection connect_edge_prolongation(const std::function< void(const bool, const unsigned int)> &slot)
unsigned int maxlevel
Definition: multigrid.h:395
VectorType vector_type
Definition: multigrid.h:177
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > matrix
Definition: multigrid.h:424
void set_maxlevel(const unsigned int)
unsigned int get_minlevel() const
void set_edge_matrices(const MGMatrixBase< VectorType > &edge_out, const MGMatrixBase< VectorType > &edge_in)
boost::signals2::connection connect_residual_step(const std::function< void(const bool, const unsigned int)> &slot)
void level_step(const unsigned int level, Cycle cycle)
void set_cycle(Cycle)
unsigned int get_maxlevel() const
const VectorType const_vector_type
Definition: multigrid.h:178
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > pre_smooth
Definition: multigrid.h:442
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > post_smooth
Definition: multigrid.h:448
boost::signals2::connection connect_post_smoother_step(const std::function< void(const bool, const unsigned int)> &slot)
boost::signals2::connection connect_coarse_solve(const std::function< void(const bool, const unsigned int)> &slot)
SmartPointer< const MGCoarseGridBase< VectorType >, Multigrid< VectorType > > coarse
Definition: multigrid.h:430
Cycle cycle_type
Definition: multigrid.h:385
boost::signals2::connection connect_prolongation(const std::function< void(const bool, const unsigned int)> &slot)
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_down
Definition: multigrid.h:470
void set_edge_in_matrix(const MGMatrixBase< VectorType > &edge_in)
MGLevelObject< VectorType > t
Definition: multigrid.h:413
void set_minlevel(const unsigned int level, bool relative=false)
void level_v_step(const unsigned int level)
Multigrid(const MGMatrixBase< VectorType > &matrix, const MGCoarseGridBase< VectorType > &coarse, const MGTransferBase< VectorType > &transfer, const MGSmootherBase< VectorType > &pre_smooth, const MGSmootherBase< VectorType > &post_smooth, const unsigned int minlevel=0, const unsigned int maxlevel=numbers::invalid_unsigned_int, Cycle cycle=v_cycle)
MGLevelObject< VectorType > defect
Definition: multigrid.h:402
void vcycle()
boost::signals2::connection connect_restriction(const std::function< void(const bool, const unsigned int)> &slot)
void reinit(const unsigned int minlevel, const unsigned int maxlevel)
boost::signals2::connection connect_pre_smoother_step(const std::function< void(const bool, const unsigned int)> &slot)
SmartPointer< Multigrid< VectorType >, PreconditionMG< dim, VectorType, TRANSFER > > multigrid
Definition: multigrid.h:626
void vmult_add(OtherVectorType &dst, const OtherVectorType &src) const
void vmult(OtherVectorType &dst, const OtherVectorType &src) const
IndexSet locally_owned_range_indices(const unsigned int block=0) const
SmartPointer< const TRANSFER, PreconditionMG< dim, VectorType, TRANSFER > > transfer
Definition: multigrid.h:632
PreconditionMG(const std::vector< const DoFHandler< dim > * > &dof_handler, Multigrid< VectorType > &mg, const TRANSFER &transfer)
std::vector< SmartPointer< const DoFHandler< dim >, PreconditionMG< dim, VectorType, TRANSFER > > > dof_handler_vector
Definition: multigrid.h:614
bool empty() const
void Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const
void Tvmult(OtherVectorType &dst, const OtherVectorType &src) const
boost::signals2::connection connect_transfer_to_mg(const std::function< void(bool)> &slot)
const Multigrid< VectorType > & get_multigrid() const
PreconditionMG(const DoFHandler< dim > &dof_handler, Multigrid< VectorType > &mg, const TRANSFER &transfer)
mg::Signals signals
Definition: multigrid.h:643
IndexSet locally_owned_domain_indices(const unsigned int block=0) const
boost::signals2::connection connect_transfer_to_global(const std::function< void(bool)> &slot)
const bool uses_dof_handler_vector
Definition: multigrid.h:638
MPI_Comm get_mpi_communicator() const
Multigrid< VectorType > & get_multigrid()
std::vector< const DoFHandler< dim > * > dof_handler_vector_raw
Definition: multigrid.h:620
Triangulation< dim, spacedim > & get_triangulation()
virtual MPI_Comm get_communicator() const override
Definition: tria_base.cc:143
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:416
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:417
unsigned int level
Definition: grid_out.cc:4606
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcNotImplemented()
#define AssertIndexRange(index, range)
Definition: exceptions.h:1732
@ matrix
Contents is actually a matrix.
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
Definition: mg.h:82
static const unsigned int invalid_unsigned_int
Definition: types.h:201
boost::signals2::signal< void(const bool before, const unsigned int level)> prolongation
Definition: multigrid.h:105
boost::signals2::signal< void(const bool before, const unsigned int level)> pre_smoother_step
Definition: multigrid.h:116
boost::signals2::signal< void(const bool before, const unsigned int level)> coarse_solve
Definition: multigrid.h:86
boost::signals2::signal< void(const bool before, const unsigned int level)> edge_prolongation
Definition: multigrid.h:139
boost::signals2::signal< void(const bool before, const unsigned int level)> post_smoother_step
Definition: multigrid.h:124
boost::signals2::signal< void(const bool before)> transfer_to_global
Definition: multigrid.h:74
boost::signals2::signal< void(const bool before, const unsigned int level)> residual_step
Definition: multigrid.h:132
boost::signals2::signal< void(const bool before, const unsigned int level)> restriction
Definition: multigrid.h:97
boost::signals2::signal< void(const bool before)> transfer_to_mg
Definition: multigrid.h:67
const ::Triangulation< dim, spacedim > & tria