Chombo + EB + MF  3.2
AMRLevelAdvectDiffuse.H
Go to the documentation of this file.
1 #ifdef CH_LANG_CC
2 /*
3  * _______ __
4  * / ___/ / ___ __ _ / / ___
5  * / /__/ _ \/ _ \/ V \/ _ \/ _ \
6  * \___/_//_/\___/_/_/_/_.__/\___/
7  * Please refer to Copyright.txt, in Chombo's root directory.
8  */
9 #endif
10 
11 #ifndef _AMRLEVELADVECTDIFFUSE_H_
12 #define _AMRLEVELADVECTDIFFUSE_H_
13 
14 #include "FArrayBox.H"
15 #include "LevelData.H"
16 #include "AMRLevel.H"
17 #include "CoarseAverage.H"
18 #include "FineInterp.H"
19 #include "LevelFluxRegister.H"
20 
21 #include "LevelAdvect.H"
22 #include "PhysIBC.H"
23 #include "ProblemDomain.H"
24 #include "IntVectSet.H"
25 #include "Vector.H"
26 #include "RealVect.H"
27 #include "LevelFluxRegister.H"
28 #include "DisjointBoxLayout.H"
29 #include "AdvectionFunctions.H"
30 #include "BCFunc.H"
31 #include "LevelTGA.H"
32 #include "AMRPoissonOp.H"
33 #include "BiCGStabSolver.H"
34 #include "RelaxSolver.H"
35 
36 #include "NamespaceHeader.H"
37 
38 /// AMRLevel for advection-diffusion
39 /**
40  */
42 {
43 public:
44 
45  /// Default constructor
47  {
48  m_isDefined = false;
49  }
50 
51  /// Full constructor. Arguments are same as in define()
53  AdvectionVelocityFunction a_advFunc,
54  BCHolder a_bcFunc,
55  const Real& a_cfl,
56  const Real& a_domainLength,
57  const Real& a_refineThresh,
58  const int& a_tagBufferSize,
59  const Real& a_initialDtMultiplier,
60  const bool& a_useLimiting,
61  const Real& a_nu)
62  {
63  define(a_gphys, a_advFunc, a_bcFunc, a_cfl, a_domainLength, a_refineThresh,
64  a_tagBufferSize, a_initialDtMultiplier, a_useLimiting, a_nu);
65  }
66 
67  void makeDiffusiveSource(LevelData<FArrayBox>& a_diffusiveSource);
68 
69  /// Defines this AMRLevelAdvectDiffuse
70  void define(/// advection physics class
71  const AdvectPhysics& a_gphys,
72  /// velocity function providing the advection velocity
73  AdvectionVelocityFunction a_advFunc,
74  /// boundary condition class for diffusion solve
75  BCHolder a_bcFunc,
76  /// CFL number
77  const Real& a_cfl,
78  /// physical length of domain
79  const Real& a_domainLength,
80  /// undivided gradient size over which a cell will be tagged for refinement
81  const Real& a_refineThresh,
82  /// number of buffer cells around each tagged cell that will also be tagged
83  const int& a_tagBufferSize,
84  /// CFL number at beginning of calculation
85  const Real& a_initialDtMultiplier,
86  /// whether to use van Leer limiting
87  const bool& a_useLimiting,
88  /// diffusion coefficient
89  const Real& a_nu);
90 
91  ///
92  virtual ~AMRLevelAdvectDiffuse();
93 
94  /// Never called: historical
95  virtual void define(AMRLevel* a_coarserLevelPtr,
96  const Box& a_problemDomain,
97  int a_level,
98  int a_refRatio)
99  {
100  MayDay::Error("never called--historical");
101  }
102 
103  /// Define new AMRLevelAdvectDiffuse from coarser
104  /**
105  */
106  virtual void define(
107  AMRLevel* a_coarserLevelPtr,
108  const ProblemDomain& a_problemDomain,
109  int a_level,
110  int a_refRatio);
111 
112  /// Advance by one timestep
113  virtual Real advance();
114 
115  /// Things to do after a timestep
116  virtual void postTimeStep();
117 
118  /// Create tags for regridding
119  virtual void tagCells(IntVectSet& a_tags) ;
120 
121  /// Create tags at initialization
122  virtual void tagCellsInit(IntVectSet& a_tags) ;
123 
124  /// Set up data on this level after regridding
125  virtual void regrid(const Vector<Box>& a_newGrids);
126 
127  /// Initialize grids
128  virtual void initialGrid(const Vector<Box>& a_newGrids);
129 
130  /// Initialize data
131  virtual void initialData();
132 
133  /// Things to do after initialization
134  virtual void postInitialize();
135 
136 #ifdef CH_USE_HDF5
137  /// Write checkpoint header
138  virtual void writeCheckpointHeader(HDF5Handle& a_handle) const;
139 
140  /// Write checkpoint data for this level
141  virtual void writeCheckpointLevel(HDF5Handle& a_handle) const;
142 
143  /// Read checkpoint header
144  virtual void readCheckpointHeader(HDF5Handle& a_handle);
145 
146  /// Read checkpoint data for this level
147  virtual void readCheckpointLevel(HDF5Handle& a_handle);
148 
149  /// Write plotfile header
150  virtual void writePlotHeader(HDF5Handle& a_handle) const;
151 
152  /// Write plotfile data for this level
153  virtual void writePlotLevel(HDF5Handle& a_handle) const;
154 #endif
155 
156  /// Returns the dt computed earlier for this level
157  virtual Real computeDt();
158 
159  /// Compute dt using initial data
160  virtual Real computeInitialDt();
161 
162  //for convergence tests
164  {
165  return m_UNew;
166  }
167 
168  //for convergence tests
170  {
171  return m_UOld;
172  }
173 protected:
174  void setSolverCoef(Real a_alpha, Real a_beta);
176  const LevelData<FArrayBox>& a_old,
177  const LevelData<FArrayBox>& a_new,
178  Real a_time, Real a_tOld, Real a_tNew);
179 
181  Vector<DisjointBoxLayout>& a_grids,
182  Vector<int>& a_refRat,
183  ProblemDomain& a_lev0Dom,
184  Real& a_lev0Dx);
185 
186  void doImplicitReflux();
187 
188  void printRefluxRHSMax(const std::string& a_preflix) ;
189 
190  Real diffusiveAdvance(LevelData<FArrayBox>& a_diffusiveSource);
191  void defineSolvers();
193  // solver to do level diffusion and flux register interaction
195  // one amrmultigrid to rule them all and in darkness bind them
197  // operator factory
199  // bottom solver
201  // static RelaxSolver<LevelData<FArrayBox> > s_botSolver;
202 
203  //embarrassing but needed function
204  void getCoarseDataPointers(LevelData<FArrayBox>** a_coarserDataOldPtr,
205  LevelData<FArrayBox>** a_coarserDataNewPtr,
206  LevelFluxRegister** a_coarserFRPtr,
207  LevelFluxRegister** a_finerFRPtr,
208  Real& a_tCoarserOld,
209  Real& a_tCoarserNew);
210 
211  void fillAdvectionVelocity();
212  // Setup menagerie of data structures
213  void levelSetup();
214 
215  // Get the next coarser level
217 
218  // Get the next finer level
220 
221  // Conserved state, U, at old and new time
224 
237 
244 
249 
251 
252 private:
253 
254  // Disallowed for all the usual reasons
255  void operator=(const AMRLevelAdvectDiffuse&);
257 };
258 
259 #include "NamespaceFooter.H"
260 
261 #endif
void interpolateInTime(LevelData< FArrayBox > &a_interp, const LevelData< FArrayBox > &a_old, const LevelData< FArrayBox > &a_new, Real a_time, Real a_tOld, Real a_tNew)
bool m_isDefined
Definition: AMRLevelAdvectDiffuse.H:226
AMRLevelAdvectDiffuse(const AdvectPhysics &a_gphys, AdvectionVelocityFunction a_advFunc, BCHolder a_bcFunc, const Real &a_cfl, const Real &a_domainLength, const Real &a_refineThresh, const int &a_tagBufferSize, const Real &a_initialDtMultiplier, const bool &a_useLimiting, const Real &a_nu)
Full constructor. Arguments are same as in define()
Definition: AMRLevelAdvectDiffuse.H:52
replaces coarse level data with an average of fine level data.
Definition: CoarseAverage.H:30
bool m_doImplicitReflux
Definition: AMRLevelAdvectDiffuse.H:247
AMRLevel for advection-diffusion.
Definition: AMRLevelAdvectDiffuse.H:41
void define(const AdvectPhysics &a_gphys, AdvectionVelocityFunction a_advFunc, BCHolder a_bcFunc, const Real &a_cfl, const Real &a_domainLength, const Real &a_refineThresh, const int &a_tagBufferSize, const Real &a_initialDtMultiplier, const bool &a_useLimiting, const Real &a_nu)
Defines this AMRLevelAdvectDiffuse.
Real(* AdvectionVelocityFunction)(const RealVect &a_point, const int &a_velComp)
Velocity function interface.
Definition: AdvectionFunctions.H:23
An irregular domain on an integer lattice.
Definition: IntVectSet.H:44
Definition: BCFunc.H:136
A class to facilitate interaction with physical boundary conditions.
Definition: ProblemDomain.H:141
LevelData< FArrayBox > & getStateNew()
Definition: AMRLevelAdvectDiffuse.H:163
LevelData< FArrayBox > m_dU
Definition: AMRLevelAdvectDiffuse.H:222
Real m_dtNew
Definition: AMRLevelAdvectDiffuse.H:240
LevelData< FArrayBox > & getStateOld()
Definition: AMRLevelAdvectDiffuse.H:169
void getHierarchyAndGrids(Vector< AMRLevelAdvectDiffuse *> &a_hierarchy, Vector< DisjointBoxLayout > &a_grids, Vector< int > &a_refRat, ProblemDomain &a_lev0Dom, Real &a_lev0Dx)
Abstract base class for time-dependent data at a level of refinement.
Definition: AMRLevel.H:47
virtual void regrid(const Vector< Box > &a_newGrids)
Set up data on this level after regridding.
virtual ~AMRLevelAdvectDiffuse()
Real m_initialDtMultiplier
Definition: AMRLevelAdvectDiffuse.H:231
static RefCountedPtr< AMRLevelOpFactory< LevelData< FArrayBox > > > s_diffuseOpFact
Definition: AMRLevelAdvectDiffuse.H:198
Advection integrator on a level.
Definition: LevelAdvect.H:30
Real m_domainLength
Definition: AMRLevelAdvectDiffuse.H:228
virtual void writeCheckpointHeader(HDF5Handle &a_handle) const
Write checkpoint header.
replaces fine level data with interpolation of coarse level data.
Definition: FineInterp.H:32
LevelFluxRegister m_fluxRegister
Definition: AMRLevelAdvectDiffuse.H:243
virtual void initialGrid(const Vector< Box > &a_newGrids)
Initialize grids.
static RefCountedPtr< AMRMultiGrid< LevelData< FArrayBox > > > s_diffuseAMRMG
Definition: AMRLevelAdvectDiffuse.H:196
bool m_hasFiner
Definition: AMRLevelAdvectDiffuse.H:246
virtual void initialData()
Initialize data.
Real m_dx
Definition: AMRLevelAdvectDiffuse.H:234
Vector< string > m_stateNames
Definition: AMRLevelAdvectDiffuse.H:225
double Real
Definition: REAL.H:33
void setSolverCoef(Real a_alpha, Real a_beta)
Real diffusiveAdvance(LevelData< FArrayBox > &a_diffusiveSource)
A BoxLayout that has a concept of disjointedness.
Definition: DisjointBoxLayout.H:30
LevelData< FArrayBox > m_UNew
Definition: AMRLevelAdvectDiffuse.H:222
virtual Real advance()
Advance by one timestep.
static void Error(const char *const a_msg=m_nullString, int m_exitCode=CH_DEFAULT_ERROR_CODE)
Print out message to cerr and exit with the specified exit code.
virtual void writeCheckpointLevel(HDF5Handle &a_handle) const
Write checkpoint data for this level.
AdvectionVelocityFunction m_advFunc
Definition: AMRLevelAdvectDiffuse.H:236
AMRLevelAdvectDiffuse * getFinerLevel() const
LevelFluxRegister-A class to encapsulate a levels worth of flux registers.
Definition: LevelFluxRegister.H:29
RefCountedPtr< AdvectPhysics > m_advPhys
Definition: AMRLevelAdvectDiffuse.H:235
A class derived from GodunovPhysics for simple advection-diffusion problems.
Definition: AdvectPhysics.H:22
void printRefluxRHSMax(const std::string &a_preflix)
A Rectangular Domain on an Integer Lattice.
Definition: Box.H:469
virtual void postTimeStep()
Things to do after a timestep.
virtual void writePlotHeader(HDF5Handle &a_handle) const
Write plotfile header.
virtual Real computeInitialDt()
Compute dt using initial data.
int m_tagBufferSize
Definition: AMRLevelAdvectDiffuse.H:230
Real m_refineThresh
Definition: AMRLevelAdvectDiffuse.H:229
BCHolder m_bcFunc
Definition: AMRLevelAdvectDiffuse.H:192
virtual void define(AMRLevel *a_coarserLevelPtr, const Box &a_problemDomain, int a_level, int a_refRatio)
Never called: historical.
Definition: AMRLevelAdvectDiffuse.H:95
virtual void tagCells(IntVectSet &a_tags)
Create tags for regridding.
virtual void writePlotLevel(HDF5Handle &a_handle) const
Write plotfile data for this level.
Real m_nu
Definition: AMRLevelAdvectDiffuse.H:233
Handle to a particular group in an HDF file.
Definition: CH_HDF5.H:294
AMRLevelAdvectDiffuse * getCoarserLevel() const
int m_numGhost
Definition: AMRLevelAdvectDiffuse.H:241
FineInterp m_fineInterp
Definition: AMRLevelAdvectDiffuse.H:238
virtual void postInitialize()
Things to do after initialization.
virtual void readCheckpointLevel(HDF5Handle &a_handle)
Read checkpoint data for this level.
AMRLevelAdvectDiffuse()
Default constructor.
Definition: AMRLevelAdvectDiffuse.H:46
virtual void tagCellsInit(IntVectSet &a_tags)
Create tags at initialization.
virtual void readCheckpointHeader(HDF5Handle &a_handle)
Read checkpoint header.
static RefCountedPtr< LevelTGA > s_diffuseLevTGA
Definition: AMRLevelAdvectDiffuse.H:194
CoarseAverage m_coarseAverage
Definition: AMRLevelAdvectDiffuse.H:239
LevelData< FArrayBox > m_UOld
Definition: AMRLevelAdvectDiffuse.H:222
void operator=(const AMRLevelAdvectDiffuse &)
virtual Real computeDt()
Returns the dt computed earlier for this level.
void getCoarseDataPointers(LevelData< FArrayBox > **a_coarserDataOldPtr, LevelData< FArrayBox > **a_coarserDataNewPtr, LevelFluxRegister **a_coarserFRPtr, LevelFluxRegister **a_finerFRPtr, Real &a_tCoarserOld, Real &a_tCoarserNew)
LevelData< FluxBox > m_advVel
Definition: AMRLevelAdvectDiffuse.H:223
bool m_hasDiffusion
Definition: AMRLevelAdvectDiffuse.H:248
static BiCGStabSolver< LevelData< FArrayBox > > s_botSolver
Definition: AMRLevelAdvectDiffuse.H:200
bool m_hasCoarser
Definition: AMRLevelAdvectDiffuse.H:245
LevelAdvect m_levelGodunov
Definition: AMRLevelAdvectDiffuse.H:242
Real m_cfl
Definition: AMRLevelAdvectDiffuse.H:227
void makeDiffusiveSource(LevelData< FArrayBox > &a_diffusiveSource)
bool m_useLimiting
Definition: AMRLevelAdvectDiffuse.H:232
Definition: BiCGStabSolver.H:26
DisjointBoxLayout m_grids
Definition: AMRLevelAdvectDiffuse.H:250