Reference documentation for deal.II version Git 08727cc441 2020-07-02 15:45:42 -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 
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 
134 template <typename VectorType>
135 class Multigrid : public Subscriptor
136 {
137 public:
141  enum Cycle
142  {
148  f_cycle
149  };
150 
153 
164  const MGCoarseGridBase<VectorType> &coarse,
165  const MGTransferBase<VectorType> & transfer,
166  const MGSmootherBase<VectorType> & pre_smooth,
167  const MGSmootherBase<VectorType> & post_smooth,
168  const unsigned int minlevel = 0,
169  const unsigned int maxlevel = numbers::invalid_unsigned_int,
170  Cycle cycle = v_cycle);
171 
175  void
176  reinit(const unsigned int minlevel, const unsigned int maxlevel);
177 
182  void
183  cycle();
184 
195  void
196  vcycle();
197 
210  void
211  set_edge_matrices(const MGMatrixBase<VectorType> &edge_out,
212  const MGMatrixBase<VectorType> &edge_in);
213 
226  void
227  set_edge_flux_matrices(const MGMatrixBase<VectorType> &edge_down,
228  const MGMatrixBase<VectorType> &edge_up);
229 
233  unsigned int
234  get_maxlevel() const;
235 
239  unsigned int
240  get_minlevel() const;
241 
247  void
248  set_maxlevel(const unsigned int);
249 
265  void
266  set_minlevel(const unsigned int level, bool relative = false);
267 
271  void set_cycle(Cycle);
272 
276  boost::signals2::connection
277  connect_coarse_solve(
278  const std::function<void(const bool, const unsigned int)> &slot);
279 
283  boost::signals2::connection
284  connect_restriction(
285  const std::function<void(const bool, const unsigned int)> &slot);
286 
290  boost::signals2::connection
291  connect_prolongation(
292  const std::function<void(const bool, const unsigned int)> &slot);
293 
297  boost::signals2::connection
298  connect_pre_smoother_step(
299  const std::function<void(const bool, const unsigned int)> &slot);
300 
304  boost::signals2::connection
305  connect_post_smoother_step(
306  const std::function<void(const bool, const unsigned int)> &slot);
307 
308 private:
313 
320  void
321  level_v_step(const unsigned int level);
322 
330  void
331  level_step(const unsigned int level, Cycle cycle);
332 
337 
341  unsigned int minlevel;
342 
346  unsigned int maxlevel;
347 
348 public:
354 
359 
360 private:
365 
370 
371 
376 
382 
388 
394 
400 
407 
415 
422 
429 
430  template <int dim, class OtherVectorType, class TRANSFER>
431  friend class PreconditionMG;
432 };
433 
434 
449 template <int dim, typename VectorType, class TRANSFER>
451 {
452 public:
457  PreconditionMG(const DoFHandler<dim> &dof_handler,
459  const TRANSFER & transfer);
460 
465  PreconditionMG(const std::vector<const DoFHandler<dim> *> &dof_handler,
467  const TRANSFER & transfer);
468 
472  bool
473  empty() const;
474 
481  template <class OtherVectorType>
482  void
483  vmult(OtherVectorType &dst, const OtherVectorType &src) const;
484 
489  template <class OtherVectorType>
490  void
491  vmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
492 
498  template <class OtherVectorType>
499  void
500  Tvmult(OtherVectorType &dst, const OtherVectorType &src) const;
501 
507  template <class OtherVectorType>
508  void
509  Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
510 
517  IndexSet
518  locally_owned_range_indices(const unsigned int block = 0) const;
519 
526  IndexSet
527  locally_owned_domain_indices(const unsigned int block = 0) const;
528 
532  const MPI_Comm &
533  get_mpi_communicator() const;
534 
538  boost::signals2::connection
539  connect_transfer_to_mg(const std::function<void(bool)> &slot);
540 
544  boost::signals2::connection
545  connect_transfer_to_global(const std::function<void(bool)> &slot);
546 
547 private:
551  std::vector<SmartPointer<const DoFHandler<dim>,
554 
559  std::vector<const DoFHandler<dim> *> dof_handler_vector_raw;
560 
566 
572 
578 
583 };
584 
587 #ifndef DOXYGEN
588 /* --------------------------- inline functions --------------------- */
589 
590 
591 template <typename VectorType>
593  const MGCoarseGridBase<VectorType> &coarse,
594  const MGTransferBase<VectorType> & transfer,
595  const MGSmootherBase<VectorType> & pre_smooth,
596  const MGSmootherBase<VectorType> &post_smooth,
597  const unsigned int min_level,
598  const unsigned int max_level,
599  Cycle cycle)
600  : cycle_type(cycle)
601  , matrix(&matrix, typeid(*this).name())
602  , coarse(&coarse, typeid(*this).name())
603  , transfer(&transfer, typeid(*this).name())
604  , pre_smooth(&pre_smooth, typeid(*this).name())
605  , post_smooth(&post_smooth, typeid(*this).name())
606  , edge_out(nullptr, typeid(*this).name())
607  , edge_in(nullptr, typeid(*this).name())
608  , edge_down(nullptr, typeid(*this).name())
609  , edge_up(nullptr, typeid(*this).name())
610 {
611  if (max_level == numbers::invalid_unsigned_int)
612  maxlevel = matrix.get_maxlevel();
613  else
614  maxlevel = max_level;
615  reinit(min_level, maxlevel);
616 }
617 
618 
619 
620 template <typename VectorType>
621 inline unsigned int
623 {
624  return maxlevel;
625 }
626 
627 
628 
629 template <typename VectorType>
630 inline unsigned int
632 {
633  return minlevel;
634 }
635 
636 
637 /* --------------------------- inline functions --------------------- */
638 
639 
640 namespace internal
641 {
642  namespace PreconditionMGImplementation
643  {
644  template <int dim,
645  typename VectorType,
646  class TRANSFER,
647  typename OtherVectorType>
648  typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
649  vmult(
650  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
652  const TRANSFER & transfer,
653  OtherVectorType & dst,
654  const OtherVectorType & src,
655  const bool uses_dof_handler_vector,
656  const typename ::mg::Signals &signals,
657  int)
658  {
659  signals.transfer_to_mg(true);
660  if (uses_dof_handler_vector)
661  transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
662  else
663  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
664  signals.transfer_to_mg(false);
665 
666  multigrid.cycle();
667 
668  signals.transfer_to_global(true);
669  if (uses_dof_handler_vector)
670  transfer.copy_from_mg(dof_handler_vector, dst, multigrid.solution);
671  else
672  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
673  signals.transfer_to_global(false);
674  }
675 
676  template <int dim,
677  typename VectorType,
678  class TRANSFER,
679  typename OtherVectorType>
680  void
681  vmult(
682  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
683  ::Multigrid<VectorType> & multigrid,
684  const TRANSFER & transfer,
685  OtherVectorType & dst,
686  const OtherVectorType & src,
687  const bool uses_dof_handler_vector,
688  const typename ::mg::Signals &signals,
689  ...)
690  {
691  (void)uses_dof_handler_vector;
692  Assert(!uses_dof_handler_vector, ExcInternalError());
693 
694  signals.transfer_to_mg(true);
695  transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
696  signals.transfer_to_mg(false);
697 
698  multigrid.cycle();
699 
700  signals.transfer_to_global(true);
701  transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
702  signals.transfer_to_global(false);
703  }
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_add(
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_add(dof_handler_vector, dst, multigrid.solution);
732  else
733  transfer.copy_from_mg_add(*dof_handler_vector[0],
734  dst,
735  multigrid.solution);
736  signals.transfer_to_global(false);
737  }
738 
739  template <int dim,
740  typename VectorType,
741  class TRANSFER,
742  typename OtherVectorType>
743  void
744  vmult_add(
745  const std::vector<const ::DoFHandler<dim> *> &dof_handler_vector,
746  ::Multigrid<VectorType> & multigrid,
747  const TRANSFER & 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_add(*dof_handler_vector[0],
765  dst,
766  multigrid.solution);
767  signals.transfer_to_global(false);
768  }
769  } // namespace PreconditionMGImplementation
770 } // namespace internal
771 
772 template <int dim, typename VectorType, class TRANSFER>
774  const DoFHandler<dim> &dof_handler,
776  const TRANSFER & transfer)
777  : dof_handler_vector(1, &dof_handler)
778  , dof_handler_vector_raw(1, &dof_handler)
779  , multigrid(&mg)
780  , transfer(&transfer)
781  , uses_dof_handler_vector(false)
782 {}
783 
784 template <int dim, typename VectorType, class TRANSFER>
786  const std::vector<const DoFHandler<dim> *> &dof_handler,
788  const TRANSFER & transfer)
789  : dof_handler_vector(dof_handler.size())
790  , dof_handler_vector_raw(dof_handler.size())
791  , multigrid(&mg)
792  , transfer(&transfer)
794 {
795  for (unsigned int i = 0; i < dof_handler.size(); ++i)
796  {
797  dof_handler_vector[i] = dof_handler[i];
798  dof_handler_vector_raw[i] = dof_handler[i];
799  }
800 }
801 
802 template <int dim, typename VectorType, class TRANSFER>
803 inline bool
805 {
806  return false;
807 }
808 
809 template <int dim, typename VectorType, class TRANSFER>
810 template <class OtherVectorType>
811 void
813  OtherVectorType & dst,
814  const OtherVectorType &src) const
815 {
816  internal::PreconditionMGImplementation::vmult(dof_handler_vector_raw,
817  *multigrid,
818  *transfer,
819  dst,
820  src,
822  this->signals,
823  0);
824 }
825 
826 
827 template <int dim, typename VectorType, class TRANSFER>
828 IndexSet
830  const unsigned int block) const
831 {
832  AssertIndexRange(block, dof_handler_vector.size());
833  return dof_handler_vector[block]->locally_owned_dofs();
834 }
835 
836 
837 template <int dim, typename VectorType, class TRANSFER>
838 IndexSet
840  const unsigned int block) const
841 {
842  AssertIndexRange(block, dof_handler_vector.size());
843  return dof_handler_vector[block]->locally_owned_dofs();
844 }
845 
846 
847 
848 template <int dim, typename VectorType, class TRANSFER>
849 const MPI_Comm &
851 {
852  // currently parallel GMG works with parallel triangulations only,
853  // so it should be a safe bet to use it to query MPI communicator:
854  const Triangulation<dim> &tria = dof_handler_vector[0]->get_triangulation();
855  const parallel::TriangulationBase<dim> *ptria =
856  dynamic_cast<const parallel::TriangulationBase<dim> *>(&tria);
857  Assert(ptria != nullptr, ExcInternalError());
858  return ptria->get_communicator();
859 }
860 
861 
862 
863 template <int dim, typename VectorType, class TRANSFER>
864 boost::signals2::connection
866  const std::function<void(bool)> &slot)
867 {
868  return this->signals.transfer_to_mg.connect(slot);
869 }
870 
871 
872 
873 template <int dim, typename VectorType, class TRANSFER>
874 boost::signals2::connection
876  const std::function<void(bool)> &slot)
877 {
878  return this->signals.transfer_to_global.connect(slot);
879 }
880 
881 
882 
883 template <int dim, typename VectorType, class TRANSFER>
884 template <class OtherVectorType>
885 void
887  OtherVectorType & dst,
888  const OtherVectorType &src) const
889 {
890  internal::PreconditionMGImplementation::vmult_add(dof_handler_vector_raw,
891  *multigrid,
892  *transfer,
893  dst,
894  src,
896  this->signals,
897  0);
898 }
899 
900 
901 template <int dim, typename VectorType, class TRANSFER>
902 template <class OtherVectorType>
903 void
905  const OtherVectorType &) const
906 {
907  Assert(false, ExcNotImplemented());
908 }
909 
910 
911 template <int dim, typename VectorType, class TRANSFER>
912 template <class OtherVectorType>
913 void
915  OtherVectorType &,
916  const OtherVectorType &) const
917 {
918  Assert(false, ExcNotImplemented());
919 }
920 
921 #endif // DOXYGEN
922 
924 
925 #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:375
void cycle()
static const unsigned int invalid_unsigned_int
Definition: types.h:191
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_down
Definition: multigrid.h:421
SmartPointer< const MGMatrixBase< VectorType > > edge_in
Definition: multigrid.h:414
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:91
#define AssertIndexRange(index, range)
Definition: exceptions.h:1628
SmartPointer< const MGCoarseGridBase< VectorType >, Multigrid< VectorType > > coarse
Definition: multigrid.h:381
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
boost::signals2::signal< void(const bool before, const unsigned int level)> restriction
Definition: multigrid.h:83
MGLevelObject< VectorType > t
Definition: multigrid.h:364
bool empty() const
const bool uses_dof_handler_vector
Definition: multigrid.h:577
mg::Signals signals
Definition: multigrid.h:312
unsigned int maxlevel
Definition: multigrid.h:346
IndexSet locally_owned_range_indices(const unsigned int block=0) const
#define Assert(cond, exc)
Definition: exceptions.h:1403
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:4355
SmartPointer< Multigrid< VectorType >, PreconditionMG< dim, VectorType, TRANSFER > > multigrid
Definition: multigrid.h:565
boost::signals2::signal< void(const bool before)> transfer_to_mg
Definition: multigrid.h:60
const MPI_Comm & get_mpi_communicator() const
SmartPointer< const MGMatrixBase< VectorType > > edge_out
Definition: multigrid.h:406
void vmult_add(OtherVectorType &dst, const OtherVectorType &src) const
The W-cycle.
Definition: multigrid.h:146
boost::signals2::signal< void(const bool before, const unsigned int level)> pre_smoother_step
Definition: multigrid.h:99
The V-cycle.
Definition: multigrid.h:144
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:341
SmartPointer< const TRANSFER, PreconditionMG< dim, VectorType, TRANSFER > > transfer
Definition: multigrid.h:571
MGLevelObject< VectorType > solution
Definition: multigrid.h:358
std::vector< SmartPointer< const DoFHandler< dim >, PreconditionMG< dim, VectorType, TRANSFER > > > dof_handler_vector
Definition: multigrid.h:553
virtual unsigned int get_maxlevel() const =0
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > post_smooth
Definition: multigrid.h:399
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:559
mg::Signals signals
Definition: multigrid.h:582
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:369
MGLevelObject< VectorType > defect
Definition: multigrid.h:353
PreconditionMG(const DoFHandler< dim > &dof_handler, Multigrid< VectorType > &mg, const TRANSFER &transfer)
SmartPointer< const MGTransferBase< VectorType >, Multigrid< VectorType > > transfer
Definition: multigrid.h:387
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_up
Definition: multigrid.h:428
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > pre_smooth
Definition: multigrid.h:393
Cycle cycle_type
Definition: multigrid.h:336
static ::ExceptionBase & ExcInternalError()