Chombo + EB + MF  3.2
MOLUtilities.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 _MOLUTILITIES_H_
12 #define _MOLUTILITIES_H_
13 
14 #include "Box.H"
15 #include "IntVectSet.H"
16 #include "Vector.H"
17 #include "REAL.H"
18 #include "FArrayBox.H"
19 #include "FluxBox.H"
20 #include "NamespaceHeader.H"
21 
22 class MOLPhysics;
23 
24 ///
25 /**
26  Utility class for higher-order methods of lines: slopes, parabolic
27  interpolants, limiters. Contains no physics-dependent methods, but
28  one of the member functions (PPMFaceLimiter) may require a pointer
29  to a physics class analysis class in order to perform limiting in
30  characteristic variables.
31 */
33 {
34 public:
35  /// Constructor
36  /**
37  */
38  MOLUtilities();
39 
40  /// Destructor
41  /**
42  */
43  ~MOLUtilities();
44 
45  /// Define the object
46  /**
47  */
48  void define(const ProblemDomain& a_domain,
49  const Real& a_dx);
50 
51  /// Compute the slope flattening coefficients
52  /**
53  Compute the slope flattening coefficients, a_flattening, using the
54  primitive variables, a_W, within a_box.
55  */
56  void computeFlattening(/// flattening coeffs, 1 component on a_box
57  FArrayBox& a_flattening,
58  /// primitive variables, on a_box grown by 3 within m_domain
59  const FArrayBox& a_W,
60  /// interval of a_W with velocity; length SpaceDim
61  const Interval& a_velInt,
62  /// index of a_W with pressure
63  const int& a_presInd,
64  /// minimum pressure to use in ratio
65  const Real& a_smallPres,
66  /// index of a_W with bulk modulus
67  const int& a_bulkModulusInd,
68  /// cell-centered box on which a_flattening lives
69  const Box& a_box);
70 
71  /// Apply the flattening to slopes.
72  /**
73  Multiply every component of a_dW by a_flat, on a_box.
74  */
75  void applyFlattening(/// slopes to be flattened, on at least a_box
76  FArrayBox& a_dW,
77  /// cell-centered flattening coefficients, on at least a_box
78  const FArrayBox& a_flat,
79  /// cell-centered box on which to flatten
80  const Box& a_box);
81 
82  /// Compute componentwise van Leer slopes.
83  /**
84  Given cell averages W, compute van Leer slopes dW on a
85  component-by-component basis.
86  */
87  void vanLeerSlopes(FArrayBox& a_dW,
88  const FArrayBox& a_W,
89  const int& a_numSlopes,
90  const bool& a_useLimiting,
91  const int& a_dir,
92  const Box& a_box);
93 
94  /// Compute fourth-order slopes.
95  /**
96  Given cell averages W and van Leer slopes dWvL, compute fourth-order
97  slopes dW4. Limiting is performed in a separate pass.
98  */
99  void fourthOrderSlopes(FArrayBox& a_dW4,
100  const FArrayBox& a_W,
101  const FArrayBox& a_dWvL,
102  const int& a_numSlopes,
103  const int& a_dir,
104  const Box& a_box);
105 
106  /// Compute slopes (dW- and dW+) using one sided differences
107  /**
108  */
109  void oneSidedDifferences(FArrayBox& a_dWMinus,
110  FArrayBox& a_dWPlus,
111  const FArrayBox& a_W,
112  const int& a_dir,
113  const Box& a_box);
114 
115  /// Compute slopes (dW (center), dW-, and dW+)
116  /**
117  a_dwCent, a_dwMinus, a_dWPlus, and a_dW all live on
118  cell-centered a_entireBox.
119 
120  For i in a_centerBox:
121  a_dwCent[i] = (a_W[i+1] - a_W[i-1])/2;
122  a_dWMinus[i] = a_W[i] - a_W[i-1];
123  a_dWPlus[i] = a_W[i+1] - a_W[i].
124 
125  For i in a_loBox, set only a_dWPlus[i] = a_W[i+1] - a_W[i]
126  and copy it to a_dWCent[i] and a_dWMinus[i].
127 
128  For i in a_hiBox, set only a_dWMinus[i] = a_W[i] - a_W[i-1]
129  and copy it to a_dWCent[i] and a_dWPlus[i].
130  */
131  void slopes(FArrayBox& a_dWCent,
132  FArrayBox& a_dWMinus,
133  FArrayBox& a_dWPlus,
134  const FArrayBox& a_W,
135  const int& a_numSlopes,
136  const int& a_dir,
137  const Box& a_loBox,
138  const Box& a_hiBox,
139  const Box& a_centerBox,
140  const Box& a_entireBox,
141  const int& a_hasLo,
142  const int& a_hasHi);
143 
144  /// van Leer slope limiter.
145  /**
146  On input, dW contains the centered, unlimited slopes, and
147  dW(Minus,Plus) contain the one-sided slopes from the minus, plus sides.
148  On output, dW contains the limited slopes.
149  slopes dW4. Limiting is performed in a separate pass.
150  */
151  void slopeLimiter(FArrayBox& a_dW,
152  const FArrayBox& a_dWMinus,
153  const FArrayBox& a_dWPlus,
154  const int& a_numSlopes,
155  const Box& a_box);
156 
157  /// extremum-preserving van Leer slope limiter.
158  /**
159  On input, dW contains the centered, unlimited slopes, and
160  dW(Minus,Plus) contain the one-sided slopes from the minus, plus sides.
161  On output, dW contains the limited slopes.
162  slopes dW4. Limiting is performed in a separate pass.
163  */
165  const FArrayBox& a_dWMinus,
166  const FArrayBox& a_dWPlus,
167  const int& a_numSlopes,
168  const Box& a_box,
169  const int& a_dir);
170 
171  /// PPM face-centered interpolant.
172  /**
173  Given the cell average a_W, compute fourth-order accurate
174  face-centered values WFace on a_box by differentiating the indefinite
175  integral. Limiting is performed in a separate pass.
176  */
177  void PPMFaceValues(FArrayBox& a_WFace,
178  const FArrayBox& a_W,
179  const int& a_numSlopes,
180  const bool& a_useLimiting,
181  const int& a_dir,
182  const Box& a_box,
183  const Real& a_time,
184  const MOLPhysics* a_physPtr = NULL);
185 
186  /// PPM Limiter.
187  /**
188  a_dWMinus and a_dWPlus are the differences between the face values
189  on the minus and plus sides of cells and the average in the cell.
190  That is,
191  a_dWMinus[i] = WFace[i - e/2] - a_W[i]
192  a_dWPlus[i] = WFace[i + e/2] - a_W[i]
193  where e is the unit vector in dimension a_dir.
194  The PPM limiter is applied to these values to obtain a monotone
195  interpolant in the cell.
196  The function returns the limited a_dWMinus and a_dWPlus on a_box.
197  petermc, 4 Sep 2008: included a_W in argument list
198 
199  If m_highOrderLimiter,
200  then need a_dWMinus and a_dWPlus on a_box,
201  and need a_W on on a_box grown by 3 in dimension a_dir.
202  Returns limited a_dWMinus and a_dWPlus on a_box.
203  */
204  void PPMLimiter(FArrayBox& a_dWMinus,
205  FArrayBox& a_dWPlus,
206  const FArrayBox& a_W,
207  const int& a_numSlopes,
208  const int& a_dir,
209  const Box& a_box);
210 
211  /// Compute face-centered velocity divergence.
212  /**
213  Returns face-centered velocity divergence on for faces in the
214  direction a_dir. The velocities are the components a_velInterval
215  of a_W.
216  */
217  void divVel(FArrayBox& a_divVel,
218  const FArrayBox& a_W,
219  const Interval& a_velInt,
220  const int& a_dir,
221  const Box& a_box);
222 
223  void divVelHO(FArrayBox& a_divVel,
224  const FArrayBox& a_W,
225  const int& a_dir,
226  const Box& a_box,
227  MOLPhysics* a_physPtr);
228 
229  /// Compute face-centered artificial viscosity flux.
230  /**
231  Increments face-centered flux in the a_dir direction with quadratic
232  artificial viscosity.
233 
234  Increments a_F on all non-boundary a_dir-faces of a_box
235  using a_U on all cells of a_box
236  and a_divVel on all non-boundary a_dir-faces of a_box.
237  */
238  void artificialViscosity(FArrayBox& a_F,
239  const FArrayBox& a_U,
240  const FArrayBox& a_divVel,
241  const Real& a_scale,
242  const int& a_dir,
243  const Box& a_box);
244 
245  /// Set whether to use high-order limiter.
246  void highOrderLimiter(bool a_highOrderLimiter);
247 
248  /// a_avgFab += a_sign * laplacian(a_cellFab) * m_dx^2 / 24, on a_cellFab.box()
249  void deconvolve(FArrayBox& a_avgFab,
250  const FArrayBox& a_cellFab,
251  int a_sign = 1);
252 
253  /// a_avgFab += a_sign * laplacian(a_cellFab) * m_dx^2 / 24, on a_box
254  void deconvolve(FArrayBox& a_avgFab,
255  const FArrayBox& a_cellFab,
256  const Box& a_box,
257  int a_sign = 1);
258 
259  /// a_avgFlux += a_sign * laplacian(a_cellFlux) * m_dx^2 / 24, on faces of a_cellFlux.box()
260  void deconvolveFace(FluxBox& a_avgFlux,
261  const FluxBox& a_cellFlux,
262  int a_sign = 1);
263 
264  /// a_avgFlux += a_sign * laplacian(a_cellFlux) * m_dx^2 / 24, on faces of a_box
265  void deconvolveFace(FluxBox& a_avgFlux,
266  const FluxBox& a_cellFlux,
267  const Box& a_box,
268  int a_sign = 1);
269 
270  /// a_faceFab += a_sign * laplacian(a_faceExtFab) * m_dx^2 / 24, on faces of a_box, which is a_dir-face-centered
271  /**
272  Convert a_faceFab from face-centered to face-averaged (a_sign = 1)
273  or vice versa (a_sign = -1) on a_dir-face-centered a_box
274  by using the second derivative from a_faceExtFab, which will
275  typically extend farther than a_faceFab does.
276  */
277  void deconvolveCenterFace(FArrayBox& a_faceFab,
278  const FArrayBox& a_faceExtFab,
279  const Box& a_box,
280  int a_dir,
281  int a_sign = 1);
282 
283 protected:
284  // Problem domain and grid spacing
287 
288  // Has this object been defined
290 
291  // Use a high-order limiter? (default false)
293 
294 private:
295  // We may allow copying and assignment later.
296  // Disallowed for all the usual reasons
297  void operator=(const MOLUtilities& a_input)
298  {
299  MayDay::Error("invalid operator");
300  }
301 
302  // Disallowed for all the usual reasons
303  MOLUtilities(const MOLUtilities& a_input)
304  {
305  MayDay::Error("invalid operator");
306  }
307 };
308 
309 #include "NamespaceFooter.H"
310 #endif
void oneSidedDifferences(FArrayBox &a_dWMinus, FArrayBox &a_dWPlus, const FArrayBox &a_W, const int &a_dir, const Box &a_box)
Compute slopes (dW- and dW+) using one sided differences.
A class to facilitate interaction with physical boundary conditions.
Definition: ProblemDomain.H:141
void divVel(FArrayBox &a_divVel, const FArrayBox &a_W, const Interval &a_velInt, const int &a_dir, const Box &a_box)
Compute face-centered velocity divergence.
bool m_highOrderLimiter
Definition: MOLUtilities.H:292
void highOrderLimiter(bool a_highOrderLimiter)
Set whether to use high-order limiter.
Real m_dx
Definition: MOLUtilities.H:286
Definition: MOLUtilities.H:32
Definition: MOLPhysics.H:37
void applyFlattening(FArrayBox &a_dW, const FArrayBox &a_flat, const Box &a_box)
Apply the flattening to slopes.
void slopes(FArrayBox &a_dWCent, FArrayBox &a_dWMinus, FArrayBox &a_dWPlus, const FArrayBox &a_W, const int &a_numSlopes, const int &a_dir, const Box &a_loBox, const Box &a_hiBox, const Box &a_centerBox, const Box &a_entireBox, const int &a_hasLo, const int &a_hasHi)
Compute slopes (dW (center), dW-, and dW+)
ProblemDomain m_domain
Definition: MOLUtilities.H:285
void fourthOrderSlopes(FArrayBox &a_dW4, const FArrayBox &a_W, const FArrayBox &a_dWvL, const int &a_numSlopes, const int &a_dir, const Box &a_box)
Compute fourth-order slopes.
~MOLUtilities()
Destructor.
A FArrayBox-like container for face-centered fluxes.
Definition: FluxBox.H:22
void slopeLimiter(FArrayBox &a_dW, const FArrayBox &a_dWMinus, const FArrayBox &a_dWPlus, const int &a_numSlopes, const Box &a_box)
van Leer slope limiter.
Structure for passing component ranges in code.
Definition: Interval.H:23
double Real
Definition: REAL.H:33
void vanLeerSlopes(FArrayBox &a_dW, const FArrayBox &a_W, const int &a_numSlopes, const bool &a_useLimiting, const int &a_dir, const Box &a_box)
Compute componentwise van Leer slopes.
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.
void deconvolveFace(FluxBox &a_avgFlux, const FluxBox &a_cellFlux, int a_sign=1)
a_avgFlux += a_sign * laplacian(a_cellFlux) * m_dx^2 / 24, on faces of a_cellFlux.box()
bool m_isDefined
Definition: MOLUtilities.H:289
void define(const ProblemDomain &a_domain, const Real &a_dx)
Define the object.
A Rectangular Domain on an Integer Lattice.
Definition: Box.H:469
void artificialViscosity(FArrayBox &a_F, const FArrayBox &a_U, const FArrayBox &a_divVel, const Real &a_scale, const int &a_dir, const Box &a_box)
Compute face-centered artificial viscosity flux.
void computeFlattening(FArrayBox &a_flattening, const FArrayBox &a_W, const Interval &a_velInt, const int &a_presInd, const Real &a_smallPres, const int &a_bulkModulusInd, const Box &a_box)
Compute the slope flattening coefficients.
Definition: FArrayBox.H:45
MOLUtilities()
Constructor.
void divVelHO(FArrayBox &a_divVel, const FArrayBox &a_W, const int &a_dir, const Box &a_box, MOLPhysics *a_physPtr)
void slopeLimiterExtPreserving(FArrayBox &a_dW, const FArrayBox &a_dWMinus, const FArrayBox &a_dWPlus, const int &a_numSlopes, const Box &a_box, const int &a_dir)
extremum-preserving van Leer slope limiter.
void deconvolveCenterFace(FArrayBox &a_faceFab, const FArrayBox &a_faceExtFab, const Box &a_box, int a_dir, int a_sign=1)
a_faceFab += a_sign * laplacian(a_faceExtFab) * m_dx^2 / 24, on faces of a_box, which is a_dir-face-c...
MOLUtilities(const MOLUtilities &a_input)
Definition: MOLUtilities.H:303
void operator=(const MOLUtilities &a_input)
Definition: MOLUtilities.H:297
void PPMFaceValues(FArrayBox &a_WFace, const FArrayBox &a_W, const int &a_numSlopes, const bool &a_useLimiting, const int &a_dir, const Box &a_box, const Real &a_time, const MOLPhysics *a_physPtr=NULL)
PPM face-centered interpolant.
void PPMLimiter(FArrayBox &a_dWMinus, FArrayBox &a_dWPlus, const FArrayBox &a_W, const int &a_numSlopes, const int &a_dir, const Box &a_box)
PPM Limiter.
void deconvolve(FArrayBox &a_avgFab, const FArrayBox &a_cellFab, int a_sign=1)
a_avgFab += a_sign * laplacian(a_cellFab) * m_dx^2 / 24, on a_cellFab.box()