Reference documentation for deal.II version Git 73c87d96ef 2021-11-30 22:54:44 +0100
\(\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 
165 template <typename VectorType>
166 class Multigrid : public Subscriptor
167 {
168 public:
172  enum Cycle
173  {
179  f_cycle
180  };
181 
184 
195  const MGCoarseGridBase<VectorType> &coarse,
196  const MGTransferBase<VectorType> & transfer,
197  const MGSmootherBase<VectorType> & pre_smooth,
198  const MGSmootherBase<VectorType> & post_smooth,
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
242  set_edge_matrices(const MGMatrixBase<VectorType> &edge_out,
243  const MGMatrixBase<VectorType> &edge_in);
244 
257  void
258  set_edge_flux_matrices(const MGMatrixBase<VectorType> &edge_down,
259  const MGMatrixBase<VectorType> &edge_up);
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 
307  boost::signals2::connection
308  connect_pre_smoother_step(
309  const std::function<void(const bool, const unsigned int)> &slot);
310 
314  boost::signals2::connection
315  connect_residual_step(
316  const std::function<void(const bool, const unsigned int)> &slot);
317 
321  boost::signals2::connection
322  connect_restriction(
323  const std::function<void(const bool, const unsigned int)> &slot);
324 
328  boost::signals2::connection
329  connect_coarse_solve(
330  const std::function<void(const bool, const unsigned int)> &slot);
331 
335  boost::signals2::connection
336  connect_prolongation(
337  const std::function<void(const bool, const unsigned int)> &slot);
338 
342  boost::signals2::connection
343  connect_edge_prolongation(
344  const std::function<void(const bool, const unsigned int)> &slot);
345 
349  boost::signals2::connection
350  connect_post_smoother_step(
351  const std::function<void(const bool, const unsigned int)> &slot);
352 
353 private:
358 
365  void
366  level_v_step(const unsigned int level);
367 
375  void
376  level_step(const unsigned int level, Cycle cycle);
377 
382 
386  unsigned int minlevel;
387 
391  unsigned int maxlevel;
392 
393 public:
399 
404 
405 private:
410 
415 
416 
421 
427 
433 
439 
445 
452 
460 
467 
474 
475  template <int dim, class OtherVectorType, class TRANSFER>
476  friend class PreconditionMG;
477 };
478 
479 
494 template <int dim, typename VectorType, class TRANSFER>
496 {
497 public:
502  PreconditionMG(const DoFHandler<dim> &dof_handler,
504  const TRANSFER & transfer);
505 
510  PreconditionMG(const std::vector<const DoFHandler<dim> *> &dof_handler,
512  const TRANSFER & transfer);
513 
517  bool
518  empty() const;
519 
526  template <class OtherVectorType>
527  void
528  vmult(OtherVectorType &dst, const OtherVectorType &src) const;
529 
534  template <class OtherVectorType>
535  void
536  vmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
537 
543  template <class OtherVectorType>
544  void
545  Tvmult(OtherVectorType &dst, const OtherVectorType &src) const;
546 
552  template <class OtherVectorType>
553  void
554  Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
555 
562  IndexSet
563  locally_owned_range_indices(const unsigned int block = 0) const;
564 
571  IndexSet
572  locally_owned_domain_indices(const unsigned int block = 0) const;
573 
577  MPI_Comm
578  get_mpi_communicator() const;
579 
583  boost::signals2::connection
584  connect_transfer_to_mg(const std::function<void(bool)> &slot);
585 
589  boost::signals2::connection
590  connect_transfer_to_global(const std::function<void(bool)> &slot);
591 
596  get_multigrid();
597 
601  const Multigrid<VectorType> &
602  get_multigrid() const;
603 
604 private:
608  std::vector<SmartPointer<const DoFHandler<dim>,
611 
616  std::vector<const DoFHandler<dim> *> dof_handler_vector_raw;
617 
623 
629 
635 
640 };
641 
644 #ifndef DOXYGEN
645 /* --------------------------- inline functions --------------------- */
646 
647 
648 template <typename VectorType>
650  const MGCoarseGridBase<VectorType> &coarse,
651  const MGTransferBase<VectorType> & transfer,
652  const MGSmootherBase<VectorType> & pre_smooth,
653  const MGSmootherBase<VectorType> &post_smooth,
654  const unsigned int min_level,
655  const unsigned int max_level,
656  Cycle cycle)
657  : cycle_type(cycle)
658  , matrix(&matrix, typeid(*this).name())
659  , coarse(&coarse, typeid(*this).name())
660  , transfer(&transfer, typeid(*this).name())
661  , pre_smooth(&pre_smooth, typeid(*this).name())
662  , post_smooth(&post_smooth, typeid(*this).name())
663  , edge_out(nullptr, typeid(*this).name())
664  , edge_in(nullptr, typeid(*this).name())
665  , edge_down(nullptr, typeid(*this).name())
666  , edge_up(nullptr, typeid(*this).name())
667 {
668  if (max_level == numbers::invalid_unsigned_int)
669  maxlevel = matrix.get_maxlevel();
670  else
671  maxlevel = max_level;
672  reinit(min_level, maxlevel);
673 }
674 
675 
676 
677 template <typename VectorType>
678 inline unsigned int
680 {
681  return maxlevel;
682 }
683 
684 
685 
686 template <typename VectorType>
687 inline unsigned int
689 {
690  return minlevel;
691 }
692 
693 
694 /* --------------------------- inline functions --------------------- */
695 
696 
697 namespace internal
698 {
699  namespace PreconditionMGImplementation
700  {
701  template <int dim,
702  typename VectorType,
703  class TRANSFER,
704  typename OtherVectorType>
705  typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
706  vmult(
707  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
709  const TRANSFER & transfer,
710  OtherVectorType & dst,
711  const OtherVectorType & src,
712  const bool uses_dof_handler_vector,
713  const typename ::mg::Signals &signals,
714  int)
715  {
716  signals.transfer_to_mg(true);
717  if (uses_dof_handler_vector)
718  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
719  else
720  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
721  signals.transfer_to_mg(false);
722 
723  multigrid.cycle();
724 
725  signals.transfer_to_global(true);
726  if (uses_dof_handler_vector)
727  transfer.copy_from_mg(dof_handler_vector, dst, multigrid.solution);
728  else
729  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
730  signals.transfer_to_global(false);
731  }
732 
733  template <int dim,
734  typename VectorType,
735  class TRANSFER,
736  typename OtherVectorType>
737  void
738  vmult(
739  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
740  ::Multigrid<VectorType> & multigrid,
741  const TRANSFER & transfer,
742  OtherVectorType & dst,
743  const OtherVectorType & src,
744  const bool uses_dof_handler_vector,
745  const typename ::mg::Signals &signals,
746  ...)
747  {
748  (void)uses_dof_handler_vector;
749  Assert(!uses_dof_handler_vector, ExcInternalError());
750 
751  signals.transfer_to_mg(true);
752  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
753  signals.transfer_to_mg(false);
754 
755  multigrid.cycle();
756 
757  signals.transfer_to_global(true);
758  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
759  signals.transfer_to_global(false);
760  }
761 
762  template <int dim,
763  typename VectorType,
764  class TRANSFER,
765  typename OtherVectorType>
766  typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
767  vmult_add(
768  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
769  ::Multigrid<VectorType> & multigrid,
770  const TRANSFER & transfer,
771  OtherVectorType & dst,
772  const OtherVectorType & src,
773  const bool uses_dof_handler_vector,
774  const typename ::mg::Signals &signals,
775  int)
776  {
777  signals.transfer_to_mg(true);
778  if (uses_dof_handler_vector)
779  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
780  else
781  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
782  signals.transfer_to_mg(false);
783 
784  multigrid.cycle();
785 
786  signals.transfer_to_global(true);
787  if (uses_dof_handler_vector)
788  transfer.copy_from_mg_add(dof_handler_vector, dst, multigrid.solution);
789  else
790  transfer.copy_from_mg_add(*dof_handler_vector[0],
791  dst,
792  multigrid.solution);
793  signals.transfer_to_global(false);
794  }
795 
796  template <int dim,
797  typename VectorType,
798  class TRANSFER,
799  typename OtherVectorType>
800  void
801  vmult_add(
802  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
803  ::Multigrid<VectorType> & multigrid,
804  const TRANSFER & transfer,
805  OtherVectorType & dst,
806  const OtherVectorType & src,
807  const bool uses_dof_handler_vector,
808  const typename ::mg::Signals &signals,
809  ...)
810  {
811  (void)uses_dof_handler_vector;
812  Assert(!uses_dof_handler_vector, ExcInternalError());
813 
814  signals.transfer_to_mg(true);
815  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
816  signals.transfer_to_mg(false);
817 
818  multigrid.cycle();
819 
820  signals.transfer_to_global(true);
821  transfer.copy_from_mg_add(*dof_handler_vector[0],
822  dst,
823  multigrid.solution);
824  signals.transfer_to_global(false);
825  }
826  } // namespace PreconditionMGImplementation
827 } // namespace internal
828 
829 template <int dim, typename VectorType, class TRANSFER>
831  const DoFHandler<dim> &dof_handler,
833  const TRANSFER & transfer)
834  : dof_handler_vector(1, &dof_handler)
835  , dof_handler_vector_raw(1, &dof_handler)
836  , multigrid(&mg)
837  , transfer(&transfer)
838  , uses_dof_handler_vector(false)
839 {}
840 
841 template <int dim, typename VectorType, class TRANSFER>
843  const std::vector<const DoFHandler<dim> *> &dof_handler,
845  const TRANSFER & transfer)
846  : dof_handler_vector(dof_handler.size())
847  , dof_handler_vector_raw(dof_handler.size())
848  , multigrid(&mg)
849  , transfer(&transfer)
851 {
852  for (unsigned int i = 0; i < dof_handler.size(); ++i)
853  {
854  dof_handler_vector[i] = dof_handler[i];
855  dof_handler_vector_raw[i] = dof_handler[i];
856  }
857 }
858 
859 template <int dim, typename VectorType, class TRANSFER>
860 inline bool
862 {
863  return false;
864 }
865 
866 template <int dim, typename VectorType, class TRANSFER>
867 template <class OtherVectorType>
868 void
870  OtherVectorType & dst,
871  const OtherVectorType &src) const
872 {
873  internal::PreconditionMGImplementation::vmult(dof_handler_vector_raw,
874  *multigrid,
875  *transfer,
876  dst,
877  src,
879  this->signals,
880  0);
881 }
882 
883 
884 template <int dim, typename VectorType, class TRANSFER>
885 IndexSet
887  const unsigned int block) const
888 {
889  AssertIndexRange(block, dof_handler_vector.size());
890  return dof_handler_vector[block]->locally_owned_dofs();
891 }
892 
893 
894 template <int dim, typename VectorType, class TRANSFER>
895 IndexSet
897  const unsigned int block) const
898 {
899  AssertIndexRange(block, dof_handler_vector.size());
900  return dof_handler_vector[block]->locally_owned_dofs();
901 }
902 
903 
904 
905 template <int dim, typename VectorType, class TRANSFER>
906 MPI_Comm
908 {
909  // currently parallel GMG works with parallel triangulations only,
910  // so it should be a safe bet to use it to query MPI communicator:
911  const Triangulation<dim> &tria = dof_handler_vector[0]->get_triangulation();
912  const parallel::TriangulationBase<dim> *ptria =
913  dynamic_cast<const parallel::TriangulationBase<dim> *>(&tria);
914  Assert(ptria != nullptr, ExcInternalError());
915  return ptria->get_communicator();
916 }
917 
918 
919 
920 template <int dim, typename VectorType, class TRANSFER>
921 boost::signals2::connection
923  const std::function<void(bool)> &slot)
924 {
925  return this->signals.transfer_to_mg.connect(slot);
926 }
927 
928 
929 
930 template <int dim, typename VectorType, class TRANSFER>
931 boost::signals2::connection
933  const std::function<void(bool)> &slot)
934 {
935  return this->signals.transfer_to_global.connect(slot);
936 }
937 
938 
939 
940 template <int dim, typename VectorType, class TRANSFER>
941 template <class OtherVectorType>
942 void
944  OtherVectorType & dst,
945  const OtherVectorType &src) const
946 {
947  internal::PreconditionMGImplementation::vmult_add(dof_handler_vector_raw,
948  *multigrid,
949  *transfer,
950  dst,
951  src,
953  this->signals,
954  0);
955 }
956 
957 
958 template <int dim, typename VectorType, class TRANSFER>
959 template <class OtherVectorType>
960 void
962  const OtherVectorType &) const
963 {
964  Assert(false, ExcNotImplemented());
965 }
966 
967 
968 template <int dim, typename VectorType, class TRANSFER>
969 template <class OtherVectorType>
970 void
972  OtherVectorType &,
973  const OtherVectorType &) const
974 {
975  Assert(false, ExcNotImplemented());
976 }
977 
978 
979 template <int dim, typename VectorType, class TRANSFER>
982 {
983  return *this->multigrid;
984 }
985 
986 
987 template <int dim, typename VectorType, class TRANSFER>
988 const Multigrid<VectorType> &
990 {
991  return *this->multigrid;
992 }
993 
994 #endif // DOXYGEN
995 
997 
998 #endif
boost::signals2::signal< void(const bool before, const unsigned int level)> coarse_solve
Definition: multigrid.h:86
boost::signals2::connection connect_transfer_to_global(const std::function< void(bool)> &slot)
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > matrix
Definition: multigrid.h:420
void cycle()
static const unsigned int invalid_unsigned_int
Definition: types.h:196
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_down
Definition: multigrid.h:466
SmartPointer< const MGMatrixBase< VectorType > > edge_in
Definition: multigrid.h:459
void vmult(OtherVectorType &dst, const OtherVectorType &src) const
MPI_Comm get_mpi_communicator() const
Contents is actually a matrix.
Definition: mg.h:81
unsigned int get_maxlevel() const
const ::Triangulation< dim, spacedim > & tria
boost::signals2::signal< void(const bool before, const unsigned int level)> prolongation
Definition: multigrid.h:105
#define AssertIndexRange(index, range)
Definition: exceptions.h:1720
SmartPointer< const MGCoarseGridBase< VectorType >, Multigrid< VectorType > > coarse
Definition: multigrid.h:426
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)> restriction
Definition: multigrid.h:97
MGLevelObject< VectorType > t
Definition: multigrid.h:409
bool empty() const
const bool uses_dof_handler_vector
Definition: multigrid.h:634
virtual MPI_Comm get_communicator() const override
Definition: tria_base.cc:142
boost::signals2::signal< void(const bool before, const unsigned int level)> residual_step
Definition: multigrid.h:132
mg::Signals signals
Definition: multigrid.h:357
unsigned int maxlevel
Definition: multigrid.h:391
IndexSet locally_owned_range_indices(const unsigned int block=0) const
#define Assert(cond, exc)
Definition: exceptions.h:1461
void Tvmult(OtherVectorType &dst, const OtherVectorType &src) const
unsigned int get_minlevel() const
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
unsigned int level
Definition: grid_out.cc:4614
SmartPointer< Multigrid< VectorType >, PreconditionMG< dim, VectorType, TRANSFER > > multigrid
Definition: multigrid.h:622
boost::signals2::signal< void(const bool before)> transfer_to_mg
Definition: multigrid.h:67
SmartPointer< const MGMatrixBase< VectorType > > edge_out
Definition: multigrid.h:451
void vmult_add(OtherVectorType &dst, const OtherVectorType &src) const
The W-cycle.
Definition: multigrid.h:177
boost::signals2::signal< void(const bool before, const unsigned int level)> pre_smoother_step
Definition: multigrid.h:116
The V-cycle.
Definition: multigrid.h:175
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:386
SmartPointer< const TRANSFER, PreconditionMG< dim, VectorType, TRANSFER > > transfer
Definition: multigrid.h:628
MGLevelObject< VectorType > solution
Definition: multigrid.h:403
Multigrid< VectorType > & get_multigrid()
std::vector< SmartPointer< const DoFHandler< dim >, PreconditionMG< dim, VectorType, TRANSFER > > > dof_handler_vector
Definition: multigrid.h:610
virtual unsigned int get_maxlevel() const =0
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > post_smooth
Definition: multigrid.h:444
IndexSet locally_owned_domain_indices(const unsigned int block=0) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
std::vector< const DoFHandler< dim > * > dof_handler_vector_raw
Definition: multigrid.h:616
boost::signals2::signal< void(const bool before, const unsigned int level)> edge_prolongation
Definition: multigrid.h:139
mg::Signals signals
Definition: multigrid.h:639
static ::ExceptionBase & ExcNotImplemented()
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 > defect2
Definition: multigrid.h:414
MGLevelObject< VectorType > defect
Definition: multigrid.h:398
PreconditionMG(const DoFHandler< dim > &dof_handler, Multigrid< VectorType > &mg, const TRANSFER &transfer)
SmartPointer< const MGTransferBase< VectorType >, Multigrid< VectorType > > transfer
Definition: multigrid.h:432
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_up
Definition: multigrid.h:473
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > pre_smooth
Definition: multigrid.h:438
Cycle cycle_type
Definition: multigrid.h:381
static ::ExceptionBase & ExcInternalError()