ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_Interpolation_WENO_Z.H
Go to the documentation of this file.
1 #ifndef ERF_INTERPOLATE_WENO_Z_H_
2 #define ERF_INTERPOLATE_WENO_Z_H_
3 
4 #include "ERF_DataStruct.H"
5 
6 /**
7  * Interpolation operators used for WENO_Z-3 scheme
8  */
9 struct WENO_Z3
10 {
11  WENO_Z3 (const amrex::Array4<const amrex::Real>& phi,
12  const amrex::Real /*upw_frac*/)
13  : m_phi(phi) {}
14 
15  AMREX_GPU_DEVICE
16  AMREX_FORCE_INLINE
17  void
18  InterpolateInX (const int& i,
19  const int& j,
20  const int& k,
21  const int& qty_index,
22  amrex::Real& val_lo,
23  amrex::Real upw_lo) const
24  {
25  // Upwinding flags
26  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
27 
28  // Data to interpolate on
29  amrex::Real sp1 = m_phi(i+1, j , k , qty_index);
30  amrex::Real s = m_phi(i , j , k , qty_index);
31  amrex::Real sm1 = m_phi(i-1, j , k , qty_index);
32  amrex::Real sm2 = m_phi(i-2, j , k , qty_index);
33 
34  // Left and right fluxes
35  amrex::Real Fhi = Evaluate(sm2,sm1,s );
36  amrex::Real Flo = Evaluate(sp1,s ,sm1);
37 
38  // Numerical flux
39  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
40  }
41 
42  AMREX_GPU_DEVICE
43  AMREX_FORCE_INLINE
44  void
45  InterpolateInY (const int& i,
46  const int& j,
47  const int& k,
48  const int& qty_index,
49  amrex::Real& val_lo,
50  amrex::Real upw_lo) const
51  {
52  // Upwinding flags
53  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
54 
55  // Data to interpolate on
56  amrex::Real sp1 = m_phi(i , j+1, k , qty_index);
57  amrex::Real s = m_phi(i , j , k , qty_index);
58  amrex::Real sm1 = m_phi(i , j-1, k , qty_index);
59  amrex::Real sm2 = m_phi(i , j-2, k , qty_index);
60 
61  // Left and right fluxes
62  amrex::Real Fhi = Evaluate(sm2,sm1,s );
63  amrex::Real Flo = Evaluate(sp1,s ,sm1);
64 
65  // Numerical flux
66  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
67  }
68 
69  AMREX_GPU_DEVICE
70  AMREX_FORCE_INLINE
71  void
72  InterpolateInZ (const int& i,
73  const int& j,
74  const int& k,
75  const int& qty_index,
76  amrex::Real& val_lo,
77  amrex::Real upw_lo) const
78  {
79  // Upwinding flags
80  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
81 
82  // Data to interpolate on
83  amrex::Real sp1 = m_phi(i , j , k+1, qty_index);
84  amrex::Real s = m_phi(i , j , k , qty_index);
85  amrex::Real sm1 = m_phi(i , j , k-1, qty_index);
86  amrex::Real sm2 = m_phi(i , j , k-2, qty_index);
87 
88  // Left and right fluxes
89  amrex::Real Fhi = Evaluate(sm2,sm1,s );
90  amrex::Real Flo = Evaluate(sp1,s ,sm1);
91 
92  // Numerical flux
93  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
94  }
95 
96  AMREX_GPU_DEVICE
97  AMREX_FORCE_INLINE
98  amrex::Real
99  Evaluate (const amrex::Real& sm1,
100  const amrex::Real& s ,
101  const amrex::Real& sp1) const
102  {
103  // Smoothing factors
104  amrex::Real b1 = (s - sm1) * (s - sm1);
105  amrex::Real b2 = (sp1 - s) * (sp1 - s);
106 
107  // Weight factors
108  amrex::Real t5 = std::abs(b2 - b1);
109  amrex::Real w1 = g1 * ( amrex::Real(1) + (t5*t5) / ((eps + b1) * (eps + b1)) );
110  amrex::Real w2 = g2 * ( amrex::Real(1) + (t5*t5) / ((eps + b2) * (eps + b2)) );
111 
112  // Weight factor norm
113  amrex::Real wsum = w1 + w2;
114 
115  // Taylor expansions
116  amrex::Real v1 = -sm1 + amrex::Real(3) * s;
117  amrex::Real v2 = s + sp1;
118 
119  // Interpolated value
120  return ( (w1 * v1 + w2 * v2) / (amrex::Real(2) * wsum) );
121  }
122 
123 private:
124  amrex::Array4<const amrex::Real> m_phi; // Quantity to interpolate
125  const amrex::Real eps=amrex::Real(1.0e-40);
126  static constexpr amrex::Real g1=(amrex::Real(1)/amrex::Real(3));
127  static constexpr amrex::Real g2=(amrex::Real(2)/amrex::Real(3));
128 };
129 
130 /**
131  * Interpolation operators used for WENO_MZQ3 scheme
132  */
133 struct WENO_MZQ3
134 {
135  WENO_MZQ3 (const amrex::Array4<const amrex::Real>& phi,
136  const amrex::Real /*upw_frac*/)
137  : m_phi(phi) {}
138 
139  AMREX_GPU_DEVICE
140  AMREX_FORCE_INLINE
141  void
142  InterpolateInX (const int& i,
143  const int& j,
144  const int& k,
145  const int& qty_index,
146  amrex::Real& val_lo,
147  amrex::Real upw_lo) const
148  {
149  // Upwinding flags
150  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
151 
152  // Data to interpolate on
153  amrex::Real sp1 = m_phi(i+1, j , k , qty_index);
154  amrex::Real s = m_phi(i , j , k , qty_index);
155  amrex::Real sm1 = m_phi(i-1, j , k , qty_index);
156  amrex::Real sm2 = m_phi(i-2, j , k , qty_index);
157 
158  // Left and right fluxes
159  amrex::Real Fhi = Evaluate(sm2,sm1,s );
160  amrex::Real Flo = Evaluate(sp1,s ,sm1);
161 
162  // Numerical flux
163  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
164  }
165 
166  AMREX_GPU_DEVICE
167  AMREX_FORCE_INLINE
168  void
169  InterpolateInY (const int& i,
170  const int& j,
171  const int& k,
172  const int& qty_index,
173  amrex::Real& val_lo,
174  amrex::Real upw_lo) const
175  {
176  // Upwinding flags
177  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
178 
179  // Data to interpolate on
180  amrex::Real sp1 = m_phi(i , j+1, k , qty_index);
181  amrex::Real s = m_phi(i , j , k , qty_index);
182  amrex::Real sm1 = m_phi(i , j-1, k , qty_index);
183  amrex::Real sm2 = m_phi(i , j-2, k , qty_index);
184 
185  // Left and right fluxes
186  amrex::Real Fhi = Evaluate(sm2,sm1,s );
187  amrex::Real Flo = Evaluate(sp1,s ,sm1);
188 
189  // Numerical flux
190  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
191  }
192 
193  AMREX_GPU_DEVICE
194  AMREX_FORCE_INLINE
195  void
196  InterpolateInZ (const int& i,
197  const int& j,
198  const int& k,
199  const int& qty_index,
200  amrex::Real& val_lo,
201  amrex::Real upw_lo) const
202  {
203  // Upwinding flags
204  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
205 
206  // Data to interpolate on
207  amrex::Real sp1 = m_phi(i , j , k+1, qty_index);
208  amrex::Real s = m_phi(i , j , k , qty_index);
209  amrex::Real sm1 = m_phi(i , j , k-1, qty_index);
210  amrex::Real sm2 = m_phi(i , j , k-2, qty_index);
211 
212  // Left and right fluxes
213  amrex::Real Fhi = Evaluate(sm2,sm1,s );
214  amrex::Real Flo = Evaluate(sp1,s ,sm1);
215 
216  // Numerical flux
217  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
218  }
219 
220  AMREX_GPU_DEVICE
221  AMREX_FORCE_INLINE
222  amrex::Real
223  Evaluate (const amrex::Real& sm1,
224  const amrex::Real& s ,
225  const amrex::Real& sp1) const
226  {
227  // Smoothing factors
228  amrex::Real b1 = (s - sm1) * (s - sm1);
229  amrex::Real b2 = (sp1 - s) * (sp1 - s);
230  amrex::Real b3 = ( (amrex::Real(13.0) / amrex::Real(12.0)) * ((sm1 - amrex::Real(2)*s + sp1)*(sm1 - amrex::Real(2)*s + sp1)) ) + ( ((sp1 - sm1)*(sp1 - sm1)) / amrex::Real(4.0) );
231 
232  // Weight factors
233  amrex::Real t5 = ( std::abs(b3 - b1) + std::abs(b3 - b2) ) / amrex::Real(32.0);
234  amrex::Real a1 = g1 * ( amrex::Real(1) + (t5*t5) / ((eps + b1) * (eps + b1)) );
235  amrex::Real a2 = g2 * ( amrex::Real(1) + (t5*t5) / ((eps + b2) * (eps + b2)) );
236  amrex::Real a3 = g3 * ( amrex::Real(1) + (t5*t5) / ((eps + b3) * (eps + b3)) );
237  amrex::Real asum = a1 + a2 + a3;
238  amrex::Real w1 = a1 / asum;
239  amrex::Real w2 = a2 / asum;
240  amrex::Real w3 = a3 / asum;
241 
242  // Taylor expansions
243  amrex::Real v1 = (-sm1 + amrex::Real(3) * s) / amrex::Real(2);
244  amrex::Real v2 = ( s + sp1) / amrex::Real(2);
245  amrex::Real v3 = (-sm1 + amrex::Real(5.0) * s + amrex::Real(2) * sp1) / amrex::Real(6.0);
246 
247  // Interpolated value
248  return ( (w3 / g3) * (v3 - g1 * v1 - g2 * v2) + w1 * v1 + w2 * v2 );
249  }
250 
251 private:
252  amrex::Array4<const amrex::Real> m_phi; // Quantity to interpolate
253  const amrex::Real eps=amrex::Real(1.0e-40);
254  static constexpr amrex::Real g1=(amrex::Real(1)/amrex::Real(3));
255  static constexpr amrex::Real g2=(amrex::Real(1)/amrex::Real(3));
256  static constexpr amrex::Real g3=(amrex::Real(1)/amrex::Real(3));
257 };
258 
259 /**
260  * Interpolation operators used for WENO_Z-5 scheme
261  */
262 struct WENO_Z5
263 {
264  WENO_Z5 (const amrex::Array4<const amrex::Real>& phi,
265  const amrex::Real /*upw_frac*/)
266  : m_phi(phi) {}
267 
268  AMREX_GPU_DEVICE
269  AMREX_FORCE_INLINE
270  void
271  InterpolateInX (const int& i,
272  const int& j,
273  const int& k,
274  const int& qty_index,
275  amrex::Real& val_lo,
276  amrex::Real upw_lo) const
277  {
278  // Upwinding flags
279  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
280 
281  // Data to interpolate on
282  amrex::Real sp2 = m_phi(i+2, j , k , qty_index);
283  amrex::Real sp1 = m_phi(i+1, j , k , qty_index);
284  amrex::Real s = m_phi(i , j , k , qty_index);
285  amrex::Real sm1 = m_phi(i-1, j , k , qty_index);
286  amrex::Real sm2 = m_phi(i-2, j , k , qty_index);
287  amrex::Real sm3 = m_phi(i-3, j , k , qty_index);
288 
289  // Left and right fluxes
290  amrex::Real Fhi = Evaluate(sm3,sm2,sm1,s ,sp1);
291  amrex::Real Flo = Evaluate(sp2,sp1,s ,sm1,sm2);
292 
293  // Numerical flux
294  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
295  }
296 
297  AMREX_GPU_DEVICE
298  AMREX_FORCE_INLINE
299  void
300  InterpolateInY (const int& i,
301  const int& j,
302  const int& k,
303  const int& qty_index,
304  amrex::Real& val_lo,
305  amrex::Real upw_lo) const
306  {
307  // Upwinding flags
308  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
309 
310  // Data to interpolate on
311  amrex::Real sp2 = m_phi(i , j+2, k , qty_index);
312  amrex::Real sp1 = m_phi(i , j+1, k , qty_index);
313  amrex::Real s = m_phi(i , j , k , qty_index);
314  amrex::Real sm1 = m_phi(i , j-1, k , qty_index);
315  amrex::Real sm2 = m_phi(i , j-2, k , qty_index);
316  amrex::Real sm3 = m_phi(i , j-3, k , qty_index);
317 
318  // Left and right fluxes
319  amrex::Real Fhi = Evaluate(sm3,sm2,sm1,s ,sp1);
320  amrex::Real Flo = Evaluate(sp2,sp1,s ,sm1,sm2);
321 
322  // Numerical flux
323  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
324  }
325 
326  AMREX_GPU_DEVICE
327  AMREX_FORCE_INLINE
328  void
329  InterpolateInZ (const int& i,
330  const int& j,
331  const int& k,
332  const int& qty_index,
333  amrex::Real& val_lo,
334  amrex::Real upw_lo) const
335  {
336  // Upwinding flags
337  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
338 
339  // Data to interpolate on
340  amrex::Real sp2 = m_phi(i , j , k+2, qty_index);
341  amrex::Real sp1 = m_phi(i , j , k+1, qty_index);
342  amrex::Real s = m_phi(i , j , k , qty_index);
343  amrex::Real sm1 = m_phi(i , j , k-1, qty_index);
344  amrex::Real sm2 = m_phi(i , j , k-2, qty_index);
345  amrex::Real sm3 = m_phi(i , j , k-3, qty_index);
346 
347  // Left and right fluxes
348  amrex::Real Fhi = Evaluate(sm3,sm2,sm1,s ,sp1);
349  amrex::Real Flo = Evaluate(sp2,sp1,s ,sm1,sm2);
350 
351  // Numerical flux
352  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
353  }
354 
355  AMREX_GPU_DEVICE
356  AMREX_FORCE_INLINE
357  amrex::Real
358  Evaluate (const amrex::Real& sm2,
359  const amrex::Real& sm1,
360  const amrex::Real& s ,
361  const amrex::Real& sp1,
362  const amrex::Real& sp2) const
363  {
364  // Smoothing factors
365  amrex::Real b1 = c1 * (sm2 - amrex::Real(2) * sm1 + s) * (sm2 - amrex::Real(2) * sm1 + s) +
366  amrex::Real(0.25) * (sm2 - amrex::Real(4.0) * sm1 + amrex::Real(3) * s) * (sm2 - amrex::Real(4.0) * sm1 + amrex::Real(3) * s);
367  amrex::Real b2 = c1 * (sm1 - amrex::Real(2) * s + sp1) * (sm1 - amrex::Real(2) * s + sp1) +
368  amrex::Real(0.25) * (sm1 - sp1) * (sm1 - sp1);
369  amrex::Real b3 = c1 * (s - amrex::Real(2) * sp1 + sp2) * (s - amrex::Real(2) * sp1 + sp2) +
370  amrex::Real(0.25) * (amrex::Real(3) * s - amrex::Real(4.0) * sp1 + sp2) * (amrex::Real(3) * s - amrex::Real(4.0) * sp1 + sp2);
371 
372  // Weight factors
373  amrex::Real t5 = std::abs(b3 - b1);
374  amrex::Real w1 = g1 * ( amrex::Real(1) + (t5*t5) / ((eps + b1) * (eps + b1)) );
375  amrex::Real w2 = g2 * ( amrex::Real(1) + (t5*t5) / ((eps + b2) * (eps + b2)) );
376  amrex::Real w3 = g3 * ( amrex::Real(1) + (t5*t5) / ((eps + b3) * (eps + b3)) );
377 
378  // Weight factor norm
379  amrex::Real wsum = w1 + w2 + w3;
380 
381  // Taylor expansions
382  amrex::Real v1 = amrex::Real(2) * sm2 - amrex::Real(7.0) * sm1 + amrex::Real(11.0) * s;
383  amrex::Real v2 = -sm1 + amrex::Real(5.0) * s + amrex::Real(2) * sp1;
384  amrex::Real v3 = amrex::Real(2) * s + amrex::Real(5.0) * sp1 - sp2;
385 
386  // Interpolated value
387  return ( (w1 * v1 + w2 * v2 + w3 * v3) / (amrex::Real(6.0) * wsum) );
388  }
389 
390 private:
391  amrex::Array4<const amrex::Real> m_phi; // Quantity to interpolate
392  const amrex::Real eps=amrex::Real(1.0e-40);
393  static constexpr amrex::Real c1=(amrex::Real(13.0)/amrex::Real(12.0));
394  static constexpr amrex::Real g1=(amrex::Real(1)/amrex::Real(10.0));
395  static constexpr amrex::Real g2=(amrex::Real(3)/amrex::Real(5.0));
396  static constexpr amrex::Real g3=(amrex::Real(3)/amrex::Real(10.0));
397 };
398 
399 /**
400  * Interpolation operators used for WENO_Z-7 scheme
401  */
402 struct WENO_Z7
403 {
404  WENO_Z7 (const amrex::Array4<const amrex::Real>& phi,
405  const amrex::Real /*upw_frac*/)
406  : m_phi(phi) {}
407 
408  AMREX_GPU_DEVICE
409  AMREX_FORCE_INLINE
410  void
411  InterpolateInX (const int& i,
412  const int& j,
413  const int& k,
414  const int& qty_index,
415  amrex::Real& val_lo,
416  amrex::Real upw_lo) const
417  {
418  // Upwinding flags
419  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
420 
421  // Data to interpolate on
422  amrex::Real sp3 = m_phi(i+3, j , k , qty_index);
423  amrex::Real sp2 = m_phi(i+2, j , k , qty_index);
424  amrex::Real sp1 = m_phi(i+1, j , k , qty_index);
425  amrex::Real s = m_phi(i , j , k , qty_index);
426  amrex::Real sm1 = m_phi(i-1, j , k , qty_index);
427  amrex::Real sm2 = m_phi(i-2, j , k , qty_index);
428  amrex::Real sm3 = m_phi(i-3, j , k , qty_index);
429  amrex::Real sm4 = m_phi(i-4, j , k , qty_index);
430 
431  // Left and right fluxes
432  amrex::Real Fhi = Evaluate(sm4,sm3,sm2,sm1,s ,sp1,sp2);
433  amrex::Real Flo = Evaluate(sp3,sp2,sp1,s ,sm1,sm2,sm3);
434 
435  // Numerical flux
436  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
437  }
438 
439  AMREX_GPU_DEVICE
440  AMREX_FORCE_INLINE
441  void
442  InterpolateInY (const int& i,
443  const int& j,
444  const int& k,
445  const int& qty_index,
446  amrex::Real& val_lo,
447  amrex::Real upw_lo) const
448  {
449  // Upwinding flags
450  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
451 
452  // Data to interpolate on
453  amrex::Real sp3 = m_phi(i , j+3, k , qty_index);
454  amrex::Real sp2 = m_phi(i , j+2, k , qty_index);
455  amrex::Real sp1 = m_phi(i , j+1, k , qty_index);
456  amrex::Real s = m_phi(i , j , k , qty_index);
457  amrex::Real sm1 = m_phi(i , j-1, k , qty_index);
458  amrex::Real sm2 = m_phi(i , j-2, k , qty_index);
459  amrex::Real sm3 = m_phi(i , j-3, k , qty_index);
460  amrex::Real sm4 = m_phi(i , j-4, k , qty_index);
461 
462  // Left and right fluxes
463  amrex::Real Fhi = Evaluate(sm4,sm3,sm2,sm1,s ,sp1,sp2);
464  amrex::Real Flo = Evaluate(sp3,sp2,sp1,s ,sm1,sm2,sm3);
465 
466  // Numerical flux
467  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
468  }
469 
470  AMREX_GPU_DEVICE
471  AMREX_FORCE_INLINE
472  void
473  InterpolateInZ (const int& i,
474  const int& j,
475  const int& k,
476  const int& qty_index,
477  amrex::Real& val_lo,
478  amrex::Real upw_lo) const
479  {
480  // Upwinding flags
481  if (upw_lo != amrex::Real(0)) upw_lo = (upw_lo > 0) ? amrex::Real(1) : -amrex::Real(1);
482 
483  // Data to interpolate on
484  amrex::Real sp3 = m_phi(i , j , k+3, qty_index);
485  amrex::Real sp2 = m_phi(i , j , k+2, qty_index);
486  amrex::Real sp1 = m_phi(i , j , k+1, qty_index);
487  amrex::Real s = m_phi(i , j , k , qty_index);
488  amrex::Real sm1 = m_phi(i , j , k-1, qty_index);
489  amrex::Real sm2 = m_phi(i , j , k-2, qty_index);
490  amrex::Real sm3 = m_phi(i , j , k-3, qty_index);
491  amrex::Real sm4 = m_phi(i , j , k-4, qty_index);
492 
493  // Left and right fluxes
494  amrex::Real Fhi = Evaluate(sm4,sm3,sm2,sm1,s ,sp1,sp2);
495  amrex::Real Flo = Evaluate(sp3,sp2,sp1,s ,sm1,sm2,sm3);
496 
497  // Numerical flux
498  val_lo = (amrex::Real(1) + upw_lo)/amrex::Real(2) * Fhi + (amrex::Real(1) - upw_lo)/amrex::Real(2) * Flo;
499  }
500 
501  AMREX_GPU_DEVICE
502  AMREX_FORCE_INLINE
503  amrex::Real
504  Evaluate (const amrex::Real& sm3,
505  const amrex::Real& sm2,
506  const amrex::Real& sm1,
507  const amrex::Real& s ,
508  const amrex::Real& sp1,
509  const amrex::Real& sp2,
510  const amrex::Real& sp3) const
511  {
512  // Smoothing factors
513  amrex::Real b1 = ( sm3 * sm3 * amrex::Real(6649.)/amrex::Real(2880.0)
514  - sm3 * sm2 * amrex::Real(2623.)/amrex::Real(160.0)
515  + sm3 * sm1 * amrex::Real(9449.)/amrex::Real(480.0)
516  - sm3 * s * amrex::Real(11389.)/amrex::Real(1440.0)
517  + sm2 * sm2 * amrex::Real(28547.)/amrex::Real(960.0)
518  - sm2 * sm1 * amrex::Real(35047.)/amrex::Real(480.0)
519  + sm2 * s * amrex::Real(14369.)/amrex::Real(480.0)
520  + sm1 * sm1 * amrex::Real(44747.)/amrex::Real(960.0)
521  - sm1 * s * amrex::Real(6383.)/amrex::Real(160.0)
522  + s * s * amrex::Real(25729.)/amrex::Real(2880.0) );
523  amrex::Real b2 = ( sm2 * sm2 * 3169/amrex::Real(2880.0)
524  - sm2 * sm1 * 3229/amrex::Real(480.0)
525  + sm2 * s * 3169/amrex::Real(480.0)
526  - sm2 * sp1 * 2989/amrex::Real(1440.0)
527  + sm1 * sm1 * 11147/amrex::Real(960.0)
528  - sm1 * s * 11767/amrex::Real(480.0)
529  + sm1 * sp1 * 1283/amrex::Real(160.0)
530  + s * s * 13667/amrex::Real(960.0)
531  - s * sp1 * 5069/amrex::Real(480.0)
532  + sp1 * sp1 * 6649/amrex::Real(2880.0) );
533  amrex::Real b3 = ( sm1 * sm1 * amrex::Real(6649.)/amrex::Real(2880.0)
534  - sm1 * s * amrex::Real(5069.)/amrex::Real(480.0)
535  + sm1 * sp1 * amrex::Real(1283.)/amrex::Real(160.0)
536  - sm1 * sp2 * amrex::Real(2989.)/amrex::Real(1440.0)
537  + s * s * amrex::Real(13667.)/amrex::Real(960.0)
538  - s * sp1 * amrex::Real(11767.)/amrex::Real(480.0)
539  + s * sp2 * amrex::Real(3169.)/amrex::Real(480.0)
540  + sp1 * sp1 * amrex::Real(11147.)/amrex::Real(960.0)
541  - sp1 * sp2 * amrex::Real(3229.)/amrex::Real(480.0)
542  + sp2 * sp2 * amrex::Real(3169.)/amrex::Real(2880.0) );
543  amrex::Real b4 = ( s * s * amrex::Real(25729.)/amrex::Real(2880.)
544  - s * sp1 * amrex::Real(6383.)/amrex::Real(160.)
545  + s * sp2 * amrex::Real(14369.)/amrex::Real(480.)
546  - s * sp3 * amrex::Real(11389.)/amrex::Real(1440.)
547  + sp1 * sp1 * amrex::Real(44747.)/amrex::Real(960.)
548  - sp1 * sp2 * amrex::Real(35047.)/amrex::Real(480.)
549  + sp1 * sp3 * amrex::Real(9449.)/amrex::Real(480.)
550  + sp2 * sp2 * amrex::Real(28547.)/amrex::Real(960.)
551  - sp2 * sp3 * amrex::Real(2623.)/amrex::Real(160.)
552  + sp3 * sp3 * amrex::Real(6649.)/amrex::Real(2880.) );
553 
554  // Weight factors
555  amrex::Real t5 = std::abs(b1 - b2 - b3 + b4);
556  amrex::Real w1 = g1 * ( amrex::Real(1) + (t5*t5) / ((eps + b1) * (eps + b1)) );
557  amrex::Real w2 = g2 * ( amrex::Real(1) + (t5*t5) / ((eps + b2) * (eps + b2)) );
558  amrex::Real w3 = g3 * ( amrex::Real(1) + (t5*t5) / ((eps + b3) * (eps + b3)) );
559  amrex::Real w4 = g4 * ( amrex::Real(1) + (t5*t5) / ((eps + b4) * (eps + b4)) );
560 
561  // Weight factor norm
562  amrex::Real wsum = w1 + w2 + w3 + w4;
563 
564  // Taylor expansions
565  amrex::Real v1 = (-amrex::Real(0.3125))*sm3 + ( amrex::Real(1.3125))*sm2 + (-amrex::Real(2.1875))*sm1 + ( amrex::Real(2.1875))*s;
566  amrex::Real v2 = ( amrex::Real(0.0625))*sm2 + (-amrex::Real(0.3125))*sm1 + ( amrex::Real(0.9375))*s + ( amrex::Real(0.3125))*sp1;
567  amrex::Real v3 = (-amrex::Real(0.0625))*sm1 + ( amrex::Real(0.5625))*s + ( amrex::Real(0.5625))*sp1 + (-amrex::Real(0.0625))*sp2;
568  amrex::Real v4 = ( amrex::Real(0.3125))*s + ( amrex::Real(0.9375))*sp1 + (-amrex::Real(0.3125))*sp2 + ( amrex::Real(0.0625))*sp3;
569 
570  // Interpolated value
571  return ( (w1 * v1 + w2 * v2 + w3 * v3 + w4 * v4) / (wsum) );
572  }
573 
574 private:
575  amrex::Array4<const amrex::Real> m_phi; // Quantity to interpolate
576  const amrex::Real eps=amrex::Real(1.0e-40);
577  static constexpr amrex::Real g1=( amrex::Real(1)/amrex::Real(64.0));
578  static constexpr amrex::Real g2=(amrex::Real(21.0)/amrex::Real(64.0));
579  static constexpr amrex::Real g3=(amrex::Real(35.0)/amrex::Real(64.0));
580  static constexpr amrex::Real g4=( amrex::Real(7.0)/amrex::Real(64.0));
581 };
582 #endif
ERF_DataStruct.H
Real
amrex::Real Real
Definition: ERF_ShocInterface.H:19
module_model_constants::a2
real(c_double), parameter a2
Definition: ERF_module_model_constants.F90:95
module_model_constants::a3
real(c_double), parameter a3
Definition: ERF_module_model_constants.F90:96
WENO_MZQ3
Definition: ERF_Interpolation_WENO_Z.H:134
WENO_MZQ3::g2
static constexpr amrex::Real g2
Definition: ERF_Interpolation_WENO_Z.H:255
WENO_MZQ3::m_phi
amrex::Array4< const amrex::Real > m_phi
Definition: ERF_Interpolation_WENO_Z.H:252
WENO_MZQ3::InterpolateInZ
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInZ(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:196
WENO_MZQ3::Evaluate
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Evaluate(const amrex::Real &sm1, const amrex::Real &s, const amrex::Real &sp1) const
Definition: ERF_Interpolation_WENO_Z.H:223
WENO_MZQ3::g3
static constexpr amrex::Real g3
Definition: ERF_Interpolation_WENO_Z.H:256
WENO_MZQ3::WENO_MZQ3
WENO_MZQ3(const amrex::Array4< const amrex::Real > &phi, const amrex::Real)
Definition: ERF_Interpolation_WENO_Z.H:135
WENO_MZQ3::eps
const amrex::Real eps
Definition: ERF_Interpolation_WENO_Z.H:253
WENO_MZQ3::g1
static constexpr amrex::Real g1
Definition: ERF_Interpolation_WENO_Z.H:254
WENO_MZQ3::InterpolateInX
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInX(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:142
WENO_MZQ3::InterpolateInY
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInY(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:169
WENO_Z3
Definition: ERF_Interpolation_WENO_Z.H:10
WENO_Z3::eps
const amrex::Real eps
Definition: ERF_Interpolation_WENO_Z.H:125
WENO_Z3::InterpolateInZ
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInZ(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:72
WENO_Z3::m_phi
amrex::Array4< const amrex::Real > m_phi
Definition: ERF_Interpolation_WENO_Z.H:124
WENO_Z3::WENO_Z3
WENO_Z3(const amrex::Array4< const amrex::Real > &phi, const amrex::Real)
Definition: ERF_Interpolation_WENO_Z.H:11
WENO_Z3::Evaluate
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Evaluate(const amrex::Real &sm1, const amrex::Real &s, const amrex::Real &sp1) const
Definition: ERF_Interpolation_WENO_Z.H:99
WENO_Z3::g2
static constexpr amrex::Real g2
Definition: ERF_Interpolation_WENO_Z.H:127
WENO_Z3::InterpolateInX
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInX(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:18
WENO_Z3::g1
static constexpr amrex::Real g1
Definition: ERF_Interpolation_WENO_Z.H:126
WENO_Z3::InterpolateInY
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInY(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:45
WENO_Z5
Definition: ERF_Interpolation_WENO_Z.H:263
WENO_Z5::Evaluate
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Evaluate(const amrex::Real &sm2, const amrex::Real &sm1, const amrex::Real &s, const amrex::Real &sp1, const amrex::Real &sp2) const
Definition: ERF_Interpolation_WENO_Z.H:358
WENO_Z5::InterpolateInZ
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInZ(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:329
WENO_Z5::g3
static constexpr amrex::Real g3
Definition: ERF_Interpolation_WENO_Z.H:396
WENO_Z5::InterpolateInY
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInY(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:300
WENO_Z5::InterpolateInX
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInX(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:271
WENO_Z5::g1
static constexpr amrex::Real g1
Definition: ERF_Interpolation_WENO_Z.H:394
WENO_Z5::m_phi
amrex::Array4< const amrex::Real > m_phi
Definition: ERF_Interpolation_WENO_Z.H:391
WENO_Z5::c1
static constexpr amrex::Real c1
Definition: ERF_Interpolation_WENO_Z.H:393
WENO_Z5::eps
const amrex::Real eps
Definition: ERF_Interpolation_WENO_Z.H:392
WENO_Z5::g2
static constexpr amrex::Real g2
Definition: ERF_Interpolation_WENO_Z.H:395
WENO_Z5::WENO_Z5
WENO_Z5(const amrex::Array4< const amrex::Real > &phi, const amrex::Real)
Definition: ERF_Interpolation_WENO_Z.H:264
WENO_Z7
Definition: ERF_Interpolation_WENO_Z.H:403
WENO_Z7::g1
static constexpr amrex::Real g1
Definition: ERF_Interpolation_WENO_Z.H:577
WENO_Z7::Evaluate
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Evaluate(const amrex::Real &sm3, const amrex::Real &sm2, const amrex::Real &sm1, const amrex::Real &s, const amrex::Real &sp1, const amrex::Real &sp2, const amrex::Real &sp3) const
Definition: ERF_Interpolation_WENO_Z.H:504
WENO_Z7::InterpolateInX
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInX(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:411
WENO_Z7::m_phi
amrex::Array4< const amrex::Real > m_phi
Definition: ERF_Interpolation_WENO_Z.H:575
WENO_Z7::eps
const amrex::Real eps
Definition: ERF_Interpolation_WENO_Z.H:576
WENO_Z7::InterpolateInZ
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInZ(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:473
WENO_Z7::WENO_Z7
WENO_Z7(const amrex::Array4< const amrex::Real > &phi, const amrex::Real)
Definition: ERF_Interpolation_WENO_Z.H:404
WENO_Z7::g4
static constexpr amrex::Real g4
Definition: ERF_Interpolation_WENO_Z.H:580
WENO_Z7::g2
static constexpr amrex::Real g2
Definition: ERF_Interpolation_WENO_Z.H:578
WENO_Z7::g3
static constexpr amrex::Real g3
Definition: ERF_Interpolation_WENO_Z.H:579
WENO_Z7::InterpolateInY
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void InterpolateInY(const int &i, const int &j, const int &k, const int &qty_index, amrex::Real &val_lo, amrex::Real upw_lo) const
Definition: ERF_Interpolation_WENO_Z.H:442