Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
multigrid.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 2018 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 
22 #include <deal.II/base/mg_level_object.h>
23 #include <deal.II/base/smartpointer.h>
24 #include <deal.II/base/subscriptor.h>
25 
26 #include <deal.II/distributed/tria.h>
27 
28 #include <deal.II/dofs/dof_handler.h>
29 
30 #include <deal.II/lac/sparse_matrix.h>
31 #include <deal.II/lac/vector.h>
32 
33 #include <deal.II/multigrid/mg_base.h>
34 
35 #include <vector>
36 
37 DEAL_II_NAMESPACE_OPEN
38 
41 
42 namespace mg
43 {
53  struct Signals
54  {
60  boost::signals2::signal<void(const bool before)> transfer_to_mg;
61 
67  boost::signals2::signal<void(const bool before)> transfer_to_global;
68 
74  boost::signals2::signal<void(const bool before, const unsigned int level)>
76 
82  boost::signals2::signal<void(const bool before, const unsigned int level)>
84 
90  boost::signals2::signal<void(const bool before, const unsigned int level)>
92 
98  boost::signals2::signal<void(const bool before, const unsigned int level)>
100 
107  boost::signals2::signal<void(const bool before, const unsigned int level)>
109  };
110 } // namespace mg
111 
136 template <typename VectorType>
137 class Multigrid : public Subscriptor
138 {
139 public:
143  enum Cycle
144  {
151  };
152 
153  using vector_type = VectorType;
154  using const_vector_type = const VectorType;
155 
171  template <int dim>
172 #ifndef _MSC_VER
173  DEAL_II_DEPRECATED
174 #endif
175  Multigrid(const DoFHandler<dim> & mg_dof_handler,
181  const unsigned int minlevel = 0,
182  const unsigned int maxlevel = numbers::invalid_unsigned_int,
183  Cycle cycle = v_cycle);
184 
199  const unsigned int minlevel = 0,
200  const unsigned int maxlevel = numbers::invalid_unsigned_int,
201  Cycle cycle = v_cycle);
202 
206  void
207  reinit(const unsigned int minlevel, const unsigned int maxlevel);
208 
213  void
214  cycle();
215 
226  void
227  vcycle();
228 
241  void
244 
257  void
260 
264  unsigned int
265  get_maxlevel() const;
266 
270  unsigned int
271  get_minlevel() const;
272 
278  void
279  set_maxlevel(const unsigned int);
280 
296  void
297  set_minlevel(const unsigned int level, bool relative = false);
298 
302  void set_cycle(Cycle);
303 
310  DEAL_II_DEPRECATED
311  void
312  set_debug(const unsigned int);
313 
317  boost::signals2::connection
319  const std::function<void(const bool, const unsigned int)> &slot);
320 
324  boost::signals2::connection
326  const std::function<void(const bool, const unsigned int)> &slot);
327 
331  boost::signals2::connection
333  const std::function<void(const bool, const unsigned int)> &slot);
334 
338  boost::signals2::connection
340  const std::function<void(const bool, const unsigned int)> &slot);
341 
345  boost::signals2::connection
347  const std::function<void(const bool, const unsigned int)> &slot);
348 
349 private:
354 
361  void
362  level_v_step(const unsigned int level);
363 
371  void
372  level_step(const unsigned int level, Cycle cycle);
373 
378 
382  unsigned int minlevel;
383 
387  unsigned int maxlevel;
388 
389 public:
395 
400 
401 private:
406 
411 
412 
417 
423 
429 
435 
441 
448 
456 
463 
470 
474  unsigned int debug;
475 
476  template <int dim, class OtherVectorType, class TRANSFER>
477  friend class PreconditionMG;
478 };
479 
480 
497 template <int dim, typename VectorType, class TRANSFER>
499 {
500 public:
505  PreconditionMG(const DoFHandler<dim> &dof_handler,
507  const TRANSFER & transfer);
508 
513  PreconditionMG(const std::vector<const DoFHandler<dim> *> &dof_handler,
515  const TRANSFER & transfer);
516 
520  bool
521  empty() const;
522 
529  template <class OtherVectorType>
530  void
531  vmult(OtherVectorType &dst, const OtherVectorType &src) const;
532 
537  template <class OtherVectorType>
538  void
539  vmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
540 
546  template <class OtherVectorType>
547  void
548  Tvmult(OtherVectorType &dst, const OtherVectorType &src) const;
549 
555  template <class OtherVectorType>
556  void
557  Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
558 
565  IndexSet
566  locally_owned_range_indices(const unsigned int block = 0) const;
567 
574  IndexSet
575  locally_owned_domain_indices(const unsigned int block = 0) const;
576 
580  MPI_Comm
581  get_mpi_communicator() const;
582 
586  boost::signals2::connection
587  connect_transfer_to_mg(const std::function<void(bool)> &slot);
588 
592  boost::signals2::connection
593  connect_transfer_to_global(const std::function<void(bool)> &slot);
594 
595 private:
599  std::vector<SmartPointer<const DoFHandler<dim>,
602 
607  std::vector<const DoFHandler<dim> *> dof_handler_vector_raw;
608 
614 
620 
626 
631 };
632 
635 #ifndef DOXYGEN
636 /* --------------------------- inline functions --------------------- */
637 
638 
639 template <typename VectorType>
640 template <int dim>
642  const MGMatrixBase<VectorType> &matrix,
643  const MGCoarseGridBase<VectorType> &coarse,
644  const MGTransferBase<VectorType> & transfer,
645  const MGSmootherBase<VectorType> & pre_smooth,
646  const MGSmootherBase<VectorType> &post_smooth,
647  const unsigned int min_level,
648  const unsigned int max_level,
649  Cycle cycle)
650  : cycle_type(cycle)
651  , minlevel(min_level)
652  , matrix(&matrix, typeid(*this).name())
653  , coarse(&coarse, typeid(*this).name())
654  , transfer(&transfer, typeid(*this).name())
655  , pre_smooth(&pre_smooth, typeid(*this).name())
656  , post_smooth(&post_smooth, typeid(*this).name())
657  , edge_down(nullptr, typeid(*this).name())
658  , edge_up(nullptr, typeid(*this).name())
659  , debug(0)
660 {
661  const unsigned int dof_handler_max_level =
662  mg_dof_handler.get_triangulation().n_global_levels() - 1;
663  if (max_level == numbers::invalid_unsigned_int)
664  maxlevel = dof_handler_max_level;
665  else
666  maxlevel = max_level;
667 
668  reinit(minlevel, maxlevel);
669 }
670 
671 
672 
673 template <typename VectorType>
675  const MGCoarseGridBase<VectorType> &coarse,
676  const MGTransferBase<VectorType> & transfer,
677  const MGSmootherBase<VectorType> & pre_smooth,
678  const MGSmootherBase<VectorType> &post_smooth,
679  const unsigned int min_level,
680  const unsigned int max_level,
681  Cycle cycle)
682  : cycle_type(cycle)
683  , matrix(&matrix, typeid(*this).name())
684  , coarse(&coarse, typeid(*this).name())
685  , transfer(&transfer, typeid(*this).name())
686  , pre_smooth(&pre_smooth, typeid(*this).name())
687  , post_smooth(&post_smooth, typeid(*this).name())
688  , edge_out(nullptr, typeid(*this).name())
689  , edge_in(nullptr, typeid(*this).name())
690  , edge_down(nullptr, typeid(*this).name())
691  , edge_up(nullptr, typeid(*this).name())
692  , debug(0)
693 {
694  if (max_level == numbers::invalid_unsigned_int)
695  maxlevel = matrix.get_maxlevel();
696  else
697  maxlevel = max_level;
698  reinit(min_level, maxlevel);
699 }
700 
701 
702 
703 template <typename VectorType>
704 inline unsigned int
706 {
707  return maxlevel;
708 }
709 
710 
711 
712 template <typename VectorType>
713 inline unsigned int
715 {
716  return minlevel;
717 }
718 
719 
720 /* --------------------------- inline functions --------------------- */
721 
722 
723 namespace internal
724 {
725  namespace PreconditionMGImplementation
726  {
727  template <int dim,
728  typename VectorType,
729  class TRANSFER,
730  typename OtherVectorType>
731  typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
732  vmult(
733  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
734  ::Multigrid<VectorType> & multigrid,
735  const TRANSFER & transfer,
736  OtherVectorType & dst,
737  const OtherVectorType & src,
738  const bool uses_dof_handler_vector,
739  const typename ::mg::Signals &signals,
740  int)
741  {
742  signals.transfer_to_mg(true);
743  if (uses_dof_handler_vector)
744  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
745  else
746  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
747  signals.transfer_to_mg(false);
748 
749  multigrid.cycle();
750 
751  signals.transfer_to_global(true);
752  if (uses_dof_handler_vector)
753  transfer.copy_from_mg(dof_handler_vector, dst, multigrid.solution);
754  else
755  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
756  signals.transfer_to_global(false);
757  }
758 
759  template <int dim,
760  typename VectorType,
761  class TRANSFER,
762  typename OtherVectorType>
763  void
764  vmult(
765  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
766  ::Multigrid<VectorType> & multigrid,
767  const TRANSFER & transfer,
768  OtherVectorType & dst,
769  const OtherVectorType & src,
770  const bool uses_dof_handler_vector,
771  const typename ::mg::Signals &signals,
772  ...)
773  {
774  (void)uses_dof_handler_vector;
775  Assert(!uses_dof_handler_vector, ExcInternalError());
776 
777  signals.transfer_to_mg(true);
778  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
779  signals.transfer_to_mg(false);
780 
781  multigrid.cycle();
782 
783  signals.transfer_to_global(true);
784  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
785  signals.transfer_to_global(false);
786  }
787 
788  template <int dim,
789  typename VectorType,
790  class TRANSFER,
791  typename OtherVectorType>
792  typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
793  vmult_add(
794  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
795  ::Multigrid<VectorType> & multigrid,
796  const TRANSFER & transfer,
797  OtherVectorType & dst,
798  const OtherVectorType & src,
799  const bool uses_dof_handler_vector,
800  const typename ::mg::Signals &signals,
801  int)
802  {
803  signals.transfer_to_mg(true);
804  if (uses_dof_handler_vector)
805  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
806  else
807  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
808  signals.transfer_to_mg(false);
809 
810  multigrid.cycle();
811 
812  signals.transfer_to_global(true);
813  if (uses_dof_handler_vector)
814  transfer.copy_from_mg_add(dof_handler_vector, dst, multigrid.solution);
815  else
816  transfer.copy_from_mg_add(*dof_handler_vector[0],
817  dst,
818  multigrid.solution);
819  signals.transfer_to_global(false);
820  }
821 
822  template <int dim,
823  typename VectorType,
824  class TRANSFER,
825  typename OtherVectorType>
826  void
827  vmult_add(
828  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
829  ::Multigrid<VectorType> & multigrid,
830  const TRANSFER & transfer,
831  OtherVectorType & dst,
832  const OtherVectorType & src,
833  const bool uses_dof_handler_vector,
834  const typename ::mg::Signals &signals,
835  ...)
836  {
837  (void)uses_dof_handler_vector;
838  Assert(!uses_dof_handler_vector, ExcInternalError());
839 
840  signals.transfer_to_mg(true);
841  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
842  signals.transfer_to_mg(false);
843 
844  multigrid.cycle();
845 
846  signals.transfer_to_global(true);
847  transfer.copy_from_mg_add(*dof_handler_vector[0],
848  dst,
849  multigrid.solution);
850  signals.transfer_to_global(false);
851  }
852  } // namespace PreconditionMGImplementation
853 } // namespace internal
854 
855 template <int dim, typename VectorType, class TRANSFER>
857  const DoFHandler<dim> &dof_handler,
859  const TRANSFER & transfer)
860  : dof_handler_vector(1, &dof_handler)
861  , dof_handler_vector_raw(1, &dof_handler)
862  , multigrid(&mg)
863  , transfer(&transfer)
864  , uses_dof_handler_vector(false)
865 {}
866 
867 template <int dim, typename VectorType, class TRANSFER>
869  const std::vector<const DoFHandler<dim> *> &dof_handler,
871  const TRANSFER & transfer)
872  : dof_handler_vector(dof_handler.size())
873  , dof_handler_vector_raw(dof_handler.size())
874  , multigrid(&mg)
875  , transfer(&transfer)
876  , uses_dof_handler_vector(true)
877 {
878  for (unsigned int i = 0; i < dof_handler.size(); ++i)
879  {
880  dof_handler_vector[i] = dof_handler[i];
881  dof_handler_vector_raw[i] = dof_handler[i];
882  }
883 }
884 
885 template <int dim, typename VectorType, class TRANSFER>
886 inline bool
888 {
889  return false;
890 }
891 
892 template <int dim, typename VectorType, class TRANSFER>
893 template <class OtherVectorType>
894 void
896  OtherVectorType & dst,
897  const OtherVectorType &src) const
898 {
899  internal::PreconditionMGImplementation::vmult(dof_handler_vector_raw,
900  *multigrid,
901  *transfer,
902  dst,
903  src,
904  uses_dof_handler_vector,
905  this->signals,
906  0);
907 }
908 
909 
910 template <int dim, typename VectorType, class TRANSFER>
911 IndexSet
913  const unsigned int block) const
914 {
915  AssertIndexRange(block, dof_handler_vector.size());
916  return dof_handler_vector[block]->locally_owned_dofs();
917 }
918 
919 
920 template <int dim, typename VectorType, class TRANSFER>
921 IndexSet
923  const unsigned int block) const
924 {
925  AssertIndexRange(block, dof_handler_vector.size());
926  return dof_handler_vector[block]->locally_owned_dofs();
927 }
928 
929 
930 
931 template <int dim, typename VectorType, class TRANSFER>
932 MPI_Comm
934 {
935  // currently parallel GMG works with distributed Triangulation only,
936  // so it should be a safe bet to use it to query MPI communicator:
937  const Triangulation<dim> &tria = dof_handler_vector[0]->get_triangulation();
939  dynamic_cast<const parallel::distributed::Triangulation<dim> *>(&tria);
940  Assert(ptria != nullptr, ExcInternalError());
941  return ptria->get_communicator();
942 }
943 
944 
945 
946 template <int dim, typename VectorType, class TRANSFER>
947 boost::signals2::connection
949  const std::function<void(bool)> &slot)
950 {
951  return this->signals.transfer_to_mg.connect(slot);
952 }
953 
954 
955 
956 template <int dim, typename VectorType, class TRANSFER>
957 boost::signals2::connection
959  const std::function<void(bool)> &slot)
960 {
961  return this->signals.transfer_to_global.connect(slot);
962 }
963 
964 
965 
966 template <int dim, typename VectorType, class TRANSFER>
967 template <class OtherVectorType>
968 void
970  OtherVectorType & dst,
971  const OtherVectorType &src) const
972 {
973  internal::PreconditionMGImplementation::vmult_add(dof_handler_vector_raw,
974  *multigrid,
975  *transfer,
976  dst,
977  src,
978  uses_dof_handler_vector,
979  this->signals,
980  0);
981 }
982 
983 
984 template <int dim, typename VectorType, class TRANSFER>
985 template <class OtherVectorType>
986 void
988  const OtherVectorType &) const
989 {
990  Assert(false, ExcNotImplemented());
991 }
992 
993 
994 template <int dim, typename VectorType, class TRANSFER>
995 template <class OtherVectorType>
996 void
998  OtherVectorType &,
999  const OtherVectorType &) const
1000 {
1001  Assert(false, ExcNotImplemented());
1002 }
1003 
1004 #endif // DOXYGEN
1005 
1006 DEAL_II_NAMESPACE_CLOSE
1007 
1008 #endif
boost::signals2::signal< void(const bool before, const unsigned int level)> coarse_solve
Definition: multigrid.h:75
boost::signals2::connection connect_transfer_to_global(const std::function< void(bool)> &slot)
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > matrix
Definition: multigrid.h:416
void vcycle()
const Triangulation< dim, spacedim > & get_triangulation() const
void cycle()
static const unsigned int invalid_unsigned_int
Definition: types.h:173
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_down
Definition: multigrid.h:462
SmartPointer< const MGMatrixBase< VectorType > > edge_in
Definition: multigrid.h:455
void vmult(OtherVectorType &dst, const OtherVectorType &src) const
MPI_Comm get_mpi_communicator() const
Contents is actually a matrix.
void set_edge_matrices(const MGMatrixBase< VectorType > &edge_out, const MGMatrixBase< VectorType > &edge_in)
Definition: mg.h:81
unsigned int get_maxlevel() const
boost::signals2::signal< void(const bool before, const unsigned int level)> prolongation
Definition: multigrid.h:91
#define AssertIndexRange(index, range)
Definition: exceptions.h:1637
SmartPointer< const MGCoarseGridBase< VectorType >, Multigrid< VectorType > > coarse
Definition: multigrid.h:422
boost::signals2::signal< void(const bool before, const unsigned int level)> post_smoother_step
Definition: multigrid.h:108
boost::signals2::signal< void(const bool before)> transfer_to_global
Definition: multigrid.h:67
void set_maxlevel(const unsigned int)
boost::signals2::signal< void(const bool before, const unsigned int level)> restriction
Definition: multigrid.h:83
MGLevelObject< VectorType > t
Definition: multigrid.h:405
boost::signals2::connection connect_coarse_solve(const std::function< void(const bool, const unsigned int)> &slot)
bool empty() const
const bool uses_dof_handler_vector
Definition: multigrid.h:625
unsigned int debug
Definition: multigrid.h:474
void set_debug(const unsigned int)
void level_v_step(const unsigned int level)
mg::Signals signals
Definition: multigrid.h:353
unsigned int maxlevel
Definition: multigrid.h:387
IndexSet locally_owned_range_indices(const unsigned int block=0) const
#define Assert(cond, exc)
Definition: exceptions.h:1407
The F-cycle.
Definition: multigrid.h:150
virtual MPI_Comm get_communicator() const
Definition: tria_base.cc:155
void Tvmult(OtherVectorType &dst, const OtherVectorType &src) const
boost::signals2::connection connect_restriction(const std::function< void(const bool, const unsigned int)> &slot)
unsigned int get_minlevel() const
SmartPointer< Multigrid< VectorType >, PreconditionMG< dim, VectorType, TRANSFER > > multigrid
Definition: multigrid.h:613
boost::signals2::signal< void(const bool before)> transfer_to_mg
Definition: multigrid.h:60
SmartPointer< const MGMatrixBase< VectorType > > edge_out
Definition: multigrid.h:447
void vmult_add(OtherVectorType &dst, const OtherVectorType &src) const
The W-cycle.
Definition: multigrid.h:148
boost::signals2::signal< void(const bool before, const unsigned int level)> pre_smoother_step
Definition: multigrid.h:99
void set_minlevel(const unsigned int level, bool relative=false)
The V-cycle.
Definition: multigrid.h:146
void set_edge_flux_matrices(const MGMatrixBase< VectorType > &edge_down, const MGMatrixBase< VectorType > &edge_up)
boost::signals2::connection connect_transfer_to_mg(const std::function< void(bool)> &slot)
void Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const
unsigned int minlevel
Definition: multigrid.h:382
SmartPointer< const TRANSFER, PreconditionMG< dim, VectorType, TRANSFER > > transfer
Definition: multigrid.h:619
MGLevelObject< VectorType > solution
Definition: multigrid.h:399
void set_cycle(Cycle)
std::vector< SmartPointer< const DoFHandler< dim >, PreconditionMG< dim, VectorType, TRANSFER > > > dof_handler_vector
Definition: multigrid.h:601
boost::signals2::connection connect_pre_smoother_step(const std::function< void(const bool, const unsigned int)> &slot)
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > post_smooth
Definition: multigrid.h:440
boost::signals2::connection connect_prolongation(const std::function< void(const bool, const unsigned int)> &slot)
IndexSet locally_owned_domain_indices(const unsigned int block=0) const
void level_step(const unsigned int level, Cycle cycle)
std::vector< const DoFHandler< dim > * > dof_handler_vector_raw
Definition: multigrid.h:607
virtual unsigned int n_global_levels() const
boost::signals2::connection connect_post_smoother_step(const std::function< void(const bool, const unsigned int)> &slot)
mg::Signals signals
Definition: multigrid.h:630
static ::ExceptionBase & ExcNotImplemented()
Multigrid(const DoFHandler< dim > &mg_dof_handler, 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 > defect2
Definition: multigrid.h:410
MGLevelObject< VectorType > defect
Definition: multigrid.h:394
PreconditionMG(const DoFHandler< dim > &dof_handler, Multigrid< VectorType > &mg, const TRANSFER &transfer)
SmartPointer< const MGTransferBase< VectorType >, Multigrid< VectorType > > transfer
Definition: multigrid.h:428
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_up
Definition: multigrid.h:469
Cycle cycle_type
Definition: multigrid.h:377
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > pre_smooth
Definition: multigrid.h:434
static ::ExceptionBase & ExcInternalError()
Triangulation< dim, spacedim > & get_triangulation()
Definition: tria.cc:13220
void reinit(const unsigned int minlevel, const unsigned int maxlevel)