MFEM  v4.6.0
Finite element discretization library
multigrid.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, 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 /// Abstract base class for Multigrid solvers
25 class MultigridBase : public Solver
26 {
27 public:
28  enum class CycleType
29  {
30  VCYCLE,
31  WCYCLE
32  };
33 
34 protected:
39 
43 
44  mutable Array2D<Vector*> X, Y, R, Z;
45  mutable int nrhs;
46 
47 public:
48  /// Constructs an empty multigrid hierarchy
49  MultigridBase();
50 
51  /// Constructs a multigrid hierarchy from the given inputs
52  /** Inputs include operators and smoothers on all levels, and ownership of
53  the given operators and smoothers */
54  MultigridBase(const Array<Operator*>& operators_,
55  const Array<Solver*>& smoothers_,
56  const Array<bool>& ownedOperators_,
57  const Array<bool>& ownedSmoothers_);
58 
59  /// Destructor
60  virtual ~MultigridBase();
61 
62  /// Adds a level to the multigrid operator hierarchy
63  /** The ownership of the operators and solvers/smoothers may be transferred
64  to the Multigrid by setting the according boolean variables */
65  void AddLevel(Operator* op, Solver* smoother, bool ownOperator,
66  bool ownSmoother);
67 
68  /// Returns the number of levels
69  int NumLevels() const { return operators.Size(); }
70 
71  /// Returns the index of the finest level
72  int GetFinestLevelIndex() const { return NumLevels() - 1; }
73 
74  /// Returns operator at given level
75  const Operator* GetOperatorAtLevel(int level) const
76  {
77  return operators[level];
78  }
80  {
81  return operators[level];
82  }
83 
84  /// Returns operator at finest level
86  {
88  }
90  {
92  }
93 
94  /// Returns smoother at given level
95  const Solver* GetSmootherAtLevel(int level) const
96  {
97  return smoothers[level];
98  }
100  {
101  return smoothers[level];
102  }
103 
104  /// Set cycle type and number of pre- and post-smoothing steps used by Mult
105  void SetCycleType(CycleType cycleType_, int preSmoothingSteps_,
106  int postSmoothingSteps_);
107 
108  /// Application of the multigrid as a preconditioner
109  virtual void Mult(const Vector& x, Vector& y) const override;
110  virtual void ArrayMult(const Array<const Vector*>& X_,
111  Array<Vector*>& Y_) const override;
112 
113  /// Not supported for multigrid
114  virtual void SetOperator(const Operator& op) override
115  {
116  MFEM_ABORT("SetOperator is not supported in Multigrid!");
117  }
118 
119 private:
120  /// Application of a multigrid cycle at particular level
121  void Cycle(int level) const;
122 
123  /// Application of a pre-/post-smoothing step at particular level
124  void SmoothingStep(int level, bool zero, bool transpose) const;
125 
126  /// Allocate or destroy temporary storage
127  void InitVectors() const;
128  void EraseVectors() const;
129 
130  /// Returns prolongation operator at given level
131  virtual const Operator* GetProlongationAtLevel(int level) const = 0;
132 };
133 
134 /// Multigrid solver class
135 class Multigrid : public MultigridBase
136 {
137 protected:
140 
141 public:
142  /// Constructs an empty multigrid hierarchy
143  Multigrid();
144 
145  /// Constructs a multigrid hierarchy from the given inputs
146  /** Inputs include operators and smoothers on all levels, prolongation
147  operators that go from coarser to finer levels, and ownership of the
148  given operators, smoothers, and prolongations */
149  Multigrid(const Array<Operator*>& operators_, const Array<Solver*>& smoothers_,
150  const Array<Operator*>& prolongations_, const Array<bool>& ownedOperators_,
151  const Array<bool>& ownedSmoothers_, const Array<bool>& ownedProlongations_);
152 
153  /// Destructor
154  virtual ~Multigrid();
155 
156 private:
157  /// Returns prolongation operator at given level
158  virtual const Operator* GetProlongationAtLevel(int level) const override
159  {
160  return prolongations[level];
161  }
162 };
163 
164 /// Geometric multigrid associated with a hierarchy of finite element spaces
166 {
167 protected:
171 
172 public:
173  /** Construct an empty multigrid object for the given finite element space
174  hierarchy @a fespaces_ */
176 
177  /// Destructor
178  virtual ~GeometricMultigrid();
179 
180  /** Form the linear system A X = B, corresponding to the operator on the
181  finest level of the geometric multigrid hierarchy */
183  Vector& B);
184 
185  /// Recover the solution of a linear system formed with FormFineLinearSystem()
186  void RecoverFineFEMSolution(const Vector& X, const Vector& b, Vector& x);
187 
188 private:
189  /// Returns prolongation operator at given level
190  virtual const Operator* GetProlongationAtLevel(int level) const override
191  {
192  return fespaces.GetProlongationAtLevel(level);
193  }
194 };
195 
196 } // namespace mfem
197 
198 #endif
Multigrid solver class.
Definition: multigrid.hpp:135
const Operator * GetOperatorAtFinestLevel() const
Returns operator at finest level.
Definition: multigrid.hpp:85
Operator * GetProlongationAtLevel(int level) const
Returns the prolongation operator from the finite element space at level to the finite element space ...
Array< bool > ownedOperators
Definition: multigrid.hpp:37
Geometric multigrid associated with a hierarchy of finite element spaces.
Definition: multigrid.hpp:165
Array2D< Vector * > Z
Definition: multigrid.hpp:44
virtual ~GeometricMultigrid()
Destructor.
Definition: multigrid.cpp:266
Array2D< Vector * > X
Definition: multigrid.hpp:44
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
Solver * GetSmootherAtLevel(int level)
Definition: multigrid.hpp:99
Array< BilinearForm * > bfs
Definition: multigrid.hpp:170
Array2D< Vector * > Y
Definition: multigrid.hpp:44
Array< Array< int > * > essentialTrueDofs
Definition: multigrid.hpp:169
Array< Operator * > prolongations
Definition: multigrid.hpp:138
Array< Operator * > operators
Definition: multigrid.hpp:35
Operator * GetOperatorAtLevel(int level)
Definition: multigrid.hpp:79
Array< Solver * > smoothers
Definition: multigrid.hpp:36
Array< bool > ownedSmoothers
Definition: multigrid.hpp:38
Array2D< Vector * > R
Definition: multigrid.hpp:44
double b
Definition: lissajous.cpp:42
const Solver * GetSmootherAtLevel(int level) const
Returns smoother at given level.
Definition: multigrid.hpp:95
int GetFinestLevelIndex() const
Returns the index of the finest level.
Definition: multigrid.hpp:72
const Operator * GetOperatorAtLevel(int level) const
Returns operator at given level.
Definition: multigrid.hpp:75
virtual void Mult(const Vector &x, Vector &y) const override
Application of the multigrid as a preconditioner.
Definition: multigrid.cpp:106
CycleType cycleType
Definition: multigrid.hpp:40
virtual void SetOperator(const Operator &op) override
Not supported for multigrid.
Definition: multigrid.hpp:114
virtual void ArrayMult(const Array< const Vector *> &X_, Array< Vector *> &Y_) const override
Operator application on a matrix: Y=A(X).
Definition: multigrid.cpp:115
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:98
Dynamic 2D array using row-major layout.
Definition: array.hpp:354
const FiniteElementSpaceHierarchy & fespaces
Definition: multigrid.hpp:168
virtual ~MultigridBase()
Destructor.
Definition: multigrid.cpp:36
Array< bool > ownedProlongations
Definition: multigrid.hpp:139
virtual ~Multigrid()
Destructor.
Definition: multigrid.cpp:250
MultigridBase()
Constructs an empty multigrid hierarchy.
Definition: multigrid.cpp:17
void RecoverFineFEMSolution(const Vector &X, const Vector &b, Vector &x)
Recover the solution of a linear system formed with FormFineLinearSystem()
Definition: multigrid.cpp:285
int NumLevels() const
Returns the number of levels.
Definition: multigrid.hpp:69
void FormFineLinearSystem(Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B)
Definition: multigrid.cpp:278
Vector data type.
Definition: vector.hpp:58
void AddLevel(Operator *op, Solver *smoother, bool ownOperator, bool ownSmoother)
Adds a level to the multigrid operator hierarchy.
Definition: multigrid.cpp:87
Multigrid()
Constructs an empty multigrid hierarchy.
Definition: multigrid.cpp:234
Base class for solvers.
Definition: operator.hpp:682
GeometricMultigrid(const FiniteElementSpaceHierarchy &fespaces_)
Definition: multigrid.cpp:262
Abstract operator.
Definition: operator.hpp:24
Operator * GetOperatorAtFinestLevel()
Definition: multigrid.hpp:89
Abstract base class for Multigrid solvers.
Definition: multigrid.hpp:25