16 #ifndef dealii_multigrid_h 17 #define dealii_multigrid_h 20 #include <deal.II/base/config.h> 22 #include <deal.II/base/mg_level_object.h> 23 #include <deal.II/base/smartpointer.h> 24 #include <deal.II/base/subscriptor.h> 26 #include <deal.II/distributed/tria.h> 28 #include <deal.II/dofs/dof_handler.h> 30 #include <deal.II/lac/sparse_matrix.h> 31 #include <deal.II/lac/vector.h> 33 #include <deal.II/multigrid/mg_base.h> 37 DEAL_II_NAMESPACE_OPEN
74 boost::signals2::signal<void(const bool before, const unsigned int level)>
82 boost::signals2::signal<void(const bool before, const unsigned int level)>
90 boost::signals2::signal<void(const bool before, const unsigned int level)>
98 boost::signals2::signal<void(const bool before, const unsigned int level)>
107 boost::signals2::signal<void(const bool before, const unsigned int level)>
136 template <
typename VectorType>
153 using vector_type = VectorType;
154 using const_vector_type =
const VectorType;
181 const unsigned int minlevel = 0,
183 Cycle cycle = v_cycle);
199 const unsigned int minlevel = 0,
201 Cycle cycle = v_cycle);
207 reinit(
const unsigned int minlevel,
const unsigned int maxlevel);
265 get_maxlevel()
const;
271 get_minlevel()
const;
279 set_maxlevel(
const unsigned int);
297 set_minlevel(
const unsigned int level,
bool relative =
false);
302 void set_cycle(
Cycle);
312 set_debug(
const unsigned int);
317 boost::signals2::connection
318 connect_coarse_solve(
319 const std::function<
void(
const bool,
const unsigned int)> &slot);
324 boost::signals2::connection
326 const std::function<
void(
const bool,
const unsigned int)> &slot);
331 boost::signals2::connection
332 connect_prolongation(
333 const std::function<
void(
const bool,
const unsigned int)> &slot);
338 boost::signals2::connection
339 connect_pre_smoother_step(
340 const std::function<
void(
const bool,
const unsigned int)> &slot);
345 boost::signals2::connection
346 connect_post_smoother_step(
347 const std::function<
void(
const bool,
const unsigned int)> &slot);
362 level_v_step(
const unsigned int level);
372 level_step(
const unsigned int level,
Cycle cycle);
476 template <
int dim,
class OtherVectorType,
class TRANSFER>
497 template <
int dim,
typename VectorType,
class TRANSFER>
507 const TRANSFER & transfer);
515 const TRANSFER & transfer);
529 template <
class OtherVectorType>
531 vmult(OtherVectorType &dst,
const OtherVectorType &src)
const;
537 template <
class OtherVectorType>
539 vmult_add(OtherVectorType &dst,
const OtherVectorType &src)
const;
546 template <
class OtherVectorType>
548 Tvmult(OtherVectorType &dst,
const OtherVectorType &src)
const;
555 template <
class OtherVectorType>
557 Tvmult_add(OtherVectorType &dst,
const OtherVectorType &src)
const;
566 locally_owned_range_indices(
const unsigned int block = 0)
const;
575 locally_owned_domain_indices(
const unsigned int block = 0)
const;
581 get_mpi_communicator()
const;
586 boost::signals2::connection
587 connect_transfer_to_mg(
const std::function<
void(
bool)> &slot);
592 boost::signals2::connection
593 connect_transfer_to_global(
const std::function<
void(
bool)> &slot);
599 std::vector<SmartPointer<const DoFHandler<dim>,
639 template <
typename VectorType>
647 const unsigned int min_level,
648 const unsigned int max_level,
651 , minlevel(min_level)
652 , matrix(&matrix, typeid(*this).name())
653 , coarse(&coarse, typeid(*this).name())
654 , transfer(&transfer, typeid(*this).name())
655 , pre_smooth(&pre_smooth, typeid(*this).name())
656 , post_smooth(&post_smooth, typeid(*this).name())
657 , edge_down(nullptr, typeid(*this).name())
658 , edge_up(nullptr, typeid(*this).name())
661 const unsigned int dof_handler_max_level =
664 maxlevel = dof_handler_max_level;
666 maxlevel = max_level;
668 reinit(minlevel, maxlevel);
673 template <
typename VectorType>
679 const unsigned int min_level,
680 const unsigned int max_level,
683 , matrix(&matrix,
typeid(*this).name())
684 , coarse(&coarse,
typeid(*this).name())
685 ,
transfer(&transfer,
typeid(*this).name())
686 , pre_smooth(&pre_smooth,
typeid(*this).name())
687 , post_smooth(&post_smooth,
typeid(*this).name())
688 , edge_out(
nullptr,
typeid(*this).name())
689 , edge_in(
nullptr,
typeid(*this).name())
690 , edge_down(
nullptr,
typeid(*this).name())
691 , edge_up(
nullptr,
typeid(*this).name())
697 maxlevel = max_level;
698 reinit(min_level, maxlevel);
703 template <
typename VectorType>
712 template <
typename VectorType>
725 namespace PreconditionMGImplementation
730 typename OtherVectorType>
731 typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
735 const TRANSFER & transfer,
736 OtherVectorType & dst,
737 const OtherVectorType & src,
739 const typename ::mg::Signals &
signals,
742 signals.transfer_to_mg(
true);
743 if (uses_dof_handler_vector)
747 signals.transfer_to_mg(
false);
751 signals.transfer_to_global(
true);
752 if (uses_dof_handler_vector)
756 signals.transfer_to_global(
false);
762 typename OtherVectorType>
767 const TRANSFER & transfer,
768 OtherVectorType & dst,
769 const OtherVectorType & src,
770 const bool uses_dof_handler_vector,
771 const typename ::mg::Signals &signals,
774 (void)uses_dof_handler_vector;
777 signals.transfer_to_mg(
true);
779 signals.transfer_to_mg(
false);
783 signals.transfer_to_global(
true);
785 signals.transfer_to_global(
false);
791 typename OtherVectorType>
792 typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
796 const TRANSFER & transfer,
797 OtherVectorType & dst,
798 const OtherVectorType & src,
799 const bool uses_dof_handler_vector,
800 const typename ::mg::Signals &signals,
803 signals.transfer_to_mg(
true);
804 if (uses_dof_handler_vector)
808 signals.transfer_to_mg(
false);
812 signals.transfer_to_global(
true);
813 if (uses_dof_handler_vector)
819 signals.transfer_to_global(
false);
825 typename OtherVectorType>
830 const TRANSFER & transfer,
831 OtherVectorType & dst,
832 const OtherVectorType & src,
833 const bool uses_dof_handler_vector,
834 const typename ::mg::Signals &signals,
837 (void)uses_dof_handler_vector;
840 signals.transfer_to_mg(
true);
842 signals.transfer_to_mg(
false);
846 signals.transfer_to_global(
true);
850 signals.transfer_to_global(
false);
855 template <
int dim,
typename VectorType,
class TRANSFER>
859 const TRANSFER & transfer)
867 template <
int dim,
typename VectorType,
class TRANSFER>
871 const TRANSFER & transfer)
878 for (
unsigned int i = 0; i < dof_handler.size(); ++i)
885 template <
int dim,
typename VectorType,
class TRANSFER>
892 template <
int dim,
typename VectorType,
class TRANSFER>
893 template <
class OtherVectorType>
896 OtherVectorType & dst,
897 const OtherVectorType &src)
const 910 template <
int dim,
typename VectorType,
class TRANSFER>
913 const unsigned int block)
const 920 template <
int dim,
typename VectorType,
class TRANSFER>
923 const unsigned int block)
const 931 template <
int dim,
typename VectorType,
class TRANSFER>
946 template <
int dim,
typename VectorType,
class TRANSFER>
947 boost::signals2::connection
949 const std::function<
void(
bool)> &slot)
956 template <
int dim,
typename VectorType,
class TRANSFER>
957 boost::signals2::connection
959 const std::function<
void(
bool)> &slot)
966 template <
int dim,
typename VectorType,
class TRANSFER>
967 template <
class OtherVectorType>
970 OtherVectorType & dst,
971 const OtherVectorType &src)
const 984 template <
int dim,
typename VectorType,
class TRANSFER>
985 template <
class OtherVectorType>
988 const OtherVectorType &)
const 994 template <
int dim,
typename VectorType,
class TRANSFER>
995 template <
class OtherVectorType>
999 const OtherVectorType &)
const 1006 DEAL_II_NAMESPACE_CLOSE
boost::signals2::signal< void(const bool before, const unsigned int level)> coarse_solve
boost::signals2::connection connect_transfer_to_global(const std::function< void(bool)> &slot)
void Tvmult(OtherVectorType &dst, const OtherVectorType &src) const
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > matrix
static const unsigned int invalid_unsigned_int
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_down
SmartPointer< const MGMatrixBase< VectorType > > edge_in
unsigned int get_minlevel() const
void vmult_add(OtherVectorType &dst, const OtherVectorType &src) const
unsigned int get_maxlevel() const
virtual unsigned int n_global_levels() const
boost::signals2::signal< void(const bool before, const unsigned int level)> prolongation
#define AssertIndexRange(index, range)
SmartPointer< const MGCoarseGridBase< VectorType >, Multigrid< VectorType > > coarse
boost::signals2::signal< void(const bool before, const unsigned int level)> post_smoother_step
boost::signals2::signal< void(const bool before)> transfer_to_global
boost::signals2::signal< void(const bool before, const unsigned int level)> restriction
MGLevelObject< VectorType > t
const Triangulation< dim, spacedim > & get_triangulation() const
const bool uses_dof_handler_vector
virtual MPI_Comm get_communicator() const
#define Assert(cond, exc)
void vmult(OtherVectorType &dst, const OtherVectorType &src) const
SmartPointer< Multigrid< VectorType >, PreconditionMG< dim, VectorType, TRANSFER > > multigrid
boost::signals2::signal< void(const bool before)> transfer_to_mg
SmartPointer< const MGMatrixBase< VectorType > > edge_out
boost::signals2::signal< void(const bool before, const unsigned int level)> pre_smoother_step
void Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const
boost::signals2::connection connect_transfer_to_mg(const std::function< void(bool)> &slot)
SmartPointer< const TRANSFER, PreconditionMG< dim, VectorType, TRANSFER > > transfer
MGLevelObject< VectorType > solution
std::vector< SmartPointer< const DoFHandler< dim >, PreconditionMG< dim, VectorType, TRANSFER > > > dof_handler_vector
virtual unsigned int get_maxlevel() const =0
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > post_smooth
MPI_Comm get_mpi_communicator() const
std::vector< const DoFHandler< dim > * > dof_handler_vector_raw
static::ExceptionBase & ExcNotImplemented()
Multigrid(const DoFHandler< dim > &mg_dof_handler, 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)
IndexSet locally_owned_range_indices(const unsigned int block=0) const
MGLevelObject< VectorType > defect2
MGLevelObject< VectorType > defect
PreconditionMG(const DoFHandler< dim > &dof_handler, Multigrid< VectorType > &mg, const TRANSFER &transfer)
SmartPointer< const MGTransferBase< VectorType >, Multigrid< VectorType > > transfer
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_up
IndexSet locally_owned_domain_indices(const unsigned int block=0) const
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > pre_smooth
static::ExceptionBase & ExcInternalError()