17 #ifndef __deal2__multigrid_h
18 #define __deal2__multigrid_h
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/subscriptor.h>
23 #include <deal.II/base/smartpointer.h>
24 #include <deal.II/dofs/dof_handler.h>
25 #include <deal.II/lac/sparse_matrix.h>
26 #include <deal.II/lac/vector.h>
27 #include <deal.II/multigrid/mg_base.h>
28 #include <deal.II/base/mg_level_object.h>
29 #include <deal.II/multigrid/mg_dof_handler.h>
67 template <
class VECTOR>
84 typedef VECTOR vector_type;
85 typedef const VECTOR const_vector_type;
336 bool relative =
false);
495 template<
int dim,
class VECTOR2,
class TRANSFER>
friend class PreconditionMG;
510 template<
int dim,
class VECTOR,
class TRANSFER>
538 template<
class VECTOR2>
539 void vmult (VECTOR2 &dst,
540 const VECTOR2 &src)
const;
548 template<
class VECTOR2>
550 const VECTOR2 &src)
const;
558 template<
class VECTOR2>
559 void Tvmult (VECTOR2 &dst,
560 const VECTOR2 &src)
const;
568 template<
class VECTOR2>
570 const VECTOR2 &src)
const;
595 template <
class VECTOR>
607 maxlevel(mg_dof_handler.get_tria().n_global_levels()-1),
608 defect(minlevel,maxlevel),
609 solution(minlevel,maxlevel),
610 t(minlevel,maxlevel),
611 defect2(minlevel,maxlevel),
612 matrix(&matrix, typeid(*this).name()),
613 coarse(&coarse, typeid(*this).name()),
614 transfer(&transfer, typeid(*this).name()),
615 pre_smooth(&pre_smooth, typeid(*this).name()),
616 post_smooth(&post_smooth, typeid(*this).name()),
617 edge_down(0, typeid(*this).name()),
618 edge_up(0, typeid(*this).name()),
624 template <
class VECTOR>
634 template <
class VECTOR>
646 template<
int dim,
class VECTOR,
class TRANSFER>
650 const TRANSFER &transfer)
652 dof_handler(&dof_handler),
657 template<
int dim,
class VECTOR,
class TRANSFER>
664 template<
int dim,
class VECTOR,
class TRANSFER>
665 template<
class VECTOR2>
669 const VECTOR2 &src)
const
671 transfer->copy_to_mg(*dof_handler,
676 transfer->copy_from_mg(*dof_handler,
678 multigrid->solution);
682 template<
int dim,
class VECTOR,
class TRANSFER>
683 template<
class VECTOR2>
687 const VECTOR2 &src)
const
689 transfer->copy_to_mg(*dof_handler,
693 transfer->copy_from_mg_add(*dof_handler,
695 multigrid->solution);
699 template<
int dim,
class VECTOR,
class TRANSFER>
700 template<
class VECTOR2>
704 const VECTOR2 &)
const
710 template<
int dim,
class VECTOR,
class TRANSFER>
711 template<
class VECTOR2>
715 const VECTOR2 &)
const
722 DEAL_II_NAMESPACE_CLOSE
SmartPointer< const MGMatrixBase< VECTOR > > edge_out
SmartPointer< const MGMatrixBase< VECTOR > > edge_in
Multigrid(const DoFHandler< dim > &mg_dof_handler, const MGMatrixBase< VECTOR > &matrix, const MGCoarseGridBase< VECTOR > &coarse, const MGTransferBase< VECTOR > &transfer, const MGSmootherBase< VECTOR > &pre_smooth, const MGSmootherBase< VECTOR > &post_smooth, Cycle cycle=v_cycle)
SmartPointer< Multigrid< VECTOR >, PreconditionMG< dim, VECTOR, TRANSFER > > multigrid
void Tvmult_add(VECTOR &dst, const VECTOR &src) const DEAL_II_DEPRECATED
SmartPointer< const MGCoarseGridBase< VECTOR >, Multigrid< VECTOR > > coarse
void set_edge_matrices(const MGMatrixBase< VECTOR > &edge_out, const MGMatrixBase< VECTOR > &edge_in)
MGLevelObject< VECTOR > solution
unsigned int get_minlevel() const
void vmult_add(VECTOR2 &dst, const VECTOR2 &src) const
void Tvmult(VECTOR &dst, const VECTOR &src) const DEAL_II_DEPRECATED
void Tvmult(VECTOR2 &dst, const VECTOR2 &src) const
void level_step(const unsigned int level, Cycle cycle)
void level_v_step(const unsigned int level)
void vmult_add(VECTOR &dst, const VECTOR &src) const DEAL_II_DEPRECATED
SmartPointer< const MGTransferBase< VECTOR >, Multigrid< VECTOR > > transfer
SmartPointer< const MGMatrixBase< VECTOR >, Multigrid< VECTOR > > edge_down
SmartPointer< const MGSmootherBase< VECTOR >, Multigrid< VECTOR > > pre_smooth
MGLevelObject< VECTOR > defect2
MGLevelObject< VECTOR > t
unsigned int get_maxlevel() const
void Tvmult_add(VECTOR2 &dst, const VECTOR2 &src) const
#define Assert(cond, exc)
SmartPointer< const MGMatrixBase< VECTOR >, Multigrid< VECTOR > > edge_up
void vmult(VECTOR &dst, const VECTOR &src) const DEAL_II_DEPRECATED
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
PreconditionMG(const DoFHandler< dim > &dof_handler, Multigrid< VECTOR > &mg, const TRANSFER &transfer)
void set_maxlevel(const unsigned int)
SmartPointer< const MGSmootherBase< VECTOR >, Multigrid< VECTOR > > post_smooth
void set_minlevel(const unsigned int level, bool relative=false)
SmartPointer< const MGMatrixBase< VECTOR >, Multigrid< VECTOR > > matrix
void reinit(const unsigned int minlevel, const unsigned int maxlevel)
SmartPointer< const TRANSFER, PreconditionMG< dim, VECTOR, TRANSFER > > transfer
void set_edge_flux_matrices(const MGMatrixBase< VECTOR > &edge_down, const MGMatrixBase< VECTOR > &edge_up)
::ExceptionBase & ExcNotImplemented()
SmartPointer< const DoFHandler< dim >, PreconditionMG< dim, VECTOR, TRANSFER > > dof_handler
void vmult(VECTOR2 &dst, const VECTOR2 &src) const
void set_debug(const unsigned int)
MGLevelObject< VECTOR > defect