ThunderEgg  1.0.0
WCycle.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_WCYCLE_H
22 #define THUNDEREGG_GMG_WCYCLE_H
23 
28 #include <ThunderEgg/GMG/Cycle.h>
30 namespace ThunderEgg::GMG {
34 template<int D>
35 class WCycle : public Cycle<D>
36 {
37 private:
38  int num_pre_sweeps = 1;
39  int num_post_sweeps = 1;
40  int num_coarse_sweeps = 1;
41  int num_mid_sweeps = 1;
42 
43 protected:
50  void visit(const Level<D>& level, const Vector<D>& f, Vector<D>& u) const override
51  {
52  if (level.coarsest()) {
53  for (int i = 0; i < num_coarse_sweeps; i++) {
54  level.getSmoother().smooth(f, u);
55  }
56  } else {
57  for (int i = 0; i < num_pre_sweeps; i++) {
58  level.getSmoother().smooth(f, u);
59  }
60 
61  Vector<D> coarser_f = this->restrict(level, f, u);
62 
63  const Level<D>& coarser_level = level.getCoarser();
64  Vector<D> coarser_u = coarser_f.getZeroClone();
65 
66  this->visit(coarser_level, coarser_f, coarser_u);
67 
68  coarser_level.getInterpolator().interpolate(coarser_u, u);
69 
70  for (int i = 0; i < num_mid_sweeps; i++) {
71  level.getSmoother().smooth(f, u);
72  }
73 
74  coarser_f = this->restrict(level, f, u);
75 
76  this->visit(coarser_level, coarser_f, coarser_u);
77 
78  coarser_level.getInterpolator().interpolate(coarser_u, u);
79 
80  for (int i = 0; i < num_post_sweeps; i++) {
81  level.getSmoother().smooth(f, u);
82  }
83  }
84  }
85 
86 public:
92  WCycle(const Level<D>& finest_level, const CycleOpts& opts)
93  : Cycle<D>(finest_level)
94  {
95  num_pre_sweeps = opts.pre_sweeps;
96  num_post_sweeps = opts.post_sweeps;
97  num_coarse_sweeps = opts.coarse_sweeps;
98  num_mid_sweeps = opts.mid_sweeps;
99  }
105  WCycle<D>* clone() const override { return new WCycle(*this); }
106 };
107 extern template class WCycle<2>;
108 extern template class WCycle<3>;
109 } // namespace ThunderEgg::GMG
110 #endif
ThunderEgg::GMG::Level::getSmoother
const Smoother< D > & getSmoother() const
Get smoother operator for this level.
Definition: Level.h:138
ThunderEgg::GMG::Level::getInterpolator
const Interpolator< D > & getInterpolator() const
Get the interpolation operator for this level.
Definition: Level.h:102
ThunderEgg::GMG::WCycle::clone
WCycle< D > * clone() const override
Get a clone of this WCycle.
Definition: WCycle.h:105
Cycle.h
Cycle class.
ThunderEgg::GMG::CycleOpts::mid_sweeps
int mid_sweeps
Number of sweeps inbetween up and down.
Definition: CycleOpts.h:46
ThunderEgg::GMG::WCycle::visit
void visit(const Level< D > &level, const Vector< D > &f, Vector< D > &u) const override
Implements W-cycle. Pre-smooth, visit coarser level, smooth, visit coarse level, and then post-smooth...
Definition: WCycle.h:50
ThunderEgg::GMG::CycleOpts::post_sweeps
int post_sweeps
Number of sweeps on up cycle.
Definition: CycleOpts.h:42
CycleOpts.h
CycleOpts struct.
ThunderEgg::GMG::Level::getCoarser
const Level & getCoarser() const
get reference to the coarser level.
Definition: Level.h:156
ThunderEgg::GMG::Cycle
Base abstract class for cycles.
Definition: Cycle.h:40
ThunderEgg::GMG::CycleOpts::pre_sweeps
int pre_sweeps
Number of sweeps on down cycle.
Definition: CycleOpts.h:38
ThunderEgg::GMG::WCycle
Implementation of a W-cycle.
Definition: WCycle.h:35
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::GMG::CycleOpts
Options for Cycle classes.
Definition: CycleOpts.h:33
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
ThunderEgg::GMG::WCycle::WCycle
WCycle(const Level< D > &finest_level, const CycleOpts &opts)
Create new W-cycle.
Definition: WCycle.h:92
ThunderEgg::GMG::CycleOpts::coarse_sweeps
int coarse_sweeps
Number of sweeps on coarse level.
Definition: CycleOpts.h:50
ThunderEgg::GMG::Level::coarsest
bool coarsest() const
Check if this level is the coarsest level.
Definition: Level.h:174