16 #ifndef dealii_multigrid_h 17 #define dealii_multigrid_h 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." 85 boost::signals2::signal<void(const bool before, const unsigned int level)>
96 boost::signals2::signal<void(const bool before, const unsigned int level)>
104 boost::signals2::signal<void(const bool before, const unsigned int level)>
115 boost::signals2::signal<void(const bool before, const unsigned int level)>
124 boost::signals2::signal<void(const bool before, const unsigned int level)>
151 template <
typename VectorType>
185 const unsigned int minlevel = 0,
187 Cycle cycle = v_cycle);
193 reinit(
const unsigned int minlevel,
const unsigned int maxlevel);
251 get_maxlevel()
const;
257 get_minlevel()
const;
265 set_maxlevel(
const unsigned int);
283 set_minlevel(
const unsigned int level,
bool relative =
false);
288 void set_cycle(
Cycle);
293 boost::signals2::connection
294 connect_coarse_solve(
295 const std::function<
void(
const bool,
const unsigned int)> &slot);
300 boost::signals2::connection
302 const std::function<
void(
const bool,
const unsigned int)> &slot);
307 boost::signals2::connection
308 connect_prolongation(
309 const std::function<
void(
const bool,
const unsigned int)> &slot);
314 boost::signals2::connection
315 connect_pre_smoother_step(
316 const std::function<
void(
const bool,
const unsigned int)> &slot);
321 boost::signals2::connection
322 connect_post_smoother_step(
323 const std::function<
void(
const bool,
const unsigned int)> &slot);
338 level_v_step(
const unsigned int level);
348 level_step(
const unsigned int level,
Cycle cycle);
447 template <
int dim,
class OtherVectorType,
class TRANSFER>
466 template <
int dim,
typename VectorType,
class TRANSFER>
476 const TRANSFER & transfer);
484 const TRANSFER & transfer);
498 template <
class OtherVectorType>
500 vmult(OtherVectorType &dst,
const OtherVectorType &src)
const;
506 template <
class OtherVectorType>
508 vmult_add(OtherVectorType &dst,
const OtherVectorType &src)
const;
515 template <
class OtherVectorType>
517 Tvmult(OtherVectorType &dst,
const OtherVectorType &src)
const;
524 template <
class OtherVectorType>
526 Tvmult_add(OtherVectorType &dst,
const OtherVectorType &src)
const;
535 locally_owned_range_indices(
const unsigned int block = 0)
const;
544 locally_owned_domain_indices(
const unsigned int block = 0)
const;
550 get_mpi_communicator()
const;
555 boost::signals2::connection
556 connect_transfer_to_mg(
const std::function<
void(
bool)> &slot);
561 boost::signals2::connection
562 connect_transfer_to_global(
const std::function<
void(
bool)> &slot);
568 std::vector<SmartPointer<const DoFHandler<dim>,
608 template <
typename VectorType>
614 const unsigned int min_level,
615 const unsigned int max_level,
618 , matrix(&matrix, typeid(*this).name())
619 , coarse(&coarse, typeid(*this).name())
620 , transfer(&transfer, typeid(*this).name())
621 , pre_smooth(&pre_smooth, typeid(*this).name())
622 , post_smooth(&post_smooth, typeid(*this).name())
623 , edge_out(nullptr, typeid(*this).name())
624 , edge_in(nullptr, typeid(*this).name())
625 , edge_down(nullptr, typeid(*this).name())
626 , edge_up(nullptr, typeid(*this).name())
631 maxlevel = max_level;
632 reinit(min_level, maxlevel);
637 template <
typename VectorType>
646 template <
typename VectorType>
659 namespace PreconditionMGImplementation
664 typename OtherVectorType>
665 typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
670 OtherVectorType & dst,
671 const OtherVectorType & src,
673 const typename ::mg::Signals &
signals,
676 signals.transfer_to_mg(
true);
677 if (uses_dof_handler_vector)
681 signals.transfer_to_mg(
false);
685 signals.transfer_to_global(
true);
686 if (uses_dof_handler_vector)
690 signals.transfer_to_global(
false);
696 typename OtherVectorType>
701 const TRANSFER & transfer,
702 OtherVectorType & dst,
703 const OtherVectorType & src,
704 const bool uses_dof_handler_vector,
705 const typename ::mg::Signals &signals,
708 (void)uses_dof_handler_vector;
711 signals.transfer_to_mg(
true);
713 signals.transfer_to_mg(
false);
717 signals.transfer_to_global(
true);
719 signals.transfer_to_global(
false);
725 typename OtherVectorType>
726 typename std::enable_if<TRANSFER::supports_dof_handler_vector>::type
730 const TRANSFER & transfer,
731 OtherVectorType & dst,
732 const OtherVectorType & src,
733 const bool uses_dof_handler_vector,
734 const typename ::mg::Signals &signals,
737 signals.transfer_to_mg(
true);
738 if (uses_dof_handler_vector)
742 signals.transfer_to_mg(
false);
746 signals.transfer_to_global(
true);
747 if (uses_dof_handler_vector)
753 signals.transfer_to_global(
false);
759 typename OtherVectorType>
764 const TRANSFER & transfer,
765 OtherVectorType & dst,
766 const OtherVectorType & src,
767 const bool uses_dof_handler_vector,
768 const typename ::mg::Signals &signals,
771 (void)uses_dof_handler_vector;
774 signals.transfer_to_mg(
true);
776 signals.transfer_to_mg(
false);
780 signals.transfer_to_global(
true);
784 signals.transfer_to_global(
false);
789 template <
int dim,
typename VectorType,
class TRANSFER>
801 template <
int dim,
typename VectorType,
class TRANSFER>
805 const TRANSFER & transfer)
812 for (
unsigned int i = 0; i < dof_handler.size(); ++i)
819 template <
int dim,
typename VectorType,
class TRANSFER>
826 template <
int dim,
typename VectorType,
class TRANSFER>
827 template <
class OtherVectorType>
830 OtherVectorType & dst,
831 const OtherVectorType &src)
const 844 template <
int dim,
typename VectorType,
class TRANSFER>
847 const unsigned int block)
const 854 template <
int dim,
typename VectorType,
class TRANSFER>
857 const unsigned int block)
const 865 template <
int dim,
typename VectorType,
class TRANSFER>
880 template <
int dim,
typename VectorType,
class TRANSFER>
881 boost::signals2::connection
883 const std::function<
void(
bool)> &slot)
890 template <
int dim,
typename VectorType,
class TRANSFER>
891 boost::signals2::connection
893 const std::function<
void(
bool)> &slot)
900 template <
int dim,
typename VectorType,
class TRANSFER>
901 template <
class OtherVectorType>
904 OtherVectorType & dst,
905 const OtherVectorType &src)
const 918 template <
int dim,
typename VectorType,
class TRANSFER>
919 template <
class OtherVectorType>
922 const OtherVectorType &)
const 928 template <
int dim,
typename VectorType,
class TRANSFER>
929 template <
class OtherVectorType>
933 const OtherVectorType &)
const 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)
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > matrix
static const unsigned int invalid_unsigned_int
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
SmartPointer< const MGMatrixBase< VectorType >, Multigrid< VectorType > > edge_down
SmartPointer< const MGMatrixBase< VectorType > > edge_in
void vmult(OtherVectorType &dst, const OtherVectorType &src) const
Contents is actually a matrix.
unsigned int get_maxlevel() 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 bool uses_dof_handler_vector
IndexSet locally_owned_range_indices(const unsigned int block=0) const
#define Assert(cond, exc)
void Tvmult(OtherVectorType &dst, const OtherVectorType &src) const
unsigned int get_minlevel() const
#define DEAL_II_NAMESPACE_CLOSE
SmartPointer< Multigrid< VectorType >, PreconditionMG< dim, VectorType, TRANSFER > > multigrid
boost::signals2::signal< void(const bool before)> transfer_to_mg
const MPI_Comm & get_mpi_communicator() const
SmartPointer< const MGMatrixBase< VectorType > > edge_out
void vmult_add(OtherVectorType &dst, const OtherVectorType &src) const
boost::signals2::signal< void(const bool before, const unsigned int level)> pre_smoother_step
virtual const MPI_Comm & get_communicator() const
boost::signals2::connection connect_transfer_to_mg(const std::function< void(bool)> &slot)
void Tvmult_add(OtherVectorType &dst, const OtherVectorType &src) const
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
IndexSet locally_owned_domain_indices(const unsigned int block=0) const
#define DEAL_II_NAMESPACE_OPEN
std::vector< const DoFHandler< dim > * > dof_handler_vector_raw
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
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
SmartPointer< const MGSmootherBase< VectorType >, Multigrid< VectorType > > pre_smooth
static ::ExceptionBase & ExcInternalError()