00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #ifndef __deal2__multigrid_h
00013 #define __deal2__multigrid_h
00014
00015
00016 #include <deal.II/base/config.h>
00017 #include <deal.II/base/subscriptor.h>
00018 #include <deal.II/base/smartpointer.h>
00019 #include <deal.II/numerics/data_out.h>
00020 #include <deal.II/lac/sparse_matrix.h>
00021 #include <deal.II/lac/vector.h>
00022 #include <deal.II/multigrid/mg_base.h>
00023 #include <deal.II/base/mg_level_object.h>
00024 #include <deal.II/multigrid/mg_dof_handler.h>
00025
00026 #include <vector>
00027
00028 DEAL_II_NAMESPACE_OPEN
00029
00032
00062 template <class VECTOR>
00063 class Multigrid : public Subscriptor
00064 {
00065 public:
00069 enum Cycle
00070 {
00072 v_cycle,
00074 w_cycle,
00076 f_cycle
00077 };
00078
00079 typedef VECTOR vector_type;
00080 typedef const VECTOR const_vector_type;
00081
00098 template <int dim>
00099 Multigrid(const MGDoFHandler<dim>& mg_dof_handler,
00100 const MGMatrixBase<VECTOR>& matrix,
00101 const MGCoarseGridBase<VECTOR>& coarse,
00102 const MGTransferBase<VECTOR>& transfer,
00103 const MGSmootherBase<VECTOR>& pre_smooth,
00104 const MGSmootherBase<VECTOR>& post_smooth,
00105 Cycle cycle = v_cycle);
00106
00114 Multigrid(const unsigned int minlevel,
00115 const unsigned int maxlevel,
00116 const MGMatrixBase<VECTOR>& matrix,
00117 const MGCoarseGridBase<VECTOR>& coarse,
00118 const MGTransferBase<VECTOR>& transfer,
00119 const MGSmootherBase<VECTOR>& pre_smooth,
00120 const MGSmootherBase<VECTOR>& post_smooth,
00121 Cycle cycle = v_cycle);
00122
00127 void reinit (const unsigned int minlevel,
00128 const unsigned int maxlevel);
00129
00137 void cycle ();
00138
00160 void vcycle ();
00161
00178 void vmult(VECTOR& dst, const VECTOR& src) const;
00179
00196 void vmult_add(VECTOR& dst, const VECTOR& src) const;
00197
00206 void Tvmult(VECTOR& dst, const VECTOR& src) const;
00207
00216 void Tvmult_add(VECTOR& dst, const VECTOR& src) const;
00217
00246 void set_edge_matrices (const MGMatrixBase<VECTOR>& edge_out,
00247 const MGMatrixBase<VECTOR>& edge_in);
00248
00276 void set_edge_flux_matrices (const MGMatrixBase<VECTOR>& edge_down,
00277 const MGMatrixBase<VECTOR>& edge_up);
00278
00283 unsigned int get_maxlevel() const;
00284
00289 unsigned int get_minlevel() const;
00290
00303 void set_maxlevel (const unsigned int);
00304
00330 void set_minlevel (const unsigned int level,
00331 bool relative = false);
00332
00336 void set_cycle(Cycle);
00337
00344 void set_debug (const unsigned int);
00345
00346 private:
00347
00360 void level_v_step (const unsigned int level);
00361
00376 void level_step (const unsigned int level, Cycle cycle);
00377
00381 Cycle cycle_type;
00382
00386 unsigned int minlevel;
00387
00391 unsigned int maxlevel;
00392
00393 public:
00400 MGLevelObject<VECTOR> defect;
00401
00406 MGLevelObject<VECTOR> solution;
00407
00408 private:
00412 MGLevelObject<VECTOR> t;
00413
00419 MGLevelObject<VECTOR> defect2;
00420
00421
00425 SmartPointer<const MGMatrixBase<VECTOR>,Multigrid<VECTOR> > matrix;
00426
00430 SmartPointer<const MGCoarseGridBase<VECTOR>,Multigrid<VECTOR> > coarse;
00431
00435 SmartPointer<const MGTransferBase<VECTOR>,Multigrid<VECTOR> > transfer;
00436
00440 SmartPointer<const MGSmootherBase<VECTOR>,Multigrid<VECTOR> > pre_smooth;
00441
00445 SmartPointer<const MGSmootherBase<VECTOR>,Multigrid<VECTOR> > post_smooth;
00446
00455 SmartPointer<const MGMatrixBase<VECTOR> > edge_out;
00456
00465 SmartPointer<const MGMatrixBase<VECTOR> > edge_in;
00466
00473 SmartPointer<const MGMatrixBase<VECTOR>,Multigrid<VECTOR> > edge_down;
00474
00481 SmartPointer<const MGMatrixBase<VECTOR>,Multigrid<VECTOR> > edge_up;
00482
00488 unsigned int debug;
00489
00490 template<int dim, class VECTOR2, class TRANSFER> friend class PreconditionMG;
00491 };
00492
00493
00505 template<int dim, class VECTOR, class TRANSFER>
00506 class PreconditionMG : public Subscriptor
00507 {
00508 public:
00515 PreconditionMG(const MGDoFHandler<dim>& mg_dof,
00516 Multigrid<VECTOR>& mg,
00517 const TRANSFER& transfer);
00518
00522 bool empty () const;
00523
00533 template<class VECTOR2>
00534 void vmult (VECTOR2 &dst,
00535 const VECTOR2 &src) const;
00536
00543 template<class VECTOR2>
00544 void vmult_add (VECTOR2 &dst,
00545 const VECTOR2 &src) const;
00546
00553 template<class VECTOR2>
00554 void Tvmult (VECTOR2 &dst,
00555 const VECTOR2 &src) const;
00556
00563 template<class VECTOR2>
00564 void Tvmult_add (VECTOR2 &dst,
00565 const VECTOR2 &src) const;
00566
00567 private:
00571 SmartPointer<const MGDoFHandler<dim>,PreconditionMG<dim,VECTOR,TRANSFER> > mg_dof_handler;
00572
00576 SmartPointer<Multigrid<VECTOR>,PreconditionMG<dim,VECTOR,TRANSFER> > multigrid;
00577
00581 SmartPointer<const TRANSFER,PreconditionMG<dim,VECTOR,TRANSFER> > transfer;
00582 };
00583
00586 #ifndef DOXYGEN
00587
00588
00589
00590 template <class VECTOR>
00591 template <int dim>
00592 Multigrid<VECTOR>::Multigrid (const MGDoFHandler<dim>& mg_dof_handler,
00593 const MGMatrixBase<VECTOR>& matrix,
00594 const MGCoarseGridBase<VECTOR>& coarse,
00595 const MGTransferBase<VECTOR>& transfer,
00596 const MGSmootherBase<VECTOR>& pre_smooth,
00597 const MGSmootherBase<VECTOR>& post_smooth,
00598 Cycle cycle)
00599 :
00600 cycle_type(cycle),
00601 minlevel(0),
00602 maxlevel(mg_dof_handler.get_tria().n_levels()-1),
00603 defect(minlevel,maxlevel),
00604 solution(minlevel,maxlevel),
00605 t(minlevel,maxlevel),
00606 defect2(minlevel,maxlevel),
00607 matrix(&matrix, typeid(*this).name()),
00608 coarse(&coarse, typeid(*this).name()),
00609 transfer(&transfer, typeid(*this).name()),
00610 pre_smooth(&pre_smooth, typeid(*this).name()),
00611 post_smooth(&post_smooth, typeid(*this).name()),
00612 edge_down(0, typeid(*this).name()),
00613 edge_up(0, typeid(*this).name()),
00614 debug(0)
00615 {}
00616
00617
00618
00619 template <class VECTOR>
00620 inline
00621 unsigned int
00622 Multigrid<VECTOR>::get_maxlevel () const
00623 {
00624 return maxlevel;
00625 }
00626
00627
00628
00629 template <class VECTOR>
00630 inline
00631 unsigned int
00632 Multigrid<VECTOR>::get_minlevel () const
00633 {
00634 return minlevel;
00635 }
00636
00637
00638
00639
00640
00641 template<int dim, class VECTOR, class TRANSFER>
00642 PreconditionMG<dim, VECTOR, TRANSFER>
00643 ::PreconditionMG(const MGDoFHandler<dim>& mg_dof_handler,
00644 Multigrid<VECTOR>& mg,
00645 const TRANSFER& transfer)
00646 :
00647 mg_dof_handler(&mg_dof_handler),
00648 multigrid(&mg),
00649 transfer(&transfer)
00650 {}
00651
00652 template<int dim, class VECTOR, class TRANSFER>
00653 inline bool
00654 PreconditionMG<dim, VECTOR, TRANSFER>::empty () const
00655 {
00656 return false;
00657 }
00658
00659 template<int dim, class VECTOR, class TRANSFER>
00660 template<class VECTOR2>
00661 void
00662 PreconditionMG<dim, VECTOR, TRANSFER>::vmult (
00663 VECTOR2& dst,
00664 const VECTOR2& src) const
00665 {
00666 transfer->copy_to_mg(*mg_dof_handler,
00667 multigrid->defect,
00668 src);
00669 multigrid->cycle();
00670
00671 transfer->copy_from_mg(*mg_dof_handler,
00672 dst,
00673 multigrid->solution);
00674 }
00675
00676
00677 template<int dim, class VECTOR, class TRANSFER>
00678 template<class VECTOR2>
00679 void
00680 PreconditionMG<dim, VECTOR, TRANSFER>::vmult_add (
00681 VECTOR2& dst,
00682 const VECTOR2& src) const
00683 {
00684 transfer->copy_to_mg(*mg_dof_handler,
00685 multigrid->defect,
00686 src);
00687 multigrid->cycle();
00688 transfer->copy_from_mg_add(*mg_dof_handler,
00689 dst,
00690 multigrid->solution);
00691 }
00692
00693
00694 template<int dim, class VECTOR, class TRANSFER>
00695 template<class VECTOR2>
00696 void
00697 PreconditionMG<dim, VECTOR, TRANSFER>::Tvmult (
00698 VECTOR2&,
00699 const VECTOR2&) const
00700 {
00701 Assert(false, ExcNotImplemented());
00702 }
00703
00704
00705 template<int dim, class VECTOR, class TRANSFER>
00706 template<class VECTOR2>
00707 void
00708 PreconditionMG<dim, VECTOR, TRANSFER>::Tvmult_add (
00709 VECTOR2&,
00710 const VECTOR2&) const
00711 {
00712 Assert(false, ExcNotImplemented());
00713 }
00714
00715 #endif // DOXYGEN
00716
00717 DEAL_II_NAMESPACE_CLOSE
00718
00719 #endif
00720