Reference documentation for deal.II version Git bed997f895 2020-09-22 11:49:20 -0400
\(\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 
81  boost::signals2::signal<void(const bool before, const unsigned int level)>
83 
89  boost::signals2::signal<void(const bool before, const unsigned int level)>
91 
97  boost::signals2::signal<void(const bool before, const unsigned int level)>
99 
105  boost::signals2::signal<void(const bool before, const unsigned int level)>
107 
114  boost::signals2::signal<void(const bool before, const unsigned int level)>
116  };
117 } // namespace mg
118 
141 template <typename VectorType>
142 class Multigrid : public Subscriptor
143 {
144 public:
148  enum Cycle
149  {
155  f_cycle
156  };
157 
160 
171  const MGCoarseGridBase<VectorType> &coarse,
172  const MGTransferBase<VectorType> & transfer,
173  const MGSmootherBase<VectorType> & pre_smooth,
174  const MGSmootherBase<VectorType> & post_smooth,
175  const unsigned int minlevel = 0,
176  const unsigned int maxlevel = numbers::invalid_unsigned_int,
177  Cycle cycle = v_cycle);
178 
182  void
183  reinit(const unsigned int minlevel, const unsigned int maxlevel);
184 
189  void
190  cycle();
191 
202  void
203  vcycle();
204 
217  void
218  set_edge_matrices(const MGMatrixBase<VectorType> &edge_out,
219  const MGMatrixBase<VectorType> &edge_in);
220 
233  void
234  set_edge_flux_matrices(const MGMatrixBase<VectorType> &edge_down,
235  const MGMatrixBase<VectorType> &edge_up);
236 
240  unsigned int
241  get_maxlevel() const;
242 
246  unsigned int
247  get_minlevel() const;
248 
254  void
255  set_maxlevel(const unsigned int);
256 
272  void
273  set_minlevel(const unsigned int level, bool relative = false);
274 
278  void set_cycle(Cycle);
279 
283  boost::signals2::connection
284  connect_coarse_solve(
285  const std::function<void(const bool, const unsigned int)> &slot);
286 
290  boost::signals2::connection
291  connect_restriction(
292  const std::function<void(const bool, const unsigned int)> &slot);
293 
297  boost::signals2::connection
298  connect_prolongation(
299  const std::function<void(const bool, const unsigned int)> &slot);
300 
304  boost::signals2::connection
305  connect_pre_smoother_step(
306  const std::function<void(const bool, const unsigned int)> &slot);
307 
311  boost::signals2::connection
312  connect_post_smoother_step(
313  const std::function<void(const bool, const unsigned int)> &slot);
314 
315 private:
320 
327  void
328  level_v_step(const unsigned int level);
329 
337  void
338  level_step(const unsigned int level, Cycle cycle);
339 
344 
348  unsigned int minlevel;
349 
353  unsigned int maxlevel;
354 
355 public:
361 
366 
367 private:
372 
377 
378 
383 
389 
395 
401 
407 
414 
422 
429 
436 
437  template <int dim, class OtherVectorType, class TRANSFER>
438  friend class PreconditionMG;
439 };
440 
441 
456 template <int dim, typename VectorType, class TRANSFER>
458 {
459 public:
464  PreconditionMG(const DoFHandler<dim> &dof_handler,
466  const TRANSFER & transfer);
467 
472  PreconditionMG(const std::vector<const DoFHandler<dim> *> &dof_handler,
474  const TRANSFER & transfer);
475 
479  bool
480  empty() const;
481 
488  template <class OtherVectorType>
489  void
490  vmult(OtherVectorType &dst, const OtherVectorType &src) const;
491 
496  template <class OtherVectorType>
497  void
498  vmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
499 
505  template <class OtherVectorType>
506  void
507  Tvmult(OtherVectorType &dst, const OtherVectorType &src) const;
508 
514  template <class OtherVectorType>
515  void
516  Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
517 
524  IndexSet
525  locally_owned_range_indices(const unsigned int block = 0) const;
526 
533  IndexSet
534  locally_owned_domain_indices(const unsigned int block = 0) const;
535 
539  const MPI_Comm &
540  get_mpi_communicator() const;
541 
545  boost::signals2::connection
546  connect_transfer_to_mg(const std::function<void(bool)> &slot);
547 
551  boost::signals2::connection
552  connect_transfer_to_global(const std::function<void(bool)> &slot);
553 
554 private:
558  std::vector<SmartPointer<const DoFHandler<dim>,
561 
566  std::vector<const DoFHandler<dim> *> dof_handler_vector_raw;
567 
573 
579 
585 
590 };
591 
594 #ifndef DOXYGEN
595 /* --------------------------- inline functions --------------------- */
596 
597 
598 template <typename VectorType>
600  const MGCoarseGridBase<VectorType> &coarse,
601  const MGTransferBase<VectorType> & transfer,
602  const MGSmootherBase<VectorType> & pre_smooth,
603  const MGSmootherBase<VectorType> &post_smooth,
604  const unsigned int min_level,
605  const unsigned int max_level,
606  Cycle cycle)
607  : cycle_type(cycle)
608  , matrix(&matrix, typeid(*this).name())
609  , coarse(&coarse, typeid(*this).name())
610  , transfer(&transfer, typeid(*this).name())
611  , pre_smooth(&pre_smooth, typeid(*this).name())
612  , post_smooth(&post_smooth, typeid(*this).name())
613  , edge_out(nullptr, typeid(*this).name())
614  , edge_in(nullptr, typeid(*this).name())
615  , edge_down(nullptr, typeid(*this).name())
616  , edge_up(nullptr, typeid(*this).name())
617 {
618  if (max_level == numbers::invalid_unsigned_int)
619  maxlevel = matrix.get_maxlevel();
620  else
621  maxlevel = max_level;
622  reinit(min_level, maxlevel);
623 }
624 
625 
626 
627 template <typename VectorType>
628 inline unsigned int
630 {
631  return maxlevel;
632 }
633 
634 
635 
636 template <typename VectorType>
637 inline unsigned int
639 {
640  return minlevel;
641 }
642 
643 
644 /* --------------------------- inline functions --------------------- */
645 
646 
647 namespace internal
648 {
649  namespace PreconditionMGImplementation
650  {
651  template <int dim,
652  typename VectorType,
653  class TRANSFER,
654  typename OtherVectorType>
655  typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
656  vmult(
657  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
659  const TRANSFER & transfer,
660  OtherVectorType & dst,
661  const OtherVectorType & src,
662  const bool uses_dof_handler_vector,
663  const typename ::mg::Signals &signals,
664  int)
665  {
666  signals.transfer_to_mg(true);
667  if (uses_dof_handler_vector)
668  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
669  else
670  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
671  signals.transfer_to_mg(false);
672 
673  multigrid.cycle();
674 
675  signals.transfer_to_global(true);
676  if (uses_dof_handler_vector)
677  transfer.copy_from_mg(dof_handler_vector, dst, multigrid.solution);
678  else
679  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
680  signals.transfer_to_global(false);
681  }
682 
683  template <int dim,
684  typename VectorType,
685  class TRANSFER,
686  typename OtherVectorType>
687  void
688  vmult(
689  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
690  ::Multigrid<VectorType> & multigrid,
691  const TRANSFER & transfer,
692  OtherVectorType & dst,
693  const OtherVectorType & src,
694  const bool uses_dof_handler_vector,
695  const typename ::mg::Signals &signals,
696  ...)
697  {
698  (void)uses_dof_handler_vector;
699  Assert(!uses_dof_handler_vector, ExcInternalError());
700 
701  signals.transfer_to_mg(true);
702  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
703  signals.transfer_to_mg(false);
704 
705  multigrid.cycle();
706 
707  signals.transfer_to_global(true);
708  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
709  signals.transfer_to_global(false);
710  }
711 
712  template <int dim,
713  typename VectorType,
714  class TRANSFER,
715  typename OtherVectorType>
716  typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
717  vmult_add(
718  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
719  ::Multigrid<VectorType> & multigrid,
720  const TRANSFER & transfer,
721  OtherVectorType & dst,
722  const OtherVectorType & src,
723  const bool uses_dof_handler_vector,
724  const typename ::mg::Signals &signals,
725  int)
726  {
727  signals.transfer_to_mg(true);
728  if (uses_dof_handler_vector)
729  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
730  else
731  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
732  signals.transfer_to_mg(false);
733 
734  multigrid.cycle();
735 
736  signals.transfer_to_global(true);
737  if (uses_dof_handler_vector)
738  transfer.copy_from_mg_add(dof_handler_vector, dst, multigrid.solution);
739  else
740  transfer.copy_from_mg_add(*dof_handler_vector[0],
741  dst,
742  multigrid.solution);
743  signals.transfer_to_global(false);
744  }
745 
746  template <int dim,
747  typename VectorType,
748  class TRANSFER,
749  typename OtherVectorType>
750  void
751  vmult_add(
752  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
753  ::Multigrid<VectorType> & multigrid,
754  const TRANSFER & transfer,
755  OtherVectorType & dst,
756  const OtherVectorType & src,
757  const bool uses_dof_handler_vector,
758  const typename ::mg::Signals &signals,
759  ...)
760  {
761  (void)uses_dof_handler_vector;
762  Assert(!uses_dof_handler_vector, ExcInternalError());
763 
764  signals.transfer_to_mg(true);
765  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
766  signals.transfer_to_mg(false);
767 
768  multigrid.cycle();
769 
770  signals.transfer_to_global(true);
771  transfer.copy_from_mg_add(*dof_handler_vector[0],
772  dst,
773  multigrid.solution);
774  signals.transfer_to_global(false);
775  }
776  } // namespace PreconditionMGImplementation
777 } // namespace internal
778 
779 template <int dim, typename VectorType, class TRANSFER>
781  const DoFHandler<dim> &dof_handler,
783  const TRANSFER & transfer)
784  : dof_handler_vector(1, &dof_handler)
785  , dof_handler_vector_raw(1, &dof_handler)
786  , multigrid(&mg)
787  , transfer(&transfer)
788  , uses_dof_handler_vector(false)
789 {}
790 
791 template <int dim, typename VectorType, class TRANSFER>
793  const std::vector<const DoFHandler<dim> *> &dof_handler,
795  const TRANSFER & transfer)
796  : dof_handler_vector(dof_handler.size())
797  , dof_handler_vector_raw(dof_handler.size())
798  , multigrid(&mg)
799  , transfer(&transfer)
801 {
802  for (unsigned int i = 0; i < dof_handler.size(); ++i)
803  {
804  dof_handler_vector[i] = dof_handler[i];
805  dof_handler_vector_raw[i] = dof_handler[i];
806  }
807 }
808 
809 template <int dim, typename VectorType, class TRANSFER>
810 inline bool
812 {
813  return false;
814 }
815 
816 template <int dim, typename VectorType, class TRANSFER>
817 template <class OtherVectorType>
818 void
820  OtherVectorType & dst,
821  const OtherVectorType &src) const
822 {
823  internal::PreconditionMGImplementation::vmult(dof_handler_vector_raw,
824  *multigrid,
825  *transfer,
826  dst,
827  src,
829  this->signals,
830  0);
831 }
832 
833 
834 template <int dim, typename VectorType, class TRANSFER>
835 IndexSet
837  const unsigned int block) const
838 {
839  AssertIndexRange(block, dof_handler_vector.size());
840  return dof_handler_vector[block]->locally_owned_dofs();
841 }
842 
843 
844 template <int dim, typename VectorType, class TRANSFER>
845 IndexSet
847  const unsigned int block) const
848 {
849  AssertIndexRange(block, dof_handler_vector.size());
850  return dof_handler_vector[block]->locally_owned_dofs();
851 }
852 
853 
854 
855 template <int dim, typename VectorType, class TRANSFER>
856 const MPI_Comm &
858 {
859  // currently parallel GMG works with parallel triangulations only,
860  // so it should be a safe bet to use it to query MPI communicator:
861  const Triangulation<dim> &tria = dof_handler_vector[0]->get_triangulation();
862  const parallel::TriangulationBase<dim> *ptria =
863  dynamic_cast<const parallel::TriangulationBase<dim> *>(&tria);
864  Assert(ptria != nullptr, ExcInternalError());
865  return ptria->get_communicator();
866 }
867 
868 
869 
870 template <int dim, typename VectorType, class TRANSFER>
871 boost::signals2::connection
873  const std::function<void(bool)> &slot)
874 {
875  return this->signals.transfer_to_mg.connect(slot);
876 }
877 
878 
879 
880 template <int dim, typename VectorType, class TRANSFER>
881 boost::signals2::connection
883  const std::function<void(bool)> &slot)
884 {
885  return this->signals.transfer_to_global.connect(slot);
886 }
887 
888 
889 
890 template <int dim, typename VectorType, class TRANSFER>
891 template <class OtherVectorType>
892 void
894  OtherVectorType & dst,
895  const OtherVectorType &src) const
896 {
897  internal::PreconditionMGImplementation::vmult_add(dof_handler_vector_raw,
898  *multigrid,
899  *transfer,
900  dst,
901  src,
903  this->signals,
904  0);
905 }
906 
907 
908 template <int dim, typename VectorType, class TRANSFER>
909 template <class OtherVectorType>
910 void
912  const OtherVectorType &) const
913 {
914  Assert(false, ExcNotImplemented());
915 }
916 
917 
918 template <int dim, typename VectorType, class TRANSFER>
919 template <class OtherVectorType>
920 void
922  OtherVectorType &,
923  const OtherVectorType &) const
924 {
925  Assert(false, ExcNotImplemented());
926 }
927 
928 #endif // DOXYGEN
929 
931 
932 #endif
boost::signals2::signal< void(const bool before, const unsigned int level)> coarse_solve
Definition: multigrid.h:82
boost::signals2::connection connect_transfer_to_global(const std::function< void(bool)> &slot)
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > matrix
Definition: multigrid.h:382
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:428
SmartPointer< const MGMatrixBase< VectorType > > edge_in
Definition: multigrid.h:421
void vmult(OtherVectorType &dst, const OtherVectorType &src) const
Contents is actually a matrix.
Definition: mg.h:81
unsigned int get_maxlevel() const
boost::signals2::signal< void(const bool before, const unsigned int level)> prolongation
Definition: multigrid.h:98
#define AssertIndexRange(index, range)
Definition: exceptions.h:1636
SmartPointer< const MGCoarseGridBase< VectorType >, Multigrid< VectorType > > coarse
Definition: multigrid.h:388
boost::signals2::signal< void(const bool before, const unsigned int level)> post_smoother_step
Definition: multigrid.h:115
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:90
MGLevelObject< VectorType > t
Definition: multigrid.h:371
bool empty() const
const bool uses_dof_handler_vector
Definition: multigrid.h:584
mg::Signals signals
Definition: multigrid.h:319
unsigned int maxlevel
Definition: multigrid.h:353
IndexSet locally_owned_range_indices(const unsigned int block=0) const
#define Assert(cond, exc)
Definition: exceptions.h:1411
void Tvmult(OtherVectorType &dst, const OtherVectorType &src) const
unsigned int get_minlevel() const
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
unsigned int level
Definition: grid_out.cc:4338
SmartPointer< Multigrid< VectorType >, PreconditionMG< dim, VectorType, TRANSFER > > multigrid
Definition: multigrid.h:572
boost::signals2::signal< void(const bool before)> transfer_to_mg
Definition: multigrid.h:67
const MPI_Comm & get_mpi_communicator() const
SmartPointer< const MGMatrixBase< VectorType > > edge_out
Definition: multigrid.h:413
void vmult_add(OtherVectorType &dst, const OtherVectorType &src) const
The W-cycle.
Definition: multigrid.h:153
boost::signals2::signal< void(const bool before, const unsigned int level)> pre_smoother_step
Definition: multigrid.h:106
The V-cycle.
Definition: multigrid.h:151
virtual const MPI_Comm & get_communicator() const
Definition: tria_base.cc:138
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:348
SmartPointer< const TRANSFER, PreconditionMG< dim, VectorType, TRANSFER > > transfer
Definition: multigrid.h:578
MGLevelObject< VectorType > solution
Definition: multigrid.h:365
std::vector< SmartPointer< const DoFHandler< dim >, PreconditionMG< dim, VectorType, TRANSFER > > > dof_handler_vector
Definition: multigrid.h:560
virtual unsigned int get_maxlevel() const =0
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > post_smooth
Definition: multigrid.h:406
IndexSet locally_owned_domain_indices(const unsigned int block=0) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
std::vector< const DoFHandler< dim > * > dof_handler_vector_raw
Definition: multigrid.h:566
mg::Signals signals
Definition: multigrid.h:589
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:376
MGLevelObject< VectorType > defect
Definition: multigrid.h:360
PreconditionMG(const DoFHandler< dim > &dof_handler, Multigrid< VectorType > &mg, const TRANSFER &transfer)
SmartPointer< const MGTransferBase< VectorType >, Multigrid< VectorType > > transfer
Definition: multigrid.h:394
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_up
Definition: multigrid.h:435
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > pre_smooth
Definition: multigrid.h:400
Cycle cycle_type
Definition: multigrid.h:343
static ::ExceptionBase & ExcInternalError()