ThunderEgg  1.0.0
Cycle.h
Go to the documentation of this file.
1 /***************************************************************************
2  * ThunderEgg, a library for solvers on adaptively refined block-structured
3  * Cartesian grids.
4  *
5  * Copyright (c) 2018-2021 Scott Aiton
6  *
7  * This program is free software: you can redistribute it and/or modify
8  * it under the terms of the GNU General Public License as published by
9  * the Free Software Foundation, either version 3 of the License, or
10  * (at your option) any later version.
11  *
12  * This program is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  * GNU General Public License for more details.
16  *
17  * You should have received a copy of the GNU General Public License
18  * along with this program. If not, see <https://www.gnu.org/licenses/>.
19  ***************************************************************************/
20 
21 #ifndef THUNDEREGG_GMG_CYCLE_H
22 #define THUNDEREGG_GMG_CYCLE_H
23 
29 #include <ThunderEgg/GMG/Level.h>
30 #include <ThunderEgg/Vector.h>
31 #include <list>
32 
33 namespace ThunderEgg::GMG {
39 template<int D>
40 class Cycle : public Operator<D>
41 {
42 private:
46  std::shared_ptr<const Level<D>> finest_level;
47 
48 protected:
57  Vector<D> restrict(const Level<D>& level, const Vector<D>& f, const Vector<D>& u) const
58  {
59  // calculate residual
60  Vector<D> r = u.getZeroClone();
61  level.getOperator().apply(u, r);
62  r.scaleThenAdd(-1, f);
63  // create vectors for coarser levels
64  return level.getRestrictor().restrict(r);
65  }
66 
74  virtual void visit(const Level<D>& level, const Vector<D>& f, Vector<D>& u) const = 0;
75 
76 public:
82  Cycle(const Level<D>& finest_level)
83  : finest_level(new Level<D>(finest_level))
84  {}
94  void apply(const Vector<D>& f, Vector<D>& u) const
95  {
96  u.setWithGhost(0);
97  visit(*finest_level, f, u);
98  }
104  const Level<D>& getFinestLevel() const { return *finest_level; }
105 };
106 extern template class Cycle<2>;
107 extern template class Cycle<3>;
108 } // namespace ThunderEgg::GMG
109 #endif
ThunderEgg::GMG::Cycle::visit
virtual void visit(const Level< D > &level, const Vector< D > &f, Vector< D > &u) const =0
Virtual visit function that needs to be implemented in derived classes.
ThunderEgg::Vector::scaleThenAdd
void scaleThenAdd(double alpha, const Vector< D > &b)
this = alpha * this + b
Definition: Vector.h:509
Vector.h
Vector class.
ThunderEgg::GMG::Cycle::getFinestLevel
const Level< D > & getFinestLevel() const
Get the finest Level.
Definition: Cycle.h:104
ThunderEgg::GMG::Level::getOperator
const Operator< D > & getOperator() const
Get the operator for this level.
Definition: Level.h:120
Level.h
Level class.
ThunderEgg::GMG::Cycle::apply
void apply(const Vector< D > &f, Vector< D > &u) const
Run one iteration of the cycle.
Definition: Cycle.h:94
ThunderEgg::GMG::Cycle::Cycle
Cycle(const Level< D > &finest_level)
Create new cycle object.
Definition: Cycle.h:82
ThunderEgg::GMG::Level::getRestrictor
const Restrictor< D > & getRestrictor() const
Get the restriction operator for this level.
Definition: Level.h:81
ThunderEgg::Vector::setWithGhost
void setWithGhost(double alpha)
set all values in the vector (including ghost cells)
Definition: Vector.h:404
ThunderEgg::GMG::Cycle
Base abstract class for cycles.
Definition: Cycle.h:40
ThunderEgg::Operator
Base class for operators.
Definition: Operator.h:37
ThunderEgg::GMG::Level
Represents a level in geometric multi-grid.
Definition: Level.h:38
ThunderEgg::Vector
Vector class for use in thunderegg.
Definition: Vector.h:42
ThunderEgg::Vector::getZeroClone
Vector< D > getZeroClone() const
Get a vector of the same length initialized to zero.
Definition: Vector.h:601
ThunderEgg::GMG::Cycle::restrict
Vector< D > restrict(const Level< D > &level, const Vector< D > &f, const Vector< D > &u) const
Prepare vectors for coarser level.
Definition: Cycle.h:57
ThunderEgg::GMG
Geometric-Multigrid classes.
Definition: Cycle.h:33