deal.II version GIT relicensing-2642-g417c8935ea 2025-02-14 17:30:00+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 <boost/signals2.hpp>
35
36#include <vector>
37
38
40
41#ifdef signals
42# error \
43 "The name 'signals' is already defined. You are most likely using the QT library \
44and using the 'signals' keyword. You can either #include the Qt headers (or any conflicting headers) \
45*after* the deal.II headers or you can define the 'QT_NO_KEYWORDS' macro and use the 'Q_SIGNALS' macro."
46#endif
47
53namespace mg
54{
64 struct Signals
65 {
71 boost::signals2::signal<void(const bool before)> transfer_to_mg;
72
78 boost::signals2::signal<void(const bool before)> transfer_to_global;
79
89 boost::signals2::signal<void(const bool before, const unsigned int level)>
91
100 boost::signals2::signal<void(const bool before, const unsigned int level)>
102
108 boost::signals2::signal<void(const bool before, const unsigned int level)>
110
119 boost::signals2::signal<void(const bool before, const unsigned int level)>
121
127 boost::signals2::signal<void(const bool before, const unsigned int level)>
129
135 boost::signals2::signal<void(const bool before, const unsigned int level)>
137
142 boost::signals2::signal<void(const bool before, const unsigned int level)>
144 };
145} // namespace mg
146
164template <typename VectorType>
166{
167public:
180
181 using vector_type = VectorType;
182 using const_vector_type = const VectorType;
183
198 const unsigned int minlevel = 0,
199 const unsigned int maxlevel = numbers::invalid_unsigned_int,
201
205 void
206 reinit(const unsigned int minlevel, const unsigned int maxlevel);
207
212 void
214
225 void
227
240 void
243
250 void
252
265 void
268
272 unsigned int
274
278 unsigned int
280
286 void
287 set_maxlevel(const unsigned int);
288
304 void
305 set_minlevel(const unsigned int level, bool relative = false);
306
311
315 boost::signals2::connection
317 const std::function<void(const bool, const unsigned int)> &slot);
318
322 boost::signals2::connection
324 const std::function<void(const bool, const unsigned int)> &slot);
325
329 boost::signals2::connection
331 const std::function<void(const bool, const unsigned int)> &slot);
332
336 boost::signals2::connection
338 const std::function<void(const bool, const unsigned int)> &slot);
339
343 boost::signals2::connection
345 const std::function<void(const bool, const unsigned int)> &slot);
346
350 boost::signals2::connection
352 const std::function<void(const bool, const unsigned int)> &slot);
353
357 boost::signals2::connection
359 const std::function<void(const bool, const unsigned int)> &slot);
360
361private:
366
373 void
374 level_v_step(const unsigned int level);
375
383 void
384 level_step(const unsigned int level, Cycle cycle);
385
390
394 unsigned int minlevel;
395
399 unsigned int maxlevel;
400
401public:
407
412
413private:
418
423
424
429
435
441
447
453
460
468
476
484
485 template <int dim, typename OtherVectorType, typename TransferType>
486 friend class PreconditionMG;
487};
488
489
504template <int dim, typename VectorType, typename TransferType>
506{
507public:
512 PreconditionMG(const DoFHandler<dim> &dof_handler,
514 const TransferType &transfer);
515
520 PreconditionMG(const std::vector<const DoFHandler<dim> *> &dof_handler,
522 const TransferType &transfer);
523
527 bool
528 empty() const;
529
536 template <typename OtherVectorType>
537 void
538 vmult(OtherVectorType &dst, const OtherVectorType &src) const;
539
544 template <typename OtherVectorType>
545 void
546 vmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
547
553 template <typename OtherVectorType>
554 void
555 Tvmult(OtherVectorType &dst, const OtherVectorType &src) const;
556
562 template <typename OtherVectorType>
563 void
564 Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const;
565
573 locally_owned_range_indices(const unsigned int block = 0) const;
574
582 locally_owned_domain_indices(const unsigned int block = 0) const;
583
589
593 boost::signals2::connection
594 connect_transfer_to_mg(const std::function<void(bool)> &slot);
595
599 boost::signals2::connection
600 connect_transfer_to_global(const std::function<void(bool)> &slot);
601
607
613
614private:
618 std::vector<ObserverPointer<const DoFHandler<dim>,
621
626 std::vector<const DoFHandler<dim> *> dof_handler_vector_raw;
627
634
638 ObserverPointer<const TransferType,
641
647
652};
653
656#ifndef DOXYGEN
657/* --------------------------- inline functions --------------------- */
658
659
660template <typename VectorType>
662 const MGCoarseGridBase<VectorType> &coarse,
663 const MGTransferBase<VectorType> &transfer,
664 const MGSmootherBase<VectorType> &pre_smooth,
665 const MGSmootherBase<VectorType> &post_smooth,
666 const unsigned int min_level,
667 const unsigned int max_level,
668 Cycle cycle)
669 : cycle_type(cycle)
670 , matrix(&matrix, typeid(*this).name())
671 , coarse(&coarse, typeid(*this).name())
672 , transfer(&transfer, typeid(*this).name())
673 , pre_smooth(&pre_smooth, typeid(*this).name())
674 , post_smooth(&post_smooth, typeid(*this).name())
675 , edge_out(nullptr, typeid(*this).name())
676 , edge_in(nullptr, typeid(*this).name())
677 , edge_down(nullptr, typeid(*this).name())
678 , edge_up(nullptr, typeid(*this).name())
679{
680 if (max_level == numbers::invalid_unsigned_int)
681 maxlevel = matrix.get_maxlevel();
682 else
683 maxlevel = max_level;
684 reinit(min_level, maxlevel);
685}
686
687
688
689template <typename VectorType>
690inline unsigned int
692{
693 return maxlevel;
694}
695
696
697
698template <typename VectorType>
699inline unsigned int
701{
702 return minlevel;
703}
704
705
706/* --------------------------- inline functions --------------------- */
707
708
709namespace internal
710{
711 namespace PreconditionMGImplementation
712 {
713 template <int dim,
714 typename VectorType,
715 typename TransferType,
716 typename OtherVectorType>
717 std::enable_if_t<TransferType::supports_dof_handler_vector>
718 vmult(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
719 Multigrid<VectorType> &multigrid,
720 const TransferType &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(dof_handler_vector, dst, multigrid.solution);
739 else
740 transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
741 signals.transfer_to_global(false);
742 }
743
744 template <int dim,
745 typename VectorType,
746 typename TransferType,
747 typename OtherVectorType>
748 void
749 vmult(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
750 Multigrid<VectorType> &multigrid,
751 const TransferType &transfer,
752 OtherVectorType &dst,
753 const OtherVectorType &src,
754 const bool uses_dof_handler_vector,
755 const typename ::mg::Signals &signals,
756 ...)
757 {
758 (void)uses_dof_handler_vector;
759 Assert(!uses_dof_handler_vector, ExcInternalError());
760
761 signals.transfer_to_mg(true);
762 transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
763 signals.transfer_to_mg(false);
764
765 multigrid.cycle();
766
767 signals.transfer_to_global(true);
768 transfer.copy_from_mg(*dof_handler_vector[0], dst, multigrid.solution);
769 signals.transfer_to_global(false);
770 }
771
772 template <int dim,
773 typename VectorType,
774 typename TransferType,
775 typename OtherVectorType>
776 std::enable_if_t<TransferType::supports_dof_handler_vector>
777 vmult_add(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
778 Multigrid<VectorType> &multigrid,
779 const TransferType &transfer,
780 OtherVectorType &dst,
781 const OtherVectorType &src,
782 const bool uses_dof_handler_vector,
783 const typename ::mg::Signals &signals,
784 int)
785 {
786 signals.transfer_to_mg(true);
787 if (uses_dof_handler_vector)
788 transfer.copy_to_mg(dof_handler_vector, multigrid.defect, src);
789 else
790 transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
791 signals.transfer_to_mg(false);
792
793 multigrid.cycle();
794
795 signals.transfer_to_global(true);
796 if (uses_dof_handler_vector)
797 transfer.copy_from_mg_add(dof_handler_vector, dst, multigrid.solution);
798 else
799 transfer.copy_from_mg_add(*dof_handler_vector[0],
800 dst,
801 multigrid.solution);
802 signals.transfer_to_global(false);
803 }
804
805 template <int dim,
806 typename VectorType,
807 typename TransferType,
808 typename OtherVectorType>
809 void
810 vmult_add(const std::vector<const DoFHandler<dim> *> &dof_handler_vector,
811 Multigrid<VectorType> &multigrid,
812 const TransferType &transfer,
813 OtherVectorType &dst,
814 const OtherVectorType &src,
815 const bool uses_dof_handler_vector,
816 const typename ::mg::Signals &signals,
817 ...)
818 {
819 (void)uses_dof_handler_vector;
820 Assert(!uses_dof_handler_vector, ExcInternalError());
821
822 signals.transfer_to_mg(true);
823 transfer.copy_to_mg(*dof_handler_vector[0], multigrid.defect, src);
824 signals.transfer_to_mg(false);
825
826 multigrid.cycle();
827
828 signals.transfer_to_global(true);
829 transfer.copy_from_mg_add(*dof_handler_vector[0],
830 dst,
831 multigrid.solution);
832 signals.transfer_to_global(false);
833 }
834 } // namespace PreconditionMGImplementation
835} // namespace internal
836
837template <int dim, typename VectorType, typename TransferType>
839 const DoFHandler<dim> &dof_handler,
841 const TransferType &transfer)
842 : dof_handler_vector(1, &dof_handler)
843 , dof_handler_vector_raw(1, &dof_handler)
844 , multigrid(&mg)
845 , transfer(&transfer)
846 , uses_dof_handler_vector(false)
847{}
848
849template <int dim, typename VectorType, typename TransferType>
851 const std::vector<const DoFHandler<dim> *> &dof_handler,
853 const TransferType &transfer)
854 : dof_handler_vector(dof_handler.size())
855 , dof_handler_vector_raw(dof_handler.size())
856 , multigrid(&mg)
857 , transfer(&transfer)
858 , uses_dof_handler_vector(true)
859{
860 for (unsigned int i = 0; i < dof_handler.size(); ++i)
861 {
862 dof_handler_vector[i] = dof_handler[i];
863 dof_handler_vector_raw[i] = dof_handler[i];
864 }
865}
866
867template <int dim, typename VectorType, typename TransferType>
868inline bool
870{
871 return false;
872}
873
874template <int dim, typename VectorType, typename TransferType>
875template <typename OtherVectorType>
876void
878 OtherVectorType &dst,
879 const OtherVectorType &src) const
880{
881 internal::PreconditionMGImplementation::vmult(dof_handler_vector_raw,
882 *multigrid,
883 *transfer,
884 dst,
885 src,
886 uses_dof_handler_vector,
887 this->signals,
888 0);
889}
890
891
892template <int dim, typename VectorType, typename TransferType>
895 const unsigned int block) const
896{
897 AssertIndexRange(block, dof_handler_vector.size());
898 return dof_handler_vector[block]->locally_owned_dofs();
899}
900
901
902template <int dim, typename VectorType, typename TransferType>
905 const unsigned int block) const
906{
907 AssertIndexRange(block, dof_handler_vector.size());
908 return dof_handler_vector[block]->locally_owned_dofs();
909}
910
911
912
913template <int dim, typename VectorType, typename TransferType>
916{
917 // currently parallel GMG works with parallel triangulations only,
918 // so it should be a safe bet to use it to query MPI communicator:
919 const Triangulation<dim> &tria = dof_handler_vector[0]->get_triangulation();
921 dynamic_cast<const parallel::TriangulationBase<dim> *>(&tria);
922 Assert(ptria != nullptr, ExcInternalError());
923 return ptria->get_mpi_communicator();
924}
925
926
927
928template <int dim, typename VectorType, typename TransferType>
929boost::signals2::connection
931 const std::function<void(bool)> &slot)
932{
933 return this->signals.transfer_to_mg.connect(slot);
934}
935
936
937
938template <int dim, typename VectorType, typename TransferType>
939boost::signals2::connection
941 const std::function<void(bool)> &slot)
942{
943 return this->signals.transfer_to_global.connect(slot);
944}
945
946
947
948template <int dim, typename VectorType, typename TransferType>
949template <typename OtherVectorType>
950void
952 OtherVectorType &dst,
953 const OtherVectorType &src) const
954{
955 internal::PreconditionMGImplementation::vmult_add(dof_handler_vector_raw,
956 *multigrid,
957 *transfer,
958 dst,
959 src,
960 uses_dof_handler_vector,
961 this->signals,
962 0);
963}
964
965
966template <int dim, typename VectorType, typename TransferType>
967template <typename OtherVectorType>
968void
970 OtherVectorType &,
971 const OtherVectorType &) const
972{
974}
975
976
977template <int dim, typename VectorType, typename TransferType>
978template <typename OtherVectorType>
979void
981 OtherVectorType &,
982 const OtherVectorType &) const
983{
985}
986
987
988template <int dim, typename VectorType, typename TransferType>
991{
992 return *this->multigrid;
993}
994
995
996template <int dim, typename VectorType, typename TransferType>
999{
1000 return *this->multigrid;
1001}
1002
1003#endif // DOXYGEN
1004
1006
1007#endif
@ v_cycle
The V-cycle.
Definition multigrid.h:174
@ w_cycle
The W-cycle.
Definition multigrid.h:176
@ f_cycle
The F-cycle.
Definition multigrid.h:178
ObserverPointer< const MGMatrixBase< VectorType > > edge_out
Definition multigrid.h:459
MGLevelObject< VectorType > defect2
Definition multigrid.h:422
void cycle()
void set_edge_flux_matrices(const MGMatrixBase< VectorType > &edge_down, const MGMatrixBase< VectorType > &edge_up)
unsigned int minlevel
Definition multigrid.h:394
mg::Signals signals
Definition multigrid.h:365
MGLevelObject< VectorType > solution
Definition multigrid.h:411
ObserverPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > post_smooth
Definition multigrid.h:452
boost::signals2::connection connect_edge_prolongation(const std::function< void(const bool, const unsigned int)> &slot)
ObserverPointer< const MGMatrixBase< VectorType > > edge_in
Definition multigrid.h:467
unsigned int maxlevel
Definition multigrid.h:399
VectorType vector_type
Definition multigrid.h:181
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
ObserverPointer< const MGCoarseGridBase< VectorType >, Multigrid< VectorType > > coarse
Definition multigrid.h:434
ObserverPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > pre_smooth
Definition multigrid.h:446
const VectorType const_vector_type
Definition multigrid.h:182
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)
ObserverPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > matrix
Definition multigrid.h:428
Cycle cycle_type
Definition multigrid.h:389
boost::signals2::connection connect_prolongation(const std::function< void(const bool, const unsigned int)> &slot)
void set_edge_in_matrix(const MGMatrixBase< VectorType > &edge_in)
MGLevelObject< VectorType > t
Definition multigrid.h:417
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:406
void vcycle()
ObserverPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_down
Definition multigrid.h:475
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)
ObserverPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_up
Definition multigrid.h:483
boost::signals2::connection connect_pre_smoother_step(const std::function< void(const bool, const unsigned int)> &slot)
ObserverPointer< const MGTransferBase< VectorType >, Multigrid< VectorType > > transfer
Definition multigrid.h:440
std::vector< const DoFHandler< dim > * > dof_handler_vector_raw
Definition multigrid.h:626
void Tvmult(OtherVectorType &dst, const OtherVectorType &src) const
IndexSet locally_owned_domain_indices(const unsigned int block=0) const
ObserverPointer< Multigrid< VectorType >, PreconditionMG< dim, VectorType, TransferType > > multigrid
Definition multigrid.h:633
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:651
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)
ObserverPointer< const TransferType, PreconditionMG< dim, VectorType, TransferType > > transfer
Definition multigrid.h:640
MPI_Comm get_mpi_communicator() const
const bool uses_dof_handler_vector
Definition multigrid.h:646
void vmult_add(OtherVectorType &dst, const OtherVectorType &src) const
bool empty() const
IndexSet locally_owned_range_indices(const unsigned int block=0) const
std::vector< ObserverPointer< const DoFHandler< dim >, PreconditionMG< dim, VectorType, TransferType > > > dof_handler_vector
Definition multigrid.h:620
PreconditionMG(const std::vector< const DoFHandler< dim > * > &dof_handler, Multigrid< VectorType > &mg, const TransferType &transfer)
Multigrid< VectorType > & get_multigrid()
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_mpi_communicator() const override
Definition tria_base.cc:160
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:518
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:519
unsigned int level
Definition grid_out.cc:4632
#define Assert(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
#define DEAL_II_NOT_IMPLEMENTED()
std::size_t size
Definition mpi.cc:739
Tpetra::Vector< Number, LO, GO, NodeType< MemorySpace > > VectorType
Definition mg.h:81
constexpr unsigned int invalid_unsigned_int
Definition types.h:232
boost::signals2::signal< void(const bool before, const unsigned int level)> prolongation
Definition multigrid.h:109
boost::signals2::signal< void(const bool before, const unsigned int level)> pre_smoother_step
Definition multigrid.h:120
boost::signals2::signal< void(const bool before, const unsigned int level)> coarse_solve
Definition multigrid.h:90
boost::signals2::signal< void(const bool before, const unsigned int level)> edge_prolongation
Definition multigrid.h:143
boost::signals2::signal< void(const bool before, const unsigned int level)> post_smoother_step
Definition multigrid.h:128
boost::signals2::signal< void(const bool before)> transfer_to_global
Definition multigrid.h:78
boost::signals2::signal< void(const bool before, const unsigned int level)> residual_step
Definition multigrid.h:136
boost::signals2::signal< void(const bool before, const unsigned int level)> restriction
Definition multigrid.h:101
boost::signals2::signal< void(const bool before)> transfer_to_mg
Definition multigrid.h:71