Reference documentation for deal.II version GIT 9042b9283b 2023-12-02 14:50:02+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 - 2023 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 
51 namespace mg
52 {
62  struct Signals
63  {
69  boost::signals2::signal<void(const bool before)> transfer_to_mg;
70 
76  boost::signals2::signal<void(const bool before)> transfer_to_global;
77 
87  boost::signals2::signal<void(const bool before, const unsigned int level)>
89 
98  boost::signals2::signal<void(const bool before, const unsigned int level)>
100 
106  boost::signals2::signal<void(const bool before, const unsigned int level)>
108 
117  boost::signals2::signal<void(const bool before, const unsigned int level)>
119 
125  boost::signals2::signal<void(const bool before, const unsigned int level)>
127 
133  boost::signals2::signal<void(const bool before, const unsigned int level)>
135 
140  boost::signals2::signal<void(const bool before, const unsigned int level)>
142  };
143 } // namespace mg
144 
162 template <typename VectorType>
163 class Multigrid : public Subscriptor
164 {
165 public:
169  enum Cycle
170  {
176  f_cycle
177  };
178 
179  using vector_type = VectorType;
180  using const_vector_type = const VectorType;
181 
196  const unsigned int minlevel = 0,
197  const unsigned int maxlevel = numbers::invalid_unsigned_int,
198  Cycle cycle = v_cycle);
199 
203  void
204  reinit(const unsigned int minlevel, const unsigned int maxlevel);
205 
210  void
211  cycle();
212 
223  void
225 
238  void
241 
248  void
250 
263  void
266 
270  unsigned int
271  get_maxlevel() const;
272 
276  unsigned int
277  get_minlevel() const;
278 
284  void
285  set_maxlevel(const unsigned int);
286 
302  void
303  set_minlevel(const unsigned int level, bool relative = false);
304 
309 
313  boost::signals2::connection
315  const std::function<void(const bool, const unsigned int)> &slot);
316 
320  boost::signals2::connection
322  const std::function<void(const bool, const unsigned int)> &slot);
323 
327  boost::signals2::connection
329  const std::function<void(const bool, const unsigned int)> &slot);
330 
334  boost::signals2::connection
336  const std::function<void(const bool, const unsigned int)> &slot);
337 
341  boost::signals2::connection
343  const std::function<void(const bool, const unsigned int)> &slot);
344 
348  boost::signals2::connection
350  const std::function<void(const bool, const unsigned int)> &slot);
351 
355  boost::signals2::connection
357  const std::function<void(const bool, const unsigned int)> &slot);
358 
359 private:
364 
371  void
372  level_v_step(const unsigned int level);
373 
381  void
382  level_step(const unsigned int level, Cycle cycle);
383 
388 
392  unsigned int minlevel;
393 
397  unsigned int maxlevel;
398 
399 public:
405 
410 
411 private:
416 
421 
422 
427 
433 
439 
445 
451 
458 
466 
473 
480 
481  template <int dim, typename OtherVectorType, typename TransferType>
482  friend class PreconditionMG;
483 };
484 
485 
500 template <int dim, typename VectorType, typename TransferType>
502 {
503 public:
508  PreconditionMG(const DoFHandler<dim> &dof_handler,
510  const TransferType &transfer);
511 
516  PreconditionMG(const std::vector<const DoFHandler<dim> *> &dof_handler,
518  const TransferType &transfer);
519 
523  bool
524  empty() const;
525 
532  template <typename OtherVectorType>
533  void
534  vmult(OtherVectorType &dst, const OtherVectorType &src) const;
535 
540  template <typename OtherVectorType>
541  void
542  vmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
543 
549  template <typename OtherVectorType>
550  void
551  Tvmult(OtherVectorType &dst, const OtherVectorType &src) const;
552 
558  template <typename OtherVectorType>
559  void
560  Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
561 
568  IndexSet
569  locally_owned_range_indices(const unsigned int block = 0) const;
570 
577  IndexSet
578  locally_owned_domain_indices(const unsigned int block = 0) const;
579 
583  MPI_Comm
585 
589  boost::signals2::connection
590  connect_transfer_to_mg(const std::function<void(bool)> &slot);
591 
595  boost::signals2::connection
596  connect_transfer_to_global(const std::function<void(bool)> &slot);
597 
603 
607  const Multigrid<VectorType> &
608  get_multigrid() const;
609 
610 private:
614  std::vector<SmartPointer<const DoFHandler<dim>,
617 
622  std::vector<const DoFHandler<dim> *> dof_handler_vector_raw;
623 
630 
634  SmartPointer<const TransferType,
637 
643 
648 };
649 
652 #ifndef DOXYGEN
653 /* --------------------------- inline functions --------------------- */
654 
655 
656 template <typename VectorType>
658  const MGCoarseGridBase<VectorType> &coarse,
659  const MGTransferBase<VectorType> &transfer,
660  const MGSmootherBase<VectorType> &pre_smooth,
661  const MGSmootherBase<VectorType> &post_smooth,
662  const unsigned int min_level,
663  const unsigned int max_level,
664  Cycle cycle)
665  : cycle_type(cycle)
666  , matrix(&matrix, typeid(*this).name())
667  , coarse(&coarse, typeid(*this).name())
668  , transfer(&transfer, typeid(*this).name())
669  , pre_smooth(&pre_smooth, typeid(*this).name())
670  , post_smooth(&post_smooth, typeid(*this).name())
671  , edge_out(nullptr, typeid(*this).name())
672  , edge_in(nullptr, typeid(*this).name())
673  , edge_down(nullptr, typeid(*this).name())
674  , edge_up(nullptr, typeid(*this).name())
675 {
676  if (max_level == numbers::invalid_unsigned_int)
677  maxlevel = matrix.get_maxlevel();
678  else
679  maxlevel = max_level;
680  reinit(min_level, maxlevel);
681 }
682 
683 
684 
685 template <typename VectorType>
686 inline unsigned int
688 {
689  return maxlevel;
690 }
691 
692 
693 
694 template <typename VectorType>
695 inline unsigned int
697 {
698  return minlevel;
699 }
700 
701 
702 /* --------------------------- inline functions --------------------- */
703 
704 
705 namespace internal
706 {
707  namespace PreconditionMGImplementation
708  {
709  template <int dim,
710  typename VectorType,
711  typename TransferType,
712  typename OtherVectorType>
713  std::enable_if_t<TransferType::supports_dof_handler_vector>
714  vmult(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
715  Multigrid<VectorType> &multigrid,
716  const TransferType &transfer,
717  OtherVectorType &dst,
718  const OtherVectorType &src,
719  const bool uses_dof_handler_vector,
720  const typename ::mg::Signals &signals,
721  int)
722  {
723  signals.transfer_to_mg(true);
724  if (uses_dof_handler_vector)
725  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
726  else
727  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
728  signals.transfer_to_mg(false);
729 
730  multigrid.cycle();
731 
732  signals.transfer_to_global(true);
733  if (uses_dof_handler_vector)
734  transfer.copy_from_mg(dof_handler_vector, dst, multigrid.solution);
735  else
736  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
737  signals.transfer_to_global(false);
738  }
739 
740  template <int dim,
741  typename VectorType,
742  typename TransferType,
743  typename OtherVectorType>
744  void
745  vmult(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
746  Multigrid<VectorType> &multigrid,
747  const TransferType &transfer,
748  OtherVectorType &dst,
749  const OtherVectorType &src,
750  const bool uses_dof_handler_vector,
751  const typename ::mg::Signals &signals,
752  ...)
753  {
754  (void)uses_dof_handler_vector;
755  Assert(!uses_dof_handler_vector, ExcInternalError());
756 
757  signals.transfer_to_mg(true);
758  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
759  signals.transfer_to_mg(false);
760 
761  multigrid.cycle();
762 
763  signals.transfer_to_global(true);
764  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
765  signals.transfer_to_global(false);
766  }
767 
768  template <int dim,
769  typename VectorType,
770  typename TransferType,
771  typename OtherVectorType>
772  std::enable_if_t<TransferType::supports_dof_handler_vector>
773  vmult_add(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
774  Multigrid<VectorType> &multigrid,
775  const TransferType &transfer,
776  OtherVectorType &dst,
777  const OtherVectorType &src,
778  const bool uses_dof_handler_vector,
779  const typename ::mg::Signals &signals,
780  int)
781  {
782  signals.transfer_to_mg(true);
783  if (uses_dof_handler_vector)
784  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
785  else
786  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
787  signals.transfer_to_mg(false);
788 
789  multigrid.cycle();
790 
791  signals.transfer_to_global(true);
792  if (uses_dof_handler_vector)
793  transfer.copy_from_mg_add(dof_handler_vector, dst, multigrid.solution);
794  else
795  transfer.copy_from_mg_add(*dof_handler_vector[0],
796  dst,
797  multigrid.solution);
798  signals.transfer_to_global(false);
799  }
800 
801  template <int dim,
802  typename VectorType,
803  typename TransferType,
804  typename OtherVectorType>
805  void
806  vmult_add(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
807  Multigrid<VectorType> &multigrid,
808  const TransferType &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, typename TransferType>
835  const DoFHandler<dim> &dof_handler,
837  const TransferType &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, typename TransferType>
847  const std::vector<const DoFHandler<dim> *> &dof_handler,
849  const TransferType &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, typename TransferType>
864 inline bool
866 {
867  return false;
868 }
869 
870 template <int dim, typename VectorType, typename TransferType>
871 template <typename 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, typename TransferType>
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, typename TransferType>
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, typename TransferType>
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, typename TransferType>
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, typename TransferType>
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, typename TransferType>
945 template <typename 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, typename TransferType>
963 template <typename OtherVectorType>
964 void
966  OtherVectorType &,
967  const OtherVectorType &) const
968 {
969  Assert(false, ExcNotImplemented());
970 }
971 
972 
973 template <int dim, typename VectorType, typename TransferType>
974 template <typename OtherVectorType>
975 void
977  OtherVectorType &,
978  const OtherVectorType &) const
979 {
980  Assert(false, ExcNotImplemented());
981 }
982 
983 
984 template <int dim, typename VectorType, typename TransferType>
987 {
988  return *this->multigrid;
989 }
990 
991 
992 template <int dim, typename VectorType, typename TransferType>
993 const Multigrid<VectorType> &
995 {
996  return *this->multigrid;
997 }
998 
999 #endif // DOXYGEN
1000 
1002 
1003 #endif
@ v_cycle
The V-cycle.
Definition: multigrid.h:172
@ w_cycle
The W-cycle.
Definition: multigrid.h:174
@ f_cycle
The F-cycle.
Definition: multigrid.h:176
MGLevelObject< VectorType > defect2
Definition: multigrid.h:420
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:438
unsigned int minlevel
Definition: multigrid.h:392
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_up
Definition: multigrid.h:479
mg::Signals signals
Definition: multigrid.h:363
SmartPointer< const MGMatrixBase< VectorType > > edge_in
Definition: multigrid.h:465
SmartPointer< const MGMatrixBase< VectorType > > edge_out
Definition: multigrid.h:457
MGLevelObject< VectorType > solution
Definition: multigrid.h:409
boost::signals2::connection connect_edge_prolongation(const std::function< void(const bool, const unsigned int)> &slot)
unsigned int maxlevel
Definition: multigrid.h:397
VectorType vector_type
Definition: multigrid.h:179
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > matrix
Definition: multigrid.h:426
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:180
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > pre_smooth
Definition: multigrid.h:444
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > post_smooth
Definition: multigrid.h:450
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:432
Cycle cycle_type
Definition: multigrid.h:387
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:472
void set_edge_in_matrix(const MGMatrixBase< VectorType > &edge_in)
MGLevelObject< VectorType > t
Definition: multigrid.h:415
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:404
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)
std::vector< const DoFHandler< dim > * > dof_handler_vector_raw
Definition: multigrid.h:622
void Tvmult(OtherVectorType &dst, const OtherVectorType &src) const
SmartPointer< const TransferType, PreconditionMG< dim, VectorType, TransferType > > transfer
Definition: multigrid.h:636
IndexSet locally_owned_domain_indices(const unsigned int block=0) const
boost::signals2::connection connect_transfer_to_global(const std::function< void(bool)> &slot)
mg::Signals signals
Definition: multigrid.h:647
boost::signals2::connection connect_transfer_to_mg(const std::function< void(bool)> &slot)
PreconditionMG(const DoFHandler< dim > &dof_handler, Multigrid< VectorType > &mg, const TransferType &transfer)
MPI_Comm get_mpi_communicator() const
const bool uses_dof_handler_vector
Definition: multigrid.h:642
void vmult_add(OtherVectorType &dst, const OtherVectorType &src) const
bool empty() const
IndexSet locally_owned_range_indices(const unsigned int block=0) const
const Multigrid< VectorType > & get_multigrid() const
SmartPointer< Multigrid< VectorType >, PreconditionMG< dim, VectorType, TransferType > > multigrid
Definition: multigrid.h:629
PreconditionMG(const std::vector< const DoFHandler< dim > * > &dof_handler, Multigrid< VectorType > &mg, const TransferType &transfer)
std::vector< SmartPointer< const DoFHandler< dim >, PreconditionMG< dim, VectorType, TransferType > > > dof_handler_vector
Definition: multigrid.h:616
void Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const
void vmult(OtherVectorType &dst, const OtherVectorType &src) const
Multigrid< VectorType > & get_multigrid()
Triangulation< dim, spacedim > & get_triangulation()
virtual MPI_Comm get_communicator() const override
Definition: tria_base.cc:160
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
unsigned int level
Definition: grid_out.cc:4617
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
Definition: exceptions.h:1631
static ::ExceptionBase & ExcNotImplemented()
#define AssertIndexRange(index, range)
Definition: exceptions.h:1888
@ 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:221
boost::signals2::signal< void(const bool before, const unsigned int level)> prolongation
Definition: multigrid.h:107
boost::signals2::signal< void(const bool before, const unsigned int level)> pre_smoother_step
Definition: multigrid.h:118
boost::signals2::signal< void(const bool before, const unsigned int level)> coarse_solve
Definition: multigrid.h:88
boost::signals2::signal< void(const bool before, const unsigned int level)> edge_prolongation
Definition: multigrid.h:141
boost::signals2::signal< void(const bool before, const unsigned int level)> post_smoother_step
Definition: multigrid.h:126
boost::signals2::signal< void(const bool before)> transfer_to_global
Definition: multigrid.h:76
boost::signals2::signal< void(const bool before, const unsigned int level)> residual_step
Definition: multigrid.h:134
boost::signals2::signal< void(const bool before, const unsigned int level)> restriction
Definition: multigrid.h:99
boost::signals2::signal< void(const bool before)> transfer_to_mg
Definition: multigrid.h:69
const ::Triangulation< dim, spacedim > & tria