Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
multigrid.h
1 // ---------------------------------------------------------------------
2 // @f$Id: multigrid.h 30253 2013-08-07 21:06:29Z bangerth @f$
3 //
4 // Copyright (C) 1999 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__multigrid_h
18 #define __deal2__multigrid_h
19 
20 
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>
30 
31 #include <vector>
32 
34 
37 
67 template <class VECTOR>
68 class Multigrid : public Subscriptor
69 {
70 public:
74  enum Cycle
75  {
82  };
83 
84  typedef VECTOR vector_type;
85  typedef const VECTOR const_vector_type;
86 
103  template <int dim>
104  Multigrid(const DoFHandler<dim> &mg_dof_handler,
110  Cycle cycle = v_cycle);
111 
119  Multigrid(const unsigned int minlevel,
120  const unsigned int maxlevel,
126  Cycle cycle = v_cycle);
127 
132  void reinit (const unsigned int minlevel,
133  const unsigned int maxlevel);
134 
142  void cycle ();
143 
165  void vcycle ();
166 
183  void vmult(VECTOR &dst, const VECTOR &src) const DEAL_II_DEPRECATED;
184 
201  void vmult_add(VECTOR &dst, const VECTOR &src) const DEAL_II_DEPRECATED;
202 
211  void Tvmult(VECTOR &dst, const VECTOR &src) const DEAL_II_DEPRECATED;
212 
221  void Tvmult_add(VECTOR &dst, const VECTOR &src) const DEAL_II_DEPRECATED;
222 
253 
283 
288  unsigned int get_maxlevel() const;
289 
294  unsigned int get_minlevel() const;
295 
308  void set_maxlevel (const unsigned int);
309 
335  void set_minlevel (const unsigned int level,
336  bool relative = false);
337 
341  void set_cycle(Cycle);
342 
349  void set_debug (const unsigned int);
350 
351 private:
352 
365  void level_v_step (const unsigned int level);
366 
381  void level_step (const unsigned int level, Cycle cycle);
382 
387 
391  unsigned int minlevel;
392 
396  unsigned int maxlevel;
397 
398 public:
406 
412 
413 private:
418 
425 
426 
431 
436 
441 
446 
451 
461 
471 
479 
487 
493  unsigned int debug;
494 
495  template<int dim, class VECTOR2, class TRANSFER> friend class PreconditionMG;
496 };
497 
498 
510 template<int dim, class VECTOR, class TRANSFER>
512 {
513 public:
522  const TRANSFER &transfer);
523 
527  bool empty () const;
528 
538  template<class VECTOR2>
539  void vmult (VECTOR2 &dst,
540  const VECTOR2 &src) const;
541 
548  template<class VECTOR2>
549  void vmult_add (VECTOR2 &dst,
550  const VECTOR2 &src) const;
551 
558  template<class VECTOR2>
559  void Tvmult (VECTOR2 &dst,
560  const VECTOR2 &src) const;
561 
568  template<class VECTOR2>
569  void Tvmult_add (VECTOR2 &dst,
570  const VECTOR2 &src) const;
571 
572 private:
577 
582 
587 };
588 
591 #ifndef DOXYGEN
592 /* --------------------------- inline functions --------------------- */
593 
594 
595 template <class VECTOR>
596 template <int dim>
597 Multigrid<VECTOR>::Multigrid (const DoFHandler<dim> &mg_dof_handler,
598  const MGMatrixBase<VECTOR> &matrix,
599  const MGCoarseGridBase<VECTOR> &coarse,
600  const MGTransferBase<VECTOR> &transfer,
601  const MGSmootherBase<VECTOR> &pre_smooth,
602  const MGSmootherBase<VECTOR> &post_smooth,
603  Cycle cycle)
604  :
605  cycle_type(cycle),
606  minlevel(0),
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()),
619  debug(0)
620 {}
621 
622 
623 
624 template <class VECTOR>
625 inline
626 unsigned int
628 {
629  return maxlevel;
630 }
631 
632 
633 
634 template <class VECTOR>
635 inline
636 unsigned int
638 {
639  return minlevel;
640 }
641 
642 
643 /* --------------------------- inline functions --------------------- */
644 
645 
646 template<int dim, class VECTOR, class TRANSFER>
648 ::PreconditionMG(const DoFHandler<dim> &dof_handler,
650  const TRANSFER &transfer)
651  :
652  dof_handler(&dof_handler),
653  multigrid(&mg),
654  transfer(&transfer)
655 {}
656 
657 template<int dim, class VECTOR, class TRANSFER>
658 inline bool
660 {
661  return false;
662 }
663 
664 template<int dim, class VECTOR, class TRANSFER>
665 template<class VECTOR2>
666 void
668  VECTOR2 &dst,
669  const VECTOR2 &src) const
670 {
671  transfer->copy_to_mg(*dof_handler,
672  multigrid->defect,
673  src);
674  multigrid->cycle();
675 
676  transfer->copy_from_mg(*dof_handler,
677  dst,
678  multigrid->solution);
679 }
680 
681 
682 template<int dim, class VECTOR, class TRANSFER>
683 template<class VECTOR2>
684 void
686  VECTOR2 &dst,
687  const VECTOR2 &src) const
688 {
689  transfer->copy_to_mg(*dof_handler,
690  multigrid->defect,
691  src);
692  multigrid->cycle();
693  transfer->copy_from_mg_add(*dof_handler,
694  dst,
695  multigrid->solution);
696 }
697 
698 
699 template<int dim, class VECTOR, class TRANSFER>
700 template<class VECTOR2>
701 void
703  VECTOR2 &,
704  const VECTOR2 &) const
705 {
706  Assert(false, ExcNotImplemented());
707 }
708 
709 
710 template<int dim, class VECTOR, class TRANSFER>
711 template<class VECTOR2>
712 void
714  VECTOR2 &,
715  const VECTOR2 &) const
716 {
717  Assert(false, ExcNotImplemented());
718 }
719 
720 #endif // DOXYGEN
721 
722 DEAL_II_NAMESPACE_CLOSE
723 
724 #endif
SmartPointer< const MGMatrixBase< VECTOR > > edge_out
Definition: multigrid.h:460
SmartPointer< const MGMatrixBase< VECTOR > > edge_in
Definition: multigrid.h:470
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
Definition: multigrid.h:581
void Tvmult_add(VECTOR &dst, const VECTOR &src) const DEAL_II_DEPRECATED
SmartPointer< const MGCoarseGridBase< VECTOR >, Multigrid< VECTOR > > coarse
Definition: multigrid.h:435
void set_edge_matrices(const MGMatrixBase< VECTOR > &edge_out, const MGMatrixBase< VECTOR > &edge_in)
MGLevelObject< VECTOR > solution
Definition: multigrid.h:411
unsigned int get_minlevel() const
void vmult_add(VECTOR2 &dst, const VECTOR2 &src) const
Definition: mg.h:82
unsigned int minlevel
Definition: multigrid.h:391
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
Definition: multigrid.h:440
SmartPointer< const MGMatrixBase< VECTOR >, Multigrid< VECTOR > > edge_down
Definition: multigrid.h:478
SmartPointer< const MGSmootherBase< VECTOR >, Multigrid< VECTOR > > pre_smooth
Definition: multigrid.h:445
MGLevelObject< VECTOR > defect2
Definition: multigrid.h:424
bool empty() const
MGLevelObject< VECTOR > t
Definition: multigrid.h:417
unsigned int get_maxlevel() const
void Tvmult_add(VECTOR2 &dst, const VECTOR2 &src) const
#define Assert(cond, exc)
Definition: exceptions.h:299
SmartPointer< const MGMatrixBase< VECTOR >, Multigrid< VECTOR > > edge_up
Definition: multigrid.h:486
void set_cycle(Cycle)
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
Definition: multigrid.h:450
void set_minlevel(const unsigned int level, bool relative=false)
SmartPointer< const MGMatrixBase< VECTOR >, Multigrid< VECTOR > > matrix
Definition: multigrid.h:430
The W-cycle.
Definition: multigrid.h:79
void reinit(const unsigned int minlevel, const unsigned int maxlevel)
unsigned int maxlevel
Definition: multigrid.h:396
unsigned int debug
Definition: multigrid.h:493
SmartPointer< const TRANSFER, PreconditionMG< dim, VECTOR, TRANSFER > > transfer
Definition: multigrid.h:586
Cycle cycle_type
Definition: multigrid.h:386
void set_edge_flux_matrices(const MGMatrixBase< VECTOR > &edge_down, const MGMatrixBase< VECTOR > &edge_up)
::ExceptionBase & ExcNotImplemented()
The V-cycle.
Definition: multigrid.h:77
The F-cycle.
Definition: multigrid.h:81
SmartPointer< const DoFHandler< dim >, PreconditionMG< dim, VECTOR, TRANSFER > > dof_handler
Definition: multigrid.h:576
void vmult(VECTOR2 &dst, const VECTOR2 &src) const
void set_debug(const unsigned int)
MGLevelObject< VECTOR > defect
Definition: multigrid.h:405