Reference documentation for deal.II version GIT relicensing-214-g6e74dec06b 2024-03-27 18:10:01+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\}}\)
Loading...
Searching...
No Matches
multigrid.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1999 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_multigrid_h
16#define dealii_multigrid_h
17
18
19#include <deal.II/base/config.h>
20
24
26
28
30#include <deal.II/lac/vector.h>
31
33
34#include <vector>
35
37
38#ifdef signals
39# error \
40 "The name 'signals' is already defined. You are most likely using the QT library \
41and using the 'signals' keyword. You can either #include the Qt headers (or any conflicting headers) \
42*after* the deal.II headers or you can define the 'QT_NO_KEYWORDS' macro and use the 'Q_SIGNALS' macro."
43#endif
44
50namespace mg
51{
61 struct Signals
62 {
68 boost::signals2::signal<void(const bool before)> transfer_to_mg;
69
75 boost::signals2::signal<void(const bool before)> transfer_to_global;
76
86 boost::signals2::signal<void(const bool before, const unsigned int level)>
88
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
116 boost::signals2::signal<void(const bool before, const unsigned int level)>
118
124 boost::signals2::signal<void(const bool before, const unsigned int level)>
126
132 boost::signals2::signal<void(const bool before, const unsigned int level)>
134
139 boost::signals2::signal<void(const bool before, const unsigned int level)>
141 };
142} // namespace mg
143
161template <typename VectorType>
162class Multigrid : public Subscriptor
163{
164public:
177
178 using vector_type = VectorType;
179 using const_vector_type = const VectorType;
180
195 const unsigned int minlevel = 0,
196 const unsigned int maxlevel = numbers::invalid_unsigned_int,
198
202 void
203 reinit(const unsigned int minlevel, const unsigned int maxlevel);
204
209 void
211
222 void
224
237 void
240
247 void
249
262 void
265
269 unsigned int
271
275 unsigned int
277
283 void
284 set_maxlevel(const unsigned int);
285
301 void
302 set_minlevel(const unsigned int level, bool relative = false);
303
308
312 boost::signals2::connection
314 const std::function<void(const bool, const unsigned int)> &slot);
315
319 boost::signals2::connection
321 const std::function<void(const bool, const unsigned int)> &slot);
322
326 boost::signals2::connection
328 const std::function<void(const bool, const unsigned int)> &slot);
329
333 boost::signals2::connection
335 const std::function<void(const bool, const unsigned int)> &slot);
336
340 boost::signals2::connection
342 const std::function<void(const bool, const unsigned int)> &slot);
343
347 boost::signals2::connection
349 const std::function<void(const bool, const unsigned int)> &slot);
350
354 boost::signals2::connection
356 const std::function<void(const bool, const unsigned int)> &slot);
357
358private:
363
370 void
371 level_v_step(const unsigned int level);
372
380 void
381 level_step(const unsigned int level, Cycle cycle);
382
387
391 unsigned int minlevel;
392
396 unsigned int maxlevel;
397
398public:
404
409
410private:
415
420
421
426
432
438
444
450
457
465
472
479
480 template <int dim, typename OtherVectorType, typename TransferType>
481 friend class PreconditionMG;
482};
483
484
499template <int dim, typename VectorType, typename TransferType>
501{
502public:
507 PreconditionMG(const DoFHandler<dim> &dof_handler,
509 const TransferType &transfer);
510
515 PreconditionMG(const std::vector<const DoFHandler<dim> *> &dof_handler,
517 const TransferType &transfer);
518
522 bool
523 empty() const;
524
531 template <typename OtherVectorType>
532 void
533 vmult(OtherVectorType &dst, const OtherVectorType &src) const;
534
539 template <typename OtherVectorType>
540 void
541 vmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
542
548 template <typename OtherVectorType>
549 void
550 Tvmult(OtherVectorType &dst, const OtherVectorType &src) const;
551
557 template <typename OtherVectorType>
558 void
559 Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
560
568 locally_owned_range_indices(const unsigned int block = 0) const;
569
577 locally_owned_domain_indices(const unsigned int block = 0) const;
578
584
588 boost::signals2::connection
589 connect_transfer_to_mg(const std::function<void(bool)> &slot);
590
594 boost::signals2::connection
595 connect_transfer_to_global(const std::function<void(bool)> &slot);
596
602
608
609private:
613 std::vector<SmartPointer<const DoFHandler<dim>,
616
621 std::vector<const DoFHandler<dim> *> dof_handler_vector_raw;
622
629
633 SmartPointer<const TransferType,
636
642
647};
648
651#ifndef DOXYGEN
652/* --------------------------- inline functions --------------------- */
653
654
655template <typename VectorType>
657 const MGCoarseGridBase<VectorType> &coarse,
658 const MGTransferBase<VectorType> &transfer,
659 const MGSmootherBase<VectorType> &pre_smooth,
660 const MGSmootherBase<VectorType> &post_smooth,
661 const unsigned int min_level,
662 const unsigned int max_level,
663 Cycle cycle)
664 : cycle_type(cycle)
665 , matrix(&matrix, typeid(*this).name())
666 , coarse(&coarse, typeid(*this).name())
667 , transfer(&transfer, typeid(*this).name())
668 , pre_smooth(&pre_smooth, typeid(*this).name())
669 , post_smooth(&post_smooth, typeid(*this).name())
670 , edge_out(nullptr, typeid(*this).name())
671 , edge_in(nullptr, typeid(*this).name())
672 , edge_down(nullptr, typeid(*this).name())
673 , edge_up(nullptr, typeid(*this).name())
674{
675 if (max_level == numbers::invalid_unsigned_int)
676 maxlevel = matrix.get_maxlevel();
677 else
678 maxlevel = max_level;
679 reinit(min_level, maxlevel);
680}
681
682
683
684template <typename VectorType>
685inline unsigned int
687{
688 return maxlevel;
689}
690
691
692
693template <typename VectorType>
694inline unsigned int
696{
697 return minlevel;
698}
699
700
701/* --------------------------- inline functions --------------------- */
702
703
704namespace internal
705{
706 namespace PreconditionMGImplementation
707 {
708 template <int dim,
709 typename VectorType,
710 typename TransferType,
711 typename OtherVectorType>
712 std::enable_if_t<TransferType::supports_dof_handler_vector>
713 vmult(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
714 Multigrid<VectorType> &multigrid,
715 const TransferType &transfer,
716 OtherVectorType &dst,
717 const OtherVectorType &src,
718 const bool uses_dof_handler_vector,
719 const typename ::mg::Signals &signals,
720 int)
721 {
722 signals.transfer_to_mg(true);
723 if (uses_dof_handler_vector)
724 transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
725 else
726 transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
727 signals.transfer_to_mg(false);
728
729 multigrid.cycle();
730
731 signals.transfer_to_global(true);
732 if (uses_dof_handler_vector)
733 transfer.copy_from_mg(dof_handler_vector, dst, multigrid.solution);
734 else
735 transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
736 signals.transfer_to_global(false);
737 }
738
739 template <int dim,
740 typename VectorType,
741 typename TransferType,
742 typename OtherVectorType>
743 void
744 vmult(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
745 Multigrid<VectorType> &multigrid,
746 const TransferType &transfer,
747 OtherVectorType &dst,
748 const OtherVectorType &src,
749 const bool uses_dof_handler_vector,
750 const typename ::mg::Signals &signals,
751 ...)
752 {
753 (void)uses_dof_handler_vector;
754 Assert(!uses_dof_handler_vector, ExcInternalError());
755
756 signals.transfer_to_mg(true);
757 transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
758 signals.transfer_to_mg(false);
759
760 multigrid.cycle();
761
762 signals.transfer_to_global(true);
763 transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
764 signals.transfer_to_global(false);
765 }
766
767 template <int dim,
768 typename VectorType,
769 typename TransferType,
770 typename OtherVectorType>
771 std::enable_if_t<TransferType::supports_dof_handler_vector>
772 vmult_add(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
773 Multigrid<VectorType> &multigrid,
774 const TransferType &transfer,
775 OtherVectorType &dst,
776 const OtherVectorType &src,
777 const bool uses_dof_handler_vector,
778 const typename ::mg::Signals &signals,
779 int)
780 {
781 signals.transfer_to_mg(true);
782 if (uses_dof_handler_vector)
783 transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
784 else
785 transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
786 signals.transfer_to_mg(false);
787
788 multigrid.cycle();
789
790 signals.transfer_to_global(true);
791 if (uses_dof_handler_vector)
792 transfer.copy_from_mg_add(dof_handler_vector, dst, multigrid.solution);
793 else
794 transfer.copy_from_mg_add(*dof_handler_vector[0],
795 dst,
796 multigrid.solution);
797 signals.transfer_to_global(false);
798 }
799
800 template <int dim,
801 typename VectorType,
802 typename TransferType,
803 typename OtherVectorType>
804 void
805 vmult_add(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
806 Multigrid<VectorType> &multigrid,
807 const TransferType &transfer,
808 OtherVectorType &dst,
809 const OtherVectorType &src,
810 const bool uses_dof_handler_vector,
811 const typename ::mg::Signals &signals,
812 ...)
813 {
814 (void)uses_dof_handler_vector;
815 Assert(!uses_dof_handler_vector, ExcInternalError());
816
817 signals.transfer_to_mg(true);
818 transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
819 signals.transfer_to_mg(false);
820
821 multigrid.cycle();
822
823 signals.transfer_to_global(true);
824 transfer.copy_from_mg_add(*dof_handler_vector[0],
825 dst,
826 multigrid.solution);
827 signals.transfer_to_global(false);
828 }
829 } // namespace PreconditionMGImplementation
830} // namespace internal
831
832template <int dim, typename VectorType, typename TransferType>
834 const DoFHandler<dim> &dof_handler,
836 const TransferType &transfer)
837 : dof_handler_vector(1, &dof_handler)
838 , dof_handler_vector_raw(1, &dof_handler)
839 , multigrid(&mg)
841 , uses_dof_handler_vector(false)
842{}
843
844template <int dim, typename VectorType, typename TransferType>
846 const std::vector<const DoFHandler<dim> *> &dof_handler,
848 const TransferType &transfer)
849 : dof_handler_vector(dof_handler.size())
850 , dof_handler_vector_raw(dof_handler.size())
851 , multigrid(&mg)
853 , uses_dof_handler_vector(true)
854{
855 for (unsigned int i = 0; i < dof_handler.size(); ++i)
856 {
857 dof_handler_vector[i] = dof_handler[i];
858 dof_handler_vector_raw[i] = dof_handler[i];
859 }
860}
861
862template <int dim, typename VectorType, typename TransferType>
863inline bool
865{
866 return false;
867}
868
869template <int dim, typename VectorType, typename TransferType>
870template <typename OtherVectorType>
871void
873 OtherVectorType &dst,
874 const OtherVectorType &src) const
875{
876 internal::PreconditionMGImplementation::vmult(dof_handler_vector_raw,
877 *multigrid,
878 *transfer,
879 dst,
880 src,
881 uses_dof_handler_vector,
882 this->signals,
883 0);
884}
885
886
887template <int dim, typename VectorType, typename TransferType>
890 const unsigned int block) const
891{
892 AssertIndexRange(block, dof_handler_vector.size());
893 return dof_handler_vector[block]->locally_owned_dofs();
894}
895
896
897template <int dim, typename VectorType, typename TransferType>
900 const unsigned int block) const
901{
902 AssertIndexRange(block, dof_handler_vector.size());
903 return dof_handler_vector[block]->locally_owned_dofs();
904}
905
906
907
908template <int dim, typename VectorType, typename TransferType>
911{
912 // currently parallel GMG works with parallel triangulations only,
913 // so it should be a safe bet to use it to query MPI communicator:
914 const Triangulation<dim> &tria = dof_handler_vector[0]->get_triangulation();
916 dynamic_cast<const parallel::TriangulationBase<dim> *>(&tria);
917 Assert(ptria != nullptr, ExcInternalError());
918 return ptria->get_communicator();
919}
920
921
922
923template <int dim, typename VectorType, typename TransferType>
924boost::signals2::connection
926 const std::function<void(bool)> &slot)
927{
928 return this->signals.transfer_to_mg.connect(slot);
929}
930
931
932
933template <int dim, typename VectorType, typename TransferType>
934boost::signals2::connection
936 const std::function<void(bool)> &slot)
937{
938 return this->signals.transfer_to_global.connect(slot);
939}
940
941
942
943template <int dim, typename VectorType, typename TransferType>
944template <typename OtherVectorType>
945void
947 OtherVectorType &dst,
948 const OtherVectorType &src) const
949{
950 internal::PreconditionMGImplementation::vmult_add(dof_handler_vector_raw,
951 *multigrid,
952 *transfer,
953 dst,
954 src,
955 uses_dof_handler_vector,
956 this->signals,
957 0);
958}
959
960
961template <int dim, typename VectorType, typename TransferType>
962template <typename OtherVectorType>
963void
965 OtherVectorType &,
966 const OtherVectorType &) const
967{
969}
970
971
972template <int dim, typename VectorType, typename TransferType>
973template <typename OtherVectorType>
974void
976 OtherVectorType &,
977 const OtherVectorType &) const
978{
980}
981
982
983template <int dim, typename VectorType, typename TransferType>
986{
987 return *this->multigrid;
988}
989
990
991template <int dim, typename VectorType, typename TransferType>
994{
995 return *this->multigrid;
996}
997
998#endif // DOXYGEN
999
1001
1002#endif
@ v_cycle
The V-cycle.
Definition multigrid.h:171
@ w_cycle
The W-cycle.
Definition multigrid.h:173
@ f_cycle
The F-cycle.
Definition multigrid.h:175
MGLevelObject< VectorType > defect2
Definition multigrid.h:419
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:437
unsigned int minlevel
Definition multigrid.h:391
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_up
Definition multigrid.h:478
mg::Signals signals
Definition multigrid.h:362
SmartPointer< const MGMatrixBase< VectorType > > edge_in
Definition multigrid.h:464
SmartPointer< const MGMatrixBase< VectorType > > edge_out
Definition multigrid.h:456
MGLevelObject< VectorType > solution
Definition multigrid.h:408
boost::signals2::connection connect_edge_prolongation(const std::function< void(const bool, const unsigned int)> &slot)
unsigned int maxlevel
Definition multigrid.h:396
VectorType vector_type
Definition multigrid.h:178
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > matrix
Definition multigrid.h:425
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:179
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > pre_smooth
Definition multigrid.h:443
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > post_smooth
Definition multigrid.h:449
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:431
Cycle cycle_type
Definition multigrid.h:386
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:471
void set_edge_in_matrix(const MGMatrixBase< VectorType > &edge_in)
MGLevelObject< VectorType > t
Definition multigrid.h:414
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:403
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:621
void Tvmult(OtherVectorType &dst, const OtherVectorType &src) const
SmartPointer< const TransferType, PreconditionMG< dim, VectorType, TransferType > > transfer
Definition multigrid.h:635
IndexSet locally_owned_domain_indices(const unsigned int block=0) const
const Multigrid< VectorType > & get_multigrid() const
boost::signals2::connection connect_transfer_to_global(const std::function< void(bool)> &slot)
mg::Signals signals
Definition multigrid.h:646
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:641
void vmult_add(OtherVectorType &dst, const OtherVectorType &src) const
bool empty() const
IndexSet locally_owned_range_indices(const unsigned int block=0) const
SmartPointer< Multigrid< VectorType >, PreconditionMG< dim, VectorType, TransferType > > multigrid
Definition multigrid.h:628
PreconditionMG(const std::vector< const DoFHandler< dim > * > &dof_handler, Multigrid< VectorType > &mg, const TransferType &transfer)
Multigrid< VectorType > & get_multigrid()
std::vector< SmartPointer< const DoFHandler< dim >, PreconditionMG< dim, VectorType, TransferType > > > dof_handler_vector
Definition multigrid.h:615
void Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const
void vmult(OtherVectorType &dst, const OtherVectorType &src) const
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:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
unsigned int level
Definition grid_out.cc:4616
#define Assert(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
#define DEAL_II_NOT_IMPLEMENTED()
const Triangulation< dim, spacedim > & tria
spacedim const DoFHandler< dim, spacedim > const FullMatrix< double > & transfer
Definition mg.h:81
static const unsigned int invalid_unsigned_int
Definition types.h:220
boost::signals2::signal< void(const bool before, const unsigned int level)> prolongation
Definition multigrid.h:106
boost::signals2::signal< void(const bool before, const unsigned int level)> pre_smoother_step
Definition multigrid.h:117
boost::signals2::signal< void(const bool before, const unsigned int level)> coarse_solve
Definition multigrid.h:87
boost::signals2::signal< void(const bool before, const unsigned int level)> edge_prolongation
Definition multigrid.h:140
boost::signals2::signal< void(const bool before, const unsigned int level)> post_smoother_step
Definition multigrid.h:125
boost::signals2::signal< void(const bool before)> transfer_to_global
Definition multigrid.h:75
boost::signals2::signal< void(const bool before, const unsigned int level)> residual_step
Definition multigrid.h:133
boost::signals2::signal< void(const bool before, const unsigned int level)> restriction
Definition multigrid.h:98
boost::signals2::signal< void(const bool before)> transfer_to_mg
Definition multigrid.h:68