chombo-discharge
CD_CdrCTU.H
Go to the documentation of this file.
1 /* chombo-discharge
2  * Copyright © 2021 SINTEF Energy Research.
3  * Please refer to Copyright.txt and LICENSE in the chombo-discharge root directory.
4  */
5 
12 #ifndef CD_CdrCTU_H
13 #define CD_CdrCTU_H
14 
15 // Our includes
16 #include <CD_CdrMultigrid.H>
17 #include <CD_CdrGodunov.H>
18 #include <CD_NamespaceHeader.H>
19 
23 class CdrCTU : public CdrMultigrid
24 {
25 public:
29  CdrCTU();
30 
34  virtual ~CdrCTU();
35 
39  virtual void
40  parseOptions() override;
41 
45  virtual void
46  parseRuntimeOptions() override;
47 
54  virtual void
55  advectToFaces(EBAMRFluxData& a_facePhi, const EBAMRCellData& a_cellPhi, const Real a_dt) override;
56 
61  virtual Real
62  computeAdvectionDt() override;
63 
64 protected:
68  enum class Limiter
69  {
70  None,
71  MinMod,
72  Superbee,
73  MonotonizedCentral,
74  };
75 
80 
84  bool m_useCTU;
85 
95  virtual void
96  computeNormalSlopes(EBCellFAB& a_normalSlopes,
97  const EBCellFAB& a_cellPhi,
98  const Box& a_cellBox,
99  const ProblemDomain& a_domain,
100  const int a_level,
101  const DataIndex& a_dit);
102 
119  virtual void
120  upwind(EBFluxFAB& a_facePhi,
121  const EBCellFAB& a_normalSlopes,
122  const EBCellFAB& a_cellPhi,
123  const EBCellFAB& a_cellVel,
124  const EBFluxFAB& a_faceVel,
125  const ProblemDomain& a_domain,
126  const Box& a_cellBox,
127  const int& a_level,
128  const DataIndex& a_dit,
129  const Real& a_dt);
130 
134  virtual void
136 
142  Real
143  minmod(const Real& dwl, const Real& dwr) const noexcept;
144 
150  Real
151  superbee(const Real& dwl, const Real& dwr) const noexcept;
152 
158  Real
159  monotonizedCentral(const Real& dwl, const Real& dwr) const noexcept;
160 };
161 
162 #include <CD_NamespaceFooter.H>
163 
164 #endif
Declaration of a class which implements CdrMultigrid using Godunov-type advection.
Extensions of CdrSolver which use EBHelmholtzOp for diffusion solves.
Class that uses a slope limited method for advection, in combination with corner transport upwind (CT...
Definition: CD_CdrCTU.H:24
virtual void upwind(EBFluxFAB &a_facePhi, const EBCellFAB &a_normalSlopes, const EBCellFAB &a_cellPhi, const EBCellFAB &a_cellVel, const EBFluxFAB &a_faceVel, const ProblemDomain &a_domain, const Box &a_cellBox, const int &a_level, const DataIndex &a_dit, const Real &a_dt)
Upwind/Riemann solve. This extrapolates the cell-centered data to the inside face centers and selects...
Definition: CD_CdrCTU.cpp:460
Real superbee(const Real &dwl, const Real &dwr) const noexcept
Superbee slope limiter.
Definition: CD_CdrCTU.cpp:774
CdrCTU()
Default constructor.
Definition: CD_CdrCTU.cpp:22
virtual void parseOptions() override
Parse class options to put object in usable state.
Definition: CD_CdrCTU.cpp:39
Real minmod(const Real &dwl, const Real &dwr) const noexcept
minmod slope function.
Definition: CD_CdrCTU.cpp:762
Limiter
Supported limiters.
Definition: CD_CdrCTU.H:69
virtual void parseSlopeLimiter()
Parse slope limiting on/off.
Definition: CD_CdrCTU.cpp:152
Real monotonizedCentral(const Real &dwl, const Real &dwr) const noexcept
Monotonized central difference slope limiter.
Definition: CD_CdrCTU.cpp:791
virtual Real computeAdvectionDt() override
Compute the largest possible advective time step (for explicit methods)
Definition: CD_CdrCTU.cpp:73
virtual void computeNormalSlopes(EBCellFAB &a_normalSlopes, const EBCellFAB &a_cellPhi, const Box &a_cellBox, const ProblemDomain &a_domain, const int a_level, const DataIndex &a_dit)
Compute slopes.
Definition: CD_CdrCTU.cpp:246
virtual void parseRuntimeOptions() override
Parse runtime options.
Definition: CD_CdrCTU.cpp:56
virtual void advectToFaces(EBAMRFluxData &a_facePhi, const EBAMRCellData &a_cellPhi, const Real a_dt) override
MUSCL advection to faces.
Definition: CD_CdrCTU.cpp:183
bool m_useCTU
If true, the discretization uses CTU. Otherwise it's DCU (donor cell upwind)
Definition: CD_CdrCTU.H:84
Limiter m_limiter
If true, slopes are limited in the normal dircetion.
Definition: CD_CdrCTU.H:79
virtual ~CdrCTU()
Destructor (does nothing)
Definition: CD_CdrCTU.cpp:33
Extension class of CdrSolver that uses multigrid for diffusion part. Can also solve for stochastic di...
Definition: CD_CdrMultigrid.H:35