MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
multigrid.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_MULTIGRID
13 #define MFEM_MULTIGRID
14 
15 #include "fespacehierarchy.hpp"
16 #include "bilinearform.hpp"
17 
18 #include "../linalg/operator.hpp"
19 #include "../linalg/handle.hpp"
20 
21 namespace mfem
22 {
23 
24 /// Multigrid solver class
25 class Multigrid : public Solver
26 {
27 public:
28  enum class CycleType
29  {
30  VCYCLE,
31  WCYCLE
32  };
33 
34 protected:
38 
42 
46 
47  mutable Array<Vector*> X;
48  mutable Array<Vector*> Y;
49  mutable Array<Vector*> R;
50  mutable Array<Vector*> Z;
51 
52 public:
53  /// Constructs an empty multigrid hierarchy.
54  Multigrid();
55 
56  /// Constructs a multigrid hierarchy from the given inputs.
57  /** Inputs include operators and smoothers on all levels, prolongation
58  operators that go from coarser to finer levels, and ownership of the
59  given operators, smoothers, and prolongations. */
60  Multigrid(const Array<Operator*>& operators_, const Array<Solver*>& smoothers_,
61  const Array<Operator*>& prolongations_, const Array<bool>& ownedOperators_,
62  const Array<bool>& ownedSmoothers_, const Array<bool>& ownedProlongations_);
63 
64  /// Destructor
65  virtual ~Multigrid();
66 
67  /// Adds a level to the multigrid operator hierarchy.
68  /** The ownership of the operators and solvers/smoothers may be transferred
69  to the Multigrid by setting the according boolean variables. */
70  void AddLevel(Operator* opr, Solver* smoother, bool ownOperator,
71  bool ownSmoother);
72 
73  /// Returns the number of levels
74  int NumLevels() const;
75 
76  /// Returns the index of the finest level
77  int GetFinestLevelIndex() const;
78 
79  /// Returns operator at given level
80  const Operator* GetOperatorAtLevel(int level) const;
81 
82  /// Returns operator at given level
83  Operator* GetOperatorAtLevel(int level);
84 
85  /// Returns operator at finest level
86  const Operator* GetOperatorAtFinestLevel() const;
87 
88  /// Returns operator at finest level
90 
91  /// Returns smoother at given level
92  Solver* GetSmootherAtLevel(int level) const;
93 
94  /// Returns smoother at given level
95  Solver* GetSmootherAtLevel(int level);
96 
97  /// Set cycle type and number of pre- and post-smoothing steps used by Mult
98  void SetCycleType(CycleType cycleType_, int preSmoothingSteps_,
99  int postSmoothingSteps_);
100 
101  /// Application of the multigrid as a preconditioner
102  virtual void Mult(const Vector& x, Vector& y) const override;
103 
104  /// Not supported for multigrid
105  virtual void SetOperator(const Operator& op) override;
106 
107 private:
108  /// Application of a smoothing step at particular level
109  void SmoothingStep(int level, bool transpose) const;
110 
111  /// Application of a multigrid cycle at particular level
112  void Cycle(int level) const;
113 
114  /// Returns prolongation operator at given level
115  virtual const Operator* GetProlongationAtLevel(int level) const;
116 };
117 
118 /// Geometric multigrid associated with a hierarchy of finite element spaces
120 {
121 protected:
125 
126 public:
127  /** Construct an empty multigrid object for the given finite element space
128  hierarchy @a fespaces_ */
130  : Multigrid(), fespaces(fespaces_) { }
131 
132  /// Destructor
133  virtual ~GeometricMultigrid();
134 
135  /** Form the linear system A X = B, corresponding to the operator on the
136  finest level of the geometric multigrid hierarchy */
138  Vector& B);
139 
140  /// Recover the solution of a linear system formed with FormFineLinearSystem()
141  void RecoverFineFEMSolution(const Vector& X, const Vector& b, Vector& x);
142 
143 private:
144  /// Returns prolongation operator at given level
145  virtual const Operator* GetProlongationAtLevel(int level) const override;
146 };
147 
148 } // namespace mfem
149 
150 #endif
Array< Vector * > Z
Definition: multigrid.hpp:50
Multigrid solver class.
Definition: multigrid.hpp:25
Array< Vector * > X
Definition: multigrid.hpp:47
Geometric multigrid associated with a hierarchy of finite element spaces.
Definition: multigrid.hpp:119
virtual ~GeometricMultigrid()
Destructor.
Definition: multigrid.cpp:229
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
Array< BilinearForm * > bfs
Definition: multigrid.hpp:124
Array< Vector * > Y
Definition: multigrid.hpp:48
Array< Array< int > * > essentialTrueDofs
Definition: multigrid.hpp:123
Array< Operator * > prolongations
Definition: multigrid.hpp:37
CycleType cycleType
Definition: multigrid.hpp:43
double b
Definition: lissajous.cpp:42
Array< Vector * > R
Definition: multigrid.hpp:49
const FiniteElementSpaceHierarchy & fespaces
Definition: multigrid.hpp:122
Array< bool > ownedProlongations
Definition: multigrid.hpp:41
virtual void SetOperator(const Operator &op) override
Not supported for multigrid.
Definition: multigrid.cpp:157
virtual ~Multigrid()
Destructor.
Definition: multigrid.cpp:51
Array< Solver * > smoothers
Definition: multigrid.hpp:36
void AddLevel(Operator *opr, Solver *smoother, bool ownOperator, bool ownSmoother)
Adds a level to the multigrid operator hierarchy.
Definition: multigrid.cpp:86
Solver * GetSmootherAtLevel(int level) const
Returns smoother at given level.
Definition: multigrid.cpp:130
const Operator * GetOperatorAtFinestLevel() const
Returns operator at finest level.
Definition: multigrid.cpp:120
virtual void Mult(const Vector &x, Vector &y) const override
Application of the multigrid as a preconditioner.
Definition: multigrid.cpp:148
void RecoverFineFEMSolution(const Vector &X, const Vector &b, Vector &x)
Recover the solution of a linear system formed with FormFineLinearSystem()
Definition: multigrid.cpp:253
int GetFinestLevelIndex() const
Returns the index of the finest level.
Definition: multigrid.cpp:108
void FormFineLinearSystem(Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B)
Definition: multigrid.cpp:246
int NumLevels() const
Returns the number of levels.
Definition: multigrid.cpp:106
const Operator * GetOperatorAtLevel(int level) const
Returns operator at given level.
Definition: multigrid.cpp:110
Array< Operator * > operators
Definition: multigrid.hpp:35
Vector data type.
Definition: vector.hpp:60
Multigrid()
Constructs an empty multigrid hierarchy.
Definition: multigrid.cpp:17
int postSmoothingSteps
Definition: multigrid.hpp:45
Array< bool > ownedSmoothers
Definition: multigrid.hpp:40
Base class for solvers.
Definition: operator.hpp:655
GeometricMultigrid(const FiniteElementSpaceHierarchy &fespaces_)
Definition: multigrid.hpp:129
Abstract operator.
Definition: operator.hpp:24
Array< bool > ownedOperators
Definition: multigrid.hpp:39
void SetCycleType(CycleType cycleType_, int preSmoothingSteps_, int postSmoothingSteps_)
Set cycle type and number of pre- and post-smoothing steps used by Mult.
Definition: multigrid.cpp:140