ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_EBMOSTStress.H
Go to the documentation of this file.
1 #ifndef ERF_EBMOSTStress_H
2 #define ERF_EBMOSTStress_H
3 
4 #include <ERF_Constants.H>
5 #include <ERF_IndexDefines.H>
6 
7 
8 /**
9  * Adiabatic with constant roughness
10  */
12 {
14  amrex::Real Qvflux)
15  {
16  mdata.surf_temp_flux = Tflux;
17  mdata.surf_moist_flux = Qvflux;
18  }
19 
20  AMREX_GPU_DEVICE
21  AMREX_FORCE_INLINE
22  void
23  iterate_flux (const int& i,
24  const int& j,
25  const int& k,
26  const int& /*max_iters*/,
27  const amrex::Array4<const amrex::Real>& zref_arr,
28  const amrex::Array4<const amrex::Real>& z0_arr,
29  const amrex::Array4<const amrex::Real>& umm_arr,
30  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
31  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
32  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
33  const amrex::Array4<amrex::Real>& u_star_arr,
34  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
35  const amrex::Array4<amrex::Real>& t_star_arr,
36  const amrex::Array4<amrex::Real>& q_star_arr,
37  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
38  const amrex::Array4<amrex::Real>& /*q_surf_arr*/,
39  const amrex::Array4<amrex::Real>& olen_arr,
40  const amrex::Array4<amrex::Real>& /*pblh_arr*/,
41  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
42  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
43  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
44  {
45  olen_arr(i,j,k) = bogus_large_value;
46  u_star_arr(i,j,k) = mdata.kappa * umm_arr(i,j,0) / std::log(zref_arr(i,j,0) / z0_arr(i,j,k));
47  t_star_arr(i,j,k) = zero;
48  q_star_arr(i,j,k) = zero;
49  }
50 
51 private:
54 };
55 
56 
57 /**
58  * Surface temperature with constant roughness
59  */
61 {
63  amrex::Real Qvflux,
64  bool cons_qflux)
65  {
66  mdata.surf_temp_flux = Tflux;
67  mdata.surf_moist_flux = Qvflux;
68  spec_qflux = cons_qflux;
69  }
70 
71  AMREX_GPU_DEVICE
72  AMREX_FORCE_INLINE
73  void
74  iterate_flux (const int& i,
75  const int& j,
76  const int& k,
77  const int& max_iters,
78  const amrex::Array4<const amrex::Real>& zref_arr,
79  const amrex::Array4<const amrex::Real>& z0_arr,
80  const amrex::Array4<const amrex::Real>& umm_arr,
81  const amrex::Array4<const amrex::Real>& tm_arr,
82  const amrex::Array4<const amrex::Real>& tvm_arr,
83  const amrex::Array4<const amrex::Real>& qvm_arr,
84  const amrex::Array4<amrex::Real>& u_star_arr,
85  const amrex::Array4<amrex::Real>& w_star_arr,
86  const amrex::Array4<amrex::Real>& t_star_arr,
87  const amrex::Array4<amrex::Real>& q_star_arr,
88  const amrex::Array4<amrex::Real>& t_surf_arr,
89  const amrex::Array4<amrex::Real>& q_surf_arr,
90  const amrex::Array4<amrex::Real>& olen_arr,
91  const amrex::Array4<amrex::Real>& pblh_arr,
92  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
93  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
94  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
95  {
96  amrex::Real Rib = zero;
97  amrex::Real zeta = zero;
98  amrex::Real zeta_old = zero;
99  amrex::Real psi_m = zero;
100  amrex::Real psi_h = zero;
101  amrex::Real num = zero;
102  amrex::Real den = zero;
103  amrex::Real zref = zref_arr(i,j,0);
104  amrex::Real z0 = z0_arr(i,j,k);
105  amrex::Real umm = std::max(umm_arr(i,j,0), WSMIN);
106  amrex::Real C = std::log(zref / z0);
107 
108  // First iteration we assume neutral (L -> inf)
109  if (u_star_arr(i,j,k) == bogus_large_value) {
110  olen_arr(i,j,k) = amrex::Real(1000);
111  }
112  zeta = zref / olen_arr(i,j,k);
113 
114  // Water vapor in atmos and surface
115  amrex::Real qv_s, qv_a;
116  if (q_surf_arr(i,j,k) > zero) {
117  qv_s = q_surf_arr(i,j,k);
118  } else {
119  // First iteration and no qv_surf was specified
120  // Use mean since there will be no flux
121  qv_s = qvm_arr(i,j,0);
122  }
123  qv_a = qvm_arr(i,j,0);
124 
125  // update w* and Umagmean from Beljaars (1995)
126  if (w_star_arr) {
127  // NOTE: Thv flux is lagged, similar to WRF
128  psi_m = sfuns.calc_psi_m2(zeta);
129  psi_h = sfuns.calc_psi_h2(zeta);
130  amrex::Real ustar = mdata.kappa * umm / (C - psi_m);
131  amrex::Real tstar = mdata.kappa * (tm_arr(i,j,0) - t_surf_arr(i,j,k)) / (C - psi_h);
133  -ustar * mdata.kappa * (qv_a - qv_s) / (C - psi_h);
134  amrex::Real tflux = -ustar*tstar*(1 + amrex::Real(0.61)*qvm_arr(i,j,0)) + amrex::Real(0.61)*tm_arr(i,j,0)*qflux;
135  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,0));
136  amrex::Real wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
137  umm = std::sqrt(umm_arr(i,j,0)*umm_arr(i,j,0) + wstar*wstar);
138  umm = std::max(umm, WSMIN);
139  }
140 
141  // Bulk Richardson number w/ moisture
142  amrex::Real thv_s = t_surf_arr(i,j,k) * (one + amrex::Real(0.61)*qv_s);
143  amrex::Real thv_a = tm_arr(i,j,0) * (one + amrex::Real(0.61)*qv_a);
144  Rib = ( (mdata.gravity * zref) / tm_arr(i,j,0) ) *
145  ( (thv_a - thv_s) / (umm * umm) );
146  Rib = std::min(std::max(Rib,-amrex::Real(4.0)),amrex::Real(4.0));
147 
148  // Fixed point iteration on zeta
149  int iter = 0;
150  do {
151  // Transfer curr to old
152  zeta_old = zeta;
153 
154  // Stability functions
155  psi_m = sfuns.calc_psi_m2(zeta_old);
156  psi_h = sfuns.calc_psi_h2(zeta_old);
157 
158  // Limiting
159  num = std::max(C - psi_m, amrex::Real(1.0));
160  den = std::max(C - psi_h, amrex::Real(1.0));
161 
162  // Update with under relaxation
163  zeta = (one - alpha) * zeta_old + alpha * Rib * num * num / den;
164 
165  ++iter;
166  } while ( (std::abs(zeta - zeta_old) > tol) && (iter <= max_iters) );
167  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
168  "Maximum number of MOST iterations reached.");
169 
170  // Populate the stored MOST arrays
171  olen_arr(i,j,k) = zref / zeta;
172  u_star_arr(i,j,k) = mdata.kappa * umm / (C - psi_m);
173  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,0) - t_surf_arr(i,j,k)) / (C - psi_h);
174  if (spec_qflux) {
175  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (C - psi_h) /
176  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,0);
177  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
178  } else {
179  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,0) - q_surf_arr(i,j,k)) / (C - psi_h);
180  }
181  }
182 
183 private:
187  const amrex::Real tol = amrex::Real(1.0e-3);
189  const amrex::Real WSMIN = amrex::Real(0.1); // minimum wind speed
190 };
191 
192 /**
193  * Surface flux with constant roughness
194  */
196 {
198  amrex::Real Qvflux,
199  bool cons_qflux)
200  {
201  mdata.surf_temp_flux = Tflux;
202  mdata.surf_moist_flux = Qvflux;
203  spec_qflux = cons_qflux;
204  }
205 
206  AMREX_GPU_DEVICE
207  AMREX_FORCE_INLINE
208  void
209  iterate_flux (const int& i,
210  const int& j,
211  const int& k,
212  const int& max_iters,
213  const amrex::Array4<const amrex::Real>& zref_arr,
214  const amrex::Array4<const amrex::Real>& z0_arr,
215  const amrex::Array4<const amrex::Real>& umm_arr,
216  const amrex::Array4<const amrex::Real>& tm_arr,
217  const amrex::Array4<const amrex::Real>& tvm_arr,
218  const amrex::Array4<const amrex::Real>& qvm_arr,
219  const amrex::Array4<amrex::Real>& u_star_arr,
220  const amrex::Array4<amrex::Real>& w_star_arr,
221  const amrex::Array4<amrex::Real>& t_star_arr,
222  const amrex::Array4<amrex::Real>& q_star_arr,
223  const amrex::Array4<amrex::Real>& t_surf_arr,
224  const amrex::Array4<amrex::Real>& q_surf_arr,
225  const amrex::Array4<amrex::Real>& olen_arr,
226  const amrex::Array4<amrex::Real>& pblh_arr,
227  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
228  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
229  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
230  {
231  int iter = 0;
232  amrex::Real ustar = zero;
233  amrex::Real wstar = zero;
234  amrex::Real tflux = zero;
235  amrex::Real qflux = zero;
236  amrex::Real zeta = zero;
237  amrex::Real psi_m = zero;
238  amrex::Real psi_h = zero;
239  amrex::Real Olen = zero;
240  amrex::Real zref = zref_arr(i,j,k);
241  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
242  if (u_star_arr(i,j,k) == bogus_large_value) {
243  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(zref / z0_arr(i,j,k));
244  } else {
245  Olen = olen_arr(i,j,k);
246  zeta = zref / Olen;
247  psi_m = sfuns.calc_psi_m(zeta);
248  psi_h = sfuns.calc_psi_h(zeta);
249  }
250  do {
251  ustar = u_star_arr(i,j,k);
252  qflux = (spec_qflux) ? mdata.surf_moist_flux :
253  -(qvm_arr(i,j,0) - q_surf_arr(i,j,k)) * ustar * mdata.kappa /
254  (std::log(zref / z0_arr(i,j,k)) - psi_h); // <w'Qv'>
255  tflux = mdata.surf_temp_flux*(1 + amrex::Real(0.61)*qvm_arr(i,j,0)) + qflux*amrex::Real(0.61)*tm_arr(i,j,0);
256  if (w_star_arr) {
257  // update w* and Umagmean
258  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,0));
259  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
260  umm = std::sqrt(umm_arr(i,j,0)*umm_arr(i,j,0) + wstar*wstar);
261  umm = std::max(umm, WSMIN);
262  }
263  Olen = -ustar * ustar * ustar * tvm_arr(i,j,0) / (mdata.kappa * mdata.gravity * tflux);
264  zeta = zref / Olen;
265  psi_m = sfuns.calc_psi_m(zeta);
266  psi_h = sfuns.calc_psi_h(zeta);
267  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(zref / z0_arr(i,j,k)) - psi_m);
268  ++iter;
269  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
270  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
271  "Maximum number of MOST iterations reached.");
272 
273  // Populate the stored MOST arrays
274  olen_arr(i,j,k) = Olen;
275  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(zref / z0_arr(i,j,k)) - psi_h) /
276  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,0);
277  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
278  if (spec_qflux) {
279  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (std::log(zref / z0_arr(i,j,k)) - psi_h) /
280  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,0);
281  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
282  } else {
283  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,0) - q_surf_arr(i,j,k)) /
284  (std::log(zref / z0_arr(i,j,k)) - psi_h);
285  }
286  }
287 
288 private:
292  const amrex::Real tol = amrex::Real(1.0e-5);
293  const amrex::Real WSMIN = amrex::Real(0.1); // minimum wind speed
294 };
295 
296 /**
297  * Moeng flux formulation
298  */
300 {
302 
303  AMREX_GPU_DEVICE
304  AMREX_FORCE_INLINE
306  compute_q_flux (const int& i,
307  const int& j,
308  const int& k,
309  const amrex::Array4<const amrex::Real>& cons_arr,
310  const amrex::Array4<const amrex::Real>& velx_arr,
311  const amrex::Array4<const amrex::Real>& vely_arr,
312  const amrex::Array4<const amrex::Real>& umm_arr,
313  const amrex::Array4<const amrex::Real>& qvm_arr,
314  const amrex::Array4<const amrex::Real>& u_star_arr,
315  const amrex::Array4<const amrex::Real>& q_star_arr,
316  const amrex::Array4<const amrex::Real>& q_surf_arr,
317  const amrex::Array4<const amrex::Real>& u_vfrac_arr,
318  const amrex::Array4<const amrex::Real>& v_vfrac_arr) const
319  {
320  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
321  amrex::Real qv = cons_arr(i,j,k,RhoQ1_comp) / rho;
322 
323  // Volume-weighted average of x-face velocities to cell center
324  amrex::Real u_vfrac_sum = u_vfrac_arr(i,j,k) + u_vfrac_arr(i+1,j,k);
325  amrex::Real velx = (u_vfrac_sum > eps) ?
326  (velx_arr(i,j,k) * u_vfrac_arr(i,j,k) + velx_arr(i+1,j,k) * u_vfrac_arr(i+1,j,k))
327  / u_vfrac_sum : zero;
328 
329  // Volume-weighted average of y-face velocities to cell center
330  amrex::Real v_vfrac_sum = v_vfrac_arr(i,j,k) + v_vfrac_arr(i,j+1,k);
331  amrex::Real vely = (v_vfrac_sum > eps) ?
332  (vely_arr(i,j,k) * v_vfrac_arr(i,j,k) + vely_arr(i,j+1,k) * v_vfrac_arr(i,j+1,k))
333  / v_vfrac_sum : zero;
334 
335  amrex::Real qv_mean = qvm_arr(i,j,0);
336  amrex::Real ustar = u_star_arr(i,j,k);
337  amrex::Real qstar = q_star_arr(i,j,k);
338  amrex::Real qv_surf = q_surf_arr(i,j,k);
339  amrex::Real wsp_mean = umm_arr(i,j,0);
340  wsp_mean = std::max(wsp_mean, WSMIN);
341 
342  amrex::Real wsp = std::sqrt(velx*velx+vely*vely);
343  amrex::Real num1 = wsp * (qv_mean-qv_surf);
344  amrex::Real num2 = wsp_mean * (qv-qv_mean);
345 
346  // NOTE: this is rho*<Qv'w'> = -K dQvdz
347  amrex::Real moflux = (std::abs(qstar) > eps) ?
348  -rho*qstar*ustar*(num1+num2)/((qv_mean-qv_surf)*wsp_mean) : zero;
349 
350  return moflux;
351  }
352 
353  AMREX_GPU_DEVICE
354  AMREX_FORCE_INLINE
356  compute_t_flux (const int& i,
357  const int& j,
358  const int& k,
359  const amrex::Array4<const amrex::Real>& cons_arr,
360  const amrex::Array4<const amrex::Real>& velx_arr,
361  const amrex::Array4<const amrex::Real>& vely_arr,
362  const amrex::Array4<const amrex::Real>& velz_arr,
363  const amrex::Array4<const amrex::Real>& umm_arr,
364  const amrex::Array4<const amrex::Real>& tm_arr,
365  const amrex::Array4<const amrex::Real>& u_star_arr,
366  const amrex::Array4<const amrex::Real>& t_star_arr,
367  const amrex::Array4<const amrex::Real>& t_surf_arr,
368  const amrex::Array4<const amrex::Real>& u_vfrac_arr,
369  const amrex::Array4<const amrex::Real>& v_vfrac_arr,
370  const amrex::Array4<const amrex::Real>& w_vfrac_arr,
371  const amrex::Array4<const amrex::Real>& bnorm_arr) const
372  {
373  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
374  amrex::Real theta = cons_arr(i,j,k,RhoTheta_comp) / rho;
375 
376  // Volume-weighted average of x-face velocities to cell center
377  amrex::Real u_vfrac_sum = u_vfrac_arr(i,j,k) + u_vfrac_arr(i+1,j,k);
378  amrex::Real velx = (u_vfrac_sum > eps) ?
379  (velx_arr(i,j,k) * u_vfrac_arr(i,j,k) + velx_arr(i+1,j,k) * u_vfrac_arr(i+1,j,k))
380  / u_vfrac_sum : zero;
381 
382  // Volume-weighted average of y-face velocities to cell center
383  amrex::Real v_vfrac_sum = v_vfrac_arr(i,j,k) + v_vfrac_arr(i,j+1,k);
384  amrex::Real vely = (v_vfrac_sum > eps) ?
385  (vely_arr(i,j,k) * v_vfrac_arr(i,j,k) + vely_arr(i,j+1,k) * v_vfrac_arr(i,j+1,k))
386  / v_vfrac_sum : zero;
387 
388  // Volume-weighted average of z-face velocities to cell center
389  amrex::Real w_vfrac_sum = w_vfrac_arr(i,j,k) + w_vfrac_arr(i,j,k+1);
390  amrex::Real velz = (w_vfrac_sum > eps) ?
391  (velz_arr(i,j,k) * w_vfrac_arr(i,j,k) + velz_arr(i,j,k+1) * w_vfrac_arr(i,j,k+1))
392  / w_vfrac_sum : zero;
393 
394  // Get boundary normal components
395  amrex::Real nx = bnorm_arr(i,j,k,0);
396  amrex::Real ny = bnorm_arr(i,j,k,1);
397  amrex::Real nz = bnorm_arr(i,j,k,2);
398 
399  // Project velocity onto tangent plane (remove normal component)
400  amrex::Real v_dot_n = velx*nx + vely*ny + velz*nz;
401  amrex::Real velx_tangent = velx - v_dot_n * nx;
402  amrex::Real vely_tangent = vely - v_dot_n * ny;
403 
404  amrex::Real theta_mean = tm_arr(i,j,0);
405  amrex::Real ustar = u_star_arr(i,j,k);
406  amrex::Real tstar = t_star_arr(i,j,k);
407  amrex::Real theta_surf = t_surf_arr(i,j,k);
408  amrex::Real wsp_mean = umm_arr(i,j,0);
409  wsp_mean = std::max(wsp_mean, WSMIN);
410 
411  // Use tangential velocity magnitude instead of Cartesian
412  amrex::Real wsp = std::sqrt(velx_tangent*velx_tangent+vely_tangent*vely_tangent);
413  amrex::Real num1 = wsp * (theta_mean-theta_surf);
414  amrex::Real num2 = wsp_mean * (theta-theta_mean);
415 
416  // NOTE: this is rho*<T'w'> = -K dTdz
417  amrex::Real moflux = (std::abs(tstar) > eps) ?
418  -rho*tstar*ustar*(num1+num2)/((theta_mean-theta_surf)*wsp_mean) : zero;
419 
420  return moflux;
421  }
422 
423  AMREX_GPU_DEVICE
424  AMREX_FORCE_INLINE
427  int j,
428  int k,
429  const amrex::Array4<const amrex::Real>& cons_arr,
430  const amrex::Array4<const amrex::Real>& velx_arr,
431  const amrex::Array4<const amrex::Real>& vely_arr,
432  const amrex::Array4<const amrex::Real>& velz_arr,
433  const amrex::Array4<const amrex::Real>& umm_arr,
434  const amrex::Array4<const amrex::Real>& um_arr,
435  const amrex::Array4<const amrex::Real>& u_star_arr,
436  const amrex::Array4<const amrex::Real>& u_vfrac_arr,
437  const amrex::Array4<const amrex::Real>& v_vfrac_arr,
438  const amrex::Array4<const amrex::Real>& w_vfrac_arr,
439  const amrex::Array4<const amrex::Real>& cc_vfrac_arr,
440  const amrex::Array4<const amrex::EBCellFlag>& cc_flag_arr,
441  const amrex::Array4<const amrex::Real>& bnorm_arr,
442  int idir = 0) const
443  {
444  amrex::Real velx, vely, rho, ustar, wsp_mean;
445  amrex::Real velx_tangent, vely_tangent;
446 
447  if (idir == 0) {
448  // x-face: average to x-face
449  velx = velx_arr(i,j,k);
450 
451  // Volume-weighted average of y-face velocities to x-face
452  amrex::Real v_vfrac_sum = v_vfrac_arr(i,j,k) + v_vfrac_arr(i,j+1,k) +
453  v_vfrac_arr(i-1,j,k) + v_vfrac_arr(i-1,j+1,k);
454  vely = (v_vfrac_sum > eps) ?
455  (vely_arr(i,j,k) * v_vfrac_arr(i,j,k) + vely_arr(i,j+1,k) * v_vfrac_arr(i,j+1,k) +
456  vely_arr(i-1,j,k) * v_vfrac_arr(i-1,j,k) + vely_arr(i-1,j+1,k) * v_vfrac_arr(i-1,j+1,k))
457  / v_vfrac_sum : zero;
458 
459  // Volume-weighted average of z-face velocities to x-face
460  amrex::Real w_vfrac_sum = w_vfrac_arr(i,j,k) + w_vfrac_arr(i,j,k+1) +
461  w_vfrac_arr(i-1,j,k) + w_vfrac_arr(i-1,j,k+1);
462  amrex::Real velz = (w_vfrac_sum > eps) ?
463  (velz_arr(i,j,k) * w_vfrac_arr(i,j,k) + velz_arr(i,j,k+1) * w_vfrac_arr(i,j,k+1) +
464  velz_arr(i-1,j,k) * w_vfrac_arr(i-1,j,k) + velz_arr(i-1,j,k+1) * w_vfrac_arr(i-1,j,k+1))
465  / w_vfrac_sum : zero;
466 
467  // Get boundary normal at x-face (already at correct staggered location)
468  amrex::Real nx = bnorm_arr(i,j,k,0);
469  amrex::Real ny = bnorm_arr(i,j,k,1);
470  amrex::Real nz = bnorm_arr(i,j,k,2);
471 
472  // Project velocity onto tangent plane
473  amrex::Real v_dot_n = velx*nx + vely*ny + velz*nz;
474  velx_tangent = velx - v_dot_n * nx;
475  vely_tangent = vely - v_dot_n * ny;
476 
477  // Volume-weighted average of cell-centered density to x-face
478  amrex::Real cc_vfrac_sum = cc_vfrac_arr(i-1,j,k) + cc_vfrac_arr(i,j,k);
479  rho = (cc_vfrac_sum > eps) ?
480  (cons_arr(i-1,j,k,Rho_comp) * cc_vfrac_arr(i-1,j,k) + cons_arr(i,j,k,Rho_comp) * cc_vfrac_arr(i,j,k))
481  / cc_vfrac_sum : zero;
482 
483  // Average cell-centered u_star and wsp_mean to x-face, using only valid (SingleValued) cells
484  bool low_valid = cc_flag_arr(i-1,j,k).isSingleValued();
485  bool high_valid = cc_flag_arr(i,j,k).isSingleValued();
486 
487  if (low_valid && high_valid) {
488  ustar = myhalf * (u_star_arr(i-1,j,k) + u_star_arr(i,j,k));
489  wsp_mean = myhalf * (umm_arr(i-1,j,0) + umm_arr(i,j,0));
490  } else if (low_valid) {
491  ustar = u_star_arr(i-1,j,k);
492  wsp_mean = umm_arr(i-1,j,0);
493  } else if (high_valid) {
494  ustar = u_star_arr(i,j,k);
495  wsp_mean = umm_arr(i,j,0);
496  } else {
497  ustar = zero;
498  wsp_mean = WSMIN;
499  }
500 
501  } else if (idir == 1) {
502  // y-face: average to y-face
503  vely = vely_arr(i,j,k);
504 
505  // Volume-weighted average of x-face velocities to y-face
506  amrex::Real u_vfrac_sum = u_vfrac_arr(i,j,k) + u_vfrac_arr(i+1,j,k) +
507  u_vfrac_arr(i,j-1,k) + u_vfrac_arr(i+1,j-1,k);
508  velx = (u_vfrac_sum > eps) ?
509  (velx_arr(i,j,k) * u_vfrac_arr(i,j,k) + velx_arr(i+1,j,k) * u_vfrac_arr(i+1,j,k) +
510  velx_arr(i,j-1,k) * u_vfrac_arr(i,j-1,k) + velx_arr(i+1,j-1,k) * u_vfrac_arr(i+1,j-1,k))
511  / u_vfrac_sum : zero;
512 
513  // Volume-weighted average of z-face velocities to y-face
514  amrex::Real w_vfrac_sum = w_vfrac_arr(i,j,k) + w_vfrac_arr(i,j,k+1) +
515  w_vfrac_arr(i,j-1,k) + w_vfrac_arr(i,j-1,k+1);
516  amrex::Real velz = (w_vfrac_sum > eps) ?
517  (velz_arr(i,j,k) * w_vfrac_arr(i,j,k) + velz_arr(i,j,k+1) * w_vfrac_arr(i,j,k+1) +
518  velz_arr(i,j-1,k) * w_vfrac_arr(i,j-1,k) + velz_arr(i,j-1,k+1) * w_vfrac_arr(i,j-1,k+1))
519  / w_vfrac_sum : zero;
520 
521  // Get boundary normal at y-face (already at correct staggered location)
522  amrex::Real nx = bnorm_arr(i,j,k,0);
523  amrex::Real ny = bnorm_arr(i,j,k,1);
524  amrex::Real nz = bnorm_arr(i,j,k,2);
525 
526  // Project velocity onto tangent plane
527  amrex::Real v_dot_n = velx*nx + vely*ny + velz*nz;
528  velx_tangent = velx - v_dot_n * nx;
529  vely_tangent = vely - v_dot_n * ny;
530 
531  // Volume-weighted average of cell-centered density to y-face
532  amrex::Real cc_vfrac_sum = cc_vfrac_arr(i,j-1,k) + cc_vfrac_arr(i,j,k);
533  rho = (cc_vfrac_sum > eps) ?
534  (cons_arr(i,j-1,k,Rho_comp) * cc_vfrac_arr(i,j-1,k) + cons_arr(i,j,k,Rho_comp) * cc_vfrac_arr(i,j,k))
535  / cc_vfrac_sum : zero;
536 
537  // Average cell-centered u_star and wsp_mean to y-face, using only valid (SingleValued) cells
538  bool low_valid = cc_flag_arr(i,j-1,k).isSingleValued();
539  bool high_valid = cc_flag_arr(i,j,k).isSingleValued();
540 
541  if (low_valid && high_valid) {
542  ustar = myhalf * (u_star_arr(i,j-1,k) + u_star_arr(i,j,k));
543  wsp_mean = myhalf * (umm_arr(i,j-1,0) + umm_arr(i,j,0));
544  } else if (low_valid) {
545  ustar = u_star_arr(i,j-1,k);
546  wsp_mean = umm_arr(i,j-1,0);
547  } else if (high_valid) {
548  ustar = u_star_arr(i,j,k);
549  wsp_mean = umm_arr(i,j,0);
550  } else {
551  ustar = zero;
552  wsp_mean = WSMIN;
553  }
554 
555  } else {
556  // z-face: average to z-face
557  // Volume-weighted average of x-face velocities to z-face
558  amrex::Real u_vfrac_sum = u_vfrac_arr(i,j,k-1) + u_vfrac_arr(i+1,j,k-1) +
559  u_vfrac_arr(i,j,k) + u_vfrac_arr(i+1,j,k);
560  velx = (u_vfrac_sum > eps) ?
561  (velx_arr(i,j,k-1) * u_vfrac_arr(i,j,k-1) + velx_arr(i+1,j,k-1) * u_vfrac_arr(i+1,j,k-1) +
562  velx_arr(i,j,k) * u_vfrac_arr(i,j,k) + velx_arr(i+1,j,k) * u_vfrac_arr(i+1,j,k))
563  / u_vfrac_sum : zero;
564 
565  // Volume-weighted average of y-face velocities to z-face
566  amrex::Real v_vfrac_sum = v_vfrac_arr(i,j,k-1) + v_vfrac_arr(i,j+1,k-1) +
567  v_vfrac_arr(i,j,k) + v_vfrac_arr(i,j+1,k);
568  vely = (v_vfrac_sum > eps) ?
569  (vely_arr(i,j,k-1) * v_vfrac_arr(i,j,k-1) + vely_arr(i,j+1,k-1) * v_vfrac_arr(i,j+1,k-1) +
570  vely_arr(i,j,k) * v_vfrac_arr(i,j,k) + vely_arr(i,j+1,k) * v_vfrac_arr(i,j+1,k))
571  / v_vfrac_sum : zero;
572 
573  // z-velocity is already at z-face
574  amrex::Real velz = velz_arr(i,j,k);
575 
576  // Get boundary normal at z-face (already at correct staggered location)
577  amrex::Real nx = bnorm_arr(i,j,k,0);
578  amrex::Real ny = bnorm_arr(i,j,k,1);
579  amrex::Real nz = bnorm_arr(i,j,k,2);
580 
581  // Project velocity onto tangent plane
582  amrex::Real v_dot_n = velx*nx + vely*ny + velz*nz;
583  velx_tangent = velx - v_dot_n * nx;
584  vely_tangent = vely - v_dot_n * ny;
585 
586  // Volume-weighted average of cell-centered density to z-face
587  amrex::Real cc_vfrac_sum = cc_vfrac_arr(i,j,k-1) + cc_vfrac_arr(i,j,k);
588  rho = (cc_vfrac_sum > eps) ?
589  (cons_arr(i,j,k-1,Rho_comp) * cc_vfrac_arr(i,j,k-1) + cons_arr(i,j,k,Rho_comp) * cc_vfrac_arr(i,j,k))
590  / cc_vfrac_sum : zero;
591 
592  // Average cell-centered u_star and wsp_mean to z-face, using only valid (SingleValued) cells
593  bool low_valid = cc_flag_arr(i,j,k-1).isSingleValued();
594  bool high_valid = cc_flag_arr(i,j,k).isSingleValued();
595 
596  if (low_valid && high_valid) {
597  ustar = myhalf * (u_star_arr(i,j,k-1) + u_star_arr(i,j,k));
598  wsp_mean = myhalf * (umm_arr(i,j,0) + umm_arr(i,j,0));
599  } else if (low_valid) {
600  ustar = u_star_arr(i,j,k-1);
601  wsp_mean = umm_arr(i,j,0);
602  } else if (high_valid) {
603  ustar = u_star_arr(i,j,k);
604  wsp_mean = umm_arr(i,j,0);
605  } else {
606  ustar = zero;
607  wsp_mean = WSMIN;
608  }
609  }
610 
611  wsp_mean = std::max(wsp_mean, WSMIN);
612  amrex::Real umean = um_arr(i,j,0);
613 
614  // Note: The surface mean shear stress is decomposed into tau_xz by
615  // multiplying the modeled shear stress (rho*ustar^2) with
616  // a factor of umean/wsp_mean for directionality; this factor
617  // modifies the denominator from what is in Moeng amrex::Real(1984.)
618  amrex::Real wsp = std::sqrt(velx_tangent*velx_tangent+vely_tangent*vely_tangent);
619  amrex::Real num1 = wsp * umean;
620  amrex::Real num2 = wsp_mean * (velx_tangent-umean);
621 
622  // NOTE: this is rho*<u'w'> = -K dudz
623  amrex::Real stressx = -rho*ustar*ustar * (num1+num2)/(wsp_mean*wsp_mean);
624 
625  return stressx;
626  }
627 
628  AMREX_GPU_DEVICE
629  AMREX_FORCE_INLINE
632  int j,
633  int k,
634  const amrex::Array4<const amrex::Real>& cons_arr,
635  const amrex::Array4<const amrex::Real>& velx_arr,
636  const amrex::Array4<const amrex::Real>& vely_arr,
637  const amrex::Array4<const amrex::Real>& velz_arr,
638  const amrex::Array4<const amrex::Real>& umm_arr,
639  const amrex::Array4<const amrex::Real>& vm_arr,
640  const amrex::Array4<const amrex::Real>& u_star_arr,
641  const amrex::Array4<const amrex::Real>& u_vfrac_arr,
642  const amrex::Array4<const amrex::Real>& v_vfrac_arr,
643  const amrex::Array4<const amrex::Real>& w_vfrac_arr,
644  const amrex::Array4<const amrex::Real>& cc_vfrac_arr,
645  const amrex::Array4<const amrex::EBCellFlag>& cc_flag_arr,
646  const amrex::Array4<const amrex::Real>& bnorm_arr,
647  int idir = 0) const
648  {
649  amrex::Real velx, vely, rho, ustar, wsp_mean;
650  amrex::Real velx_tangent, vely_tangent;
651 
652  if (idir == 0) {
653  // x-face: average from cells (i-1) and (i)
654  // Volume-weighted average of y-face velocities to x-face
655  amrex::Real v_vfrac_sum = v_vfrac_arr(i,j,k) + v_vfrac_arr(i,j+1,k) +
656  v_vfrac_arr(i-1,j,k) + v_vfrac_arr(i-1,j+1,k);
657  vely = (v_vfrac_sum > eps) ?
658  (vely_arr(i,j,k) * v_vfrac_arr(i,j,k) + vely_arr(i,j+1,k) * v_vfrac_arr(i,j+1,k) +
659  vely_arr(i-1,j,k) * v_vfrac_arr(i-1,j,k) + vely_arr(i-1,j+1,k) * v_vfrac_arr(i-1,j+1,k))
660  / v_vfrac_sum : zero;
661 
662  velx = velx_arr(i,j,k);
663 
664  // Volume-weighted average of z-face velocities to x-face
665  amrex::Real w_vfrac_sum = w_vfrac_arr(i,j,k) + w_vfrac_arr(i,j,k+1) +
666  w_vfrac_arr(i-1,j,k) + w_vfrac_arr(i-1,j,k+1);
667  amrex::Real velz = (w_vfrac_sum > eps) ?
668  (velz_arr(i,j,k) * w_vfrac_arr(i,j,k) + velz_arr(i,j,k+1) * w_vfrac_arr(i,j,k+1) +
669  velz_arr(i-1,j,k) * w_vfrac_arr(i-1,j,k) + velz_arr(i-1,j,k+1) * w_vfrac_arr(i-1,j,k+1))
670  / w_vfrac_sum : zero;
671 
672  // Get boundary normal at x-face (already at correct staggered location)
673  amrex::Real nx = bnorm_arr(i,j,k,0);
674  amrex::Real ny = bnorm_arr(i,j,k,1);
675  amrex::Real nz = bnorm_arr(i,j,k,2);
676 
677  // Project velocity onto tangent plane
678  amrex::Real v_dot_n = velx*nx + vely*ny + velz*nz;
679  velx_tangent = velx - v_dot_n * nx;
680  vely_tangent = vely - v_dot_n * ny;
681 
682  // Volume-weighted average of cell-centered density to x-face
683  amrex::Real cc_vfrac_sum = cc_vfrac_arr(i-1,j,k) + cc_vfrac_arr(i,j,k);
684  rho = (cc_vfrac_sum > eps) ?
685  (cons_arr(i-1,j,k,Rho_comp) * cc_vfrac_arr(i-1,j,k) + cons_arr(i,j,k,Rho_comp) * cc_vfrac_arr(i,j,k))
686  / cc_vfrac_sum : zero;
687 
688  // Average cell-centered u_star and wsp_mean to x-face, using only valid (SingleValued) cells
689  bool low_valid = cc_flag_arr(i-1,j,k).isSingleValued();
690  bool high_valid = cc_flag_arr(i,j,k).isSingleValued();
691 
692  if (low_valid && high_valid) {
693  ustar = myhalf * (u_star_arr(i-1,j,k) + u_star_arr(i,j,k));
694  wsp_mean = myhalf * (umm_arr(i-1,j,0) + umm_arr(i,j,0));
695  } else if (low_valid) {
696  ustar = u_star_arr(i-1,j,k);
697  wsp_mean = umm_arr(i-1,j,0);
698  } else if (high_valid) {
699  ustar = u_star_arr(i,j,k);
700  wsp_mean = umm_arr(i,j,0);
701  } else {
702  ustar = zero;
703  wsp_mean = WSMIN;
704  }
705 
706  } else if (idir == 1) {
707  // y-face: average from cells (j-1) and (j)
708  // Volume-weighted average of x-face velocities to y-face
709  amrex::Real u_vfrac_sum = u_vfrac_arr(i,j,k) + u_vfrac_arr(i+1,j,k) +
710  u_vfrac_arr(i,j-1,k) + u_vfrac_arr(i+1,j-1,k);
711  velx = (u_vfrac_sum > eps) ?
712  (velx_arr(i,j,k) * u_vfrac_arr(i,j,k) + velx_arr(i+1,j,k) * u_vfrac_arr(i+1,j,k) +
713  velx_arr(i,j-1,k) * u_vfrac_arr(i,j-1,k) + velx_arr(i+1,j-1,k) * u_vfrac_arr(i+1,j-1,k))
714  / u_vfrac_sum : zero;
715 
716  vely = vely_arr(i,j,k);
717 
718  // Volume-weighted average of z-face velocities to y-face
719  amrex::Real w_vfrac_sum = w_vfrac_arr(i,j,k) + w_vfrac_arr(i,j,k+1) +
720  w_vfrac_arr(i,j-1,k) + w_vfrac_arr(i,j-1,k+1);
721  amrex::Real velz = (w_vfrac_sum > eps) ?
722  (velz_arr(i,j,k) * w_vfrac_arr(i,j,k) + velz_arr(i,j,k+1) * w_vfrac_arr(i,j,k+1) +
723  velz_arr(i,j-1,k) * w_vfrac_arr(i,j-1,k) + velz_arr(i,j-1,k+1) * w_vfrac_arr(i,j-1,k+1))
724  / w_vfrac_sum : zero;
725 
726  // Get boundary normal at y-face (already at correct staggered location)
727  amrex::Real nx = bnorm_arr(i,j,k,0);
728  amrex::Real ny = bnorm_arr(i,j,k,1);
729  amrex::Real nz = bnorm_arr(i,j,k,2);
730 
731  // Project velocity onto tangent plane
732  amrex::Real v_dot_n = velx*nx + vely*ny + velz*nz;
733  velx_tangent = velx - v_dot_n * nx;
734  vely_tangent = vely - v_dot_n * ny;
735 
736  // Volume-weighted average of cell-centered density to y-face
737  amrex::Real cc_vfrac_sum = cc_vfrac_arr(i,j-1,k) + cc_vfrac_arr(i,j,k);
738  rho = (cc_vfrac_sum > eps) ?
739  (cons_arr(i,j-1,k,Rho_comp) * cc_vfrac_arr(i,j-1,k) + cons_arr(i,j,k,Rho_comp) * cc_vfrac_arr(i,j,k))
740  / cc_vfrac_sum : zero;
741 
742  // Average cell-centered u_star and wsp_mean to y-face, using only valid (SingleValued) cells
743  bool low_valid = cc_flag_arr(i,j-1,k).isSingleValued();
744  bool high_valid = cc_flag_arr(i,j,k).isSingleValued();
745 
746  if (low_valid && high_valid) {
747  ustar = myhalf * (u_star_arr(i,j-1,k) + u_star_arr(i,j,k));
748  wsp_mean = myhalf * (umm_arr(i,j-1,0) + umm_arr(i,j,0));
749  } else if (low_valid) {
750  ustar = u_star_arr(i,j-1,k);
751  wsp_mean = umm_arr(i,j-1,0);
752  } else if (high_valid) {
753  ustar = u_star_arr(i,j,k);
754  wsp_mean = umm_arr(i,j,0);
755  } else {
756  ustar = zero;
757  wsp_mean = WSMIN;
758  }
759 
760  } else {
761  // z-face: average from cells (k-1) and (k)
762  // Volume-weighted average of x-face velocities to z-face
763  amrex::Real u_vfrac_sum = u_vfrac_arr(i,j,k-1) + u_vfrac_arr(i+1,j,k-1) +
764  u_vfrac_arr(i,j,k) + u_vfrac_arr(i+1,j,k);
765  velx = (u_vfrac_sum > eps) ?
766  (velx_arr(i,j,k-1) * u_vfrac_arr(i,j,k-1) + velx_arr(i+1,j,k-1) * u_vfrac_arr(i+1,j,k-1) +
767  velx_arr(i,j,k) * u_vfrac_arr(i,j,k) + velx_arr(i+1,j,k) * u_vfrac_arr(i+1,j,k))
768  / u_vfrac_sum : zero;
769 
770  // Volume-weighted average of y-face velocities to z-face
771  amrex::Real v_vfrac_sum = v_vfrac_arr(i,j,k-1) + v_vfrac_arr(i,j+1,k-1) +
772  v_vfrac_arr(i,j,k) + v_vfrac_arr(i,j+1,k);
773  vely = (v_vfrac_sum > eps) ?
774  (vely_arr(i,j,k-1) * v_vfrac_arr(i,j,k-1) + vely_arr(i,j+1,k-1) * v_vfrac_arr(i,j+1,k-1) +
775  vely_arr(i,j,k) * v_vfrac_arr(i,j,k) + vely_arr(i,j+1,k) * v_vfrac_arr(i,j+1,k))
776  / v_vfrac_sum : zero;
777 
778  // z-velocity is already at z-face
779  amrex::Real velz = velz_arr(i,j,k);
780 
781  // Get boundary normal at z-face (already at correct staggered location)
782  amrex::Real nx = bnorm_arr(i,j,k,0);
783  amrex::Real ny = bnorm_arr(i,j,k,1);
784  amrex::Real nz = bnorm_arr(i,j,k,2);
785 
786  // Project velocity onto tangent plane
787  amrex::Real v_dot_n = velx*nx + vely*ny + velz*nz;
788  velx_tangent = velx - v_dot_n * nx;
789  vely_tangent = vely - v_dot_n * ny;
790 
791  // Volume-weighted average of cell-centered density to z-face
792  amrex::Real cc_vfrac_sum = cc_vfrac_arr(i,j,k-1) + cc_vfrac_arr(i,j,k);
793  rho = (cc_vfrac_sum > eps) ?
794  (cons_arr(i,j,k-1,Rho_comp) * cc_vfrac_arr(i,j,k-1) + cons_arr(i,j,k,Rho_comp) * cc_vfrac_arr(i,j,k))
795  / cc_vfrac_sum : zero;
796 
797  // Average cell-centered u_star and wsp_mean to z-face, using only valid (SingleValued) cells
798  bool low_valid = cc_flag_arr(i,j,k-1).isSingleValued();
799  bool high_valid = cc_flag_arr(i,j,k).isSingleValued();
800 
801  if (low_valid && high_valid) {
802  ustar = myhalf * (u_star_arr(i,j,k-1) + u_star_arr(i,j,k));
803  wsp_mean = myhalf * (umm_arr(i,j,0) + umm_arr(i,j,0));
804  } else if (low_valid) {
805  ustar = u_star_arr(i,j,k-1);
806  wsp_mean = umm_arr(i,j,0);
807  } else if (high_valid) {
808  ustar = u_star_arr(i,j,k);
809  wsp_mean = umm_arr(i,j,0);
810  } else {
811  ustar = zero;
812  wsp_mean = WSMIN;
813  }
814  }
815 
816  wsp_mean = std::max(wsp_mean, WSMIN);
817  amrex::Real vmean = vm_arr(i,j,0);
818 
819  // Note: The surface mean shear stress is decomposed into tau_yz by
820  // multiplying the modeled shear stress (rho*ustar^2) with
821  // a factor of vmean/wsp_mean for directionality; this factor
822  // modifies the denominator from what is in Moeng amrex::Real(1984.)
823  amrex::Real wsp = std::sqrt(velx_tangent*velx_tangent+vely_tangent*vely_tangent);
824  amrex::Real num1 = wsp * vmean;
825  amrex::Real num2 = wsp_mean * (vely_tangent-vmean);
826 
827  // NOTE: this is rho*<v'w'> = -K dvdz
828  amrex::Real stressy = -rho*ustar*ustar * (num1+num2)/(wsp_mean*wsp_mean);
829 
830  return stressy;
831  }
832 
833 private:
834 #ifdef AMREX_USE_FLOAT
835  const amrex::Real eps = amrex::Real(1e-6);
836 #else
837  const amrex::Real eps = amrex::Real(1e-12);
838 #endif
839  const amrex::Real WSMIN = amrex::Real(0.1); // minimum wind speed
840 };
841 #endif
constexpr amrex::Real bogus_large_value
Definition: ERF_Constants.H:26
constexpr amrex::Real one
Definition: ERF_Constants.H:9
constexpr amrex::Real zero
Definition: ERF_Constants.H:8
constexpr amrex::Real myhalf
Definition: ERF_Constants.H:13
@ num
Definition: ERF_DataStruct.H:24
#define Rho_comp
Definition: ERF_IndexDefines.H:36
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:37
#define RhoQ1_comp
Definition: ERF_IndexDefines.H:42
Real z0
Definition: ERF_InitCustomPertVels_ScalarAdvDiff.H:8
rho
Definition: ERF_InitCustomPert_Bubble.H:106
amrex::Real Real
Definition: ERF_ShocInterface.H:19
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real calc_wstar(const amrex::Real &ust, const amrex::Real &tst, const amrex::Real &qst, const amrex::Real &pblh, const amrex::Real &th, const amrex::Real &thv, const amrex::Real &qv=amrex::Real(0))
Definition: ERF_Wstar.H:13
@ theta
Definition: ERF_MM5.H:20
@ qv
Definition: ERF_Kessler.H:29
@ den
Definition: ERF_AdvanceWSM6.cpp:109
Definition: ERF_EBMOSTStress.H:12
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void iterate_flux(const int &i, const int &j, const int &k, const int &, const amrex::Array4< const amrex::Real > &zref_arr, const amrex::Array4< const amrex::Real > &z0_arr, const amrex::Array4< const amrex::Real > &umm_arr, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< amrex::Real > &u_star_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &t_star_arr, const amrex::Array4< amrex::Real > &q_star_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_EBMOSTStress.H:23
adiabatic_eb(amrex::Real Tflux, amrex::Real Qvflux)
Definition: ERF_EBMOSTStress.H:13
most_data mdata
Definition: ERF_EBMOSTStress.H:52
similarity_funs sfuns
Definition: ERF_EBMOSTStress.H:53
Definition: ERF_EBMOSTStress.H:300
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_u_flux(int i, int j, int k, const amrex::Array4< const amrex::Real > &cons_arr, const amrex::Array4< const amrex::Real > &velx_arr, const amrex::Array4< const amrex::Real > &vely_arr, const amrex::Array4< const amrex::Real > &velz_arr, const amrex::Array4< const amrex::Real > &umm_arr, const amrex::Array4< const amrex::Real > &um_arr, const amrex::Array4< const amrex::Real > &u_star_arr, const amrex::Array4< const amrex::Real > &u_vfrac_arr, const amrex::Array4< const amrex::Real > &v_vfrac_arr, const amrex::Array4< const amrex::Real > &w_vfrac_arr, const amrex::Array4< const amrex::Real > &cc_vfrac_arr, const amrex::Array4< const amrex::EBCellFlag > &cc_flag_arr, const amrex::Array4< const amrex::Real > &bnorm_arr, int idir=0) const
Definition: ERF_EBMOSTStress.H:426
const amrex::Real WSMIN
Definition: ERF_EBMOSTStress.H:839
moeng_flux_eb()
Definition: ERF_EBMOSTStress.H:301
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_t_flux(const int &i, const int &j, const int &k, const amrex::Array4< const amrex::Real > &cons_arr, const amrex::Array4< const amrex::Real > &velx_arr, const amrex::Array4< const amrex::Real > &vely_arr, const amrex::Array4< const amrex::Real > &velz_arr, const amrex::Array4< const amrex::Real > &umm_arr, const amrex::Array4< const amrex::Real > &tm_arr, const amrex::Array4< const amrex::Real > &u_star_arr, const amrex::Array4< const amrex::Real > &t_star_arr, const amrex::Array4< const amrex::Real > &t_surf_arr, const amrex::Array4< const amrex::Real > &u_vfrac_arr, const amrex::Array4< const amrex::Real > &v_vfrac_arr, const amrex::Array4< const amrex::Real > &w_vfrac_arr, const amrex::Array4< const amrex::Real > &bnorm_arr) const
Definition: ERF_EBMOSTStress.H:356
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_v_flux(int i, int j, int k, const amrex::Array4< const amrex::Real > &cons_arr, const amrex::Array4< const amrex::Real > &velx_arr, const amrex::Array4< const amrex::Real > &vely_arr, const amrex::Array4< const amrex::Real > &velz_arr, const amrex::Array4< const amrex::Real > &umm_arr, const amrex::Array4< const amrex::Real > &vm_arr, const amrex::Array4< const amrex::Real > &u_star_arr, const amrex::Array4< const amrex::Real > &u_vfrac_arr, const amrex::Array4< const amrex::Real > &v_vfrac_arr, const amrex::Array4< const amrex::Real > &w_vfrac_arr, const amrex::Array4< const amrex::Real > &cc_vfrac_arr, const amrex::Array4< const amrex::EBCellFlag > &cc_flag_arr, const amrex::Array4< const amrex::Real > &bnorm_arr, int idir=0) const
Definition: ERF_EBMOSTStress.H:631
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_q_flux(const int &i, const int &j, const int &k, const amrex::Array4< const amrex::Real > &cons_arr, const amrex::Array4< const amrex::Real > &velx_arr, const amrex::Array4< const amrex::Real > &vely_arr, const amrex::Array4< const amrex::Real > &umm_arr, const amrex::Array4< const amrex::Real > &qvm_arr, const amrex::Array4< const amrex::Real > &u_star_arr, const amrex::Array4< const amrex::Real > &q_star_arr, const amrex::Array4< const amrex::Real > &q_surf_arr, const amrex::Array4< const amrex::Real > &u_vfrac_arr, const amrex::Array4< const amrex::Real > &v_vfrac_arr) const
Definition: ERF_EBMOSTStress.H:306
const amrex::Real eps
Definition: ERF_EBMOSTStress.H:837
Definition: ERF_MOSTStress.H:13
amrex::Real surf_moist_flux
Moisture flux.
Definition: ERF_MOSTStress.H:19
amrex::Real kappa
von Karman constant
Definition: ERF_MOSTStress.H:16
amrex::Real gravity
Acceleration due to gravity (m/s^2)
Definition: ERF_MOSTStress.H:17
const amrex::Real Bjr_beta
Definition: ERF_MOSTStress.H:32
amrex::Real surf_temp_flux
Heat flux TODO: decide whether this is <θ'w'> or <θv'w'> under moist conditions.
Definition: ERF_MOSTStress.H:18
Definition: ERF_MOSTStress.H:40
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real calc_psi_m(amrex::Real zeta) const
Definition: ERF_MOSTStress.H:90
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real calc_psi_h2(amrex::Real zeta) const
Definition: ERF_MOSTStress.H:67
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real calc_psi_h(amrex::Real zeta) const
Definition: ERF_MOSTStress.H:104
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real calc_psi_m2(amrex::Real zeta) const
Definition: ERF_MOSTStress.H:47
Definition: ERF_EBMOSTStress.H:196
const amrex::Real WSMIN
Definition: ERF_EBMOSTStress.H:293
similarity_funs sfuns
Definition: ERF_EBMOSTStress.H:291
const amrex::Real tol
Definition: ERF_EBMOSTStress.H:292
most_data mdata
Definition: ERF_EBMOSTStress.H:289
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void iterate_flux(const int &i, const int &j, const int &k, const int &max_iters, const amrex::Array4< const amrex::Real > &zref_arr, const amrex::Array4< const amrex::Real > &z0_arr, const amrex::Array4< const amrex::Real > &umm_arr, const amrex::Array4< const amrex::Real > &tm_arr, const amrex::Array4< const amrex::Real > &tvm_arr, const amrex::Array4< const amrex::Real > &qvm_arr, const amrex::Array4< amrex::Real > &u_star_arr, const amrex::Array4< amrex::Real > &w_star_arr, const amrex::Array4< amrex::Real > &t_star_arr, const amrex::Array4< amrex::Real > &q_star_arr, const amrex::Array4< amrex::Real > &t_surf_arr, const amrex::Array4< amrex::Real > &q_surf_arr, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_EBMOSTStress.H:209
bool spec_qflux
Definition: ERF_EBMOSTStress.H:290
surface_flux_eb(amrex::Real Tflux, amrex::Real Qvflux, bool cons_qflux)
Definition: ERF_EBMOSTStress.H:197
Definition: ERF_EBMOSTStress.H:61
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void iterate_flux(const int &i, const int &j, const int &k, const int &max_iters, const amrex::Array4< const amrex::Real > &zref_arr, const amrex::Array4< const amrex::Real > &z0_arr, const amrex::Array4< const amrex::Real > &umm_arr, const amrex::Array4< const amrex::Real > &tm_arr, const amrex::Array4< const amrex::Real > &tvm_arr, const amrex::Array4< const amrex::Real > &qvm_arr, const amrex::Array4< amrex::Real > &u_star_arr, const amrex::Array4< amrex::Real > &w_star_arr, const amrex::Array4< amrex::Real > &t_star_arr, const amrex::Array4< amrex::Real > &q_star_arr, const amrex::Array4< amrex::Real > &t_surf_arr, const amrex::Array4< amrex::Real > &q_surf_arr, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_EBMOSTStress.H:74
const amrex::Real tol
Definition: ERF_EBMOSTStress.H:187
const amrex::Real alpha
Definition: ERF_EBMOSTStress.H:188
const amrex::Real WSMIN
Definition: ERF_EBMOSTStress.H:189
similarity_funs sfuns
Definition: ERF_EBMOSTStress.H:186
surface_temp_eb(amrex::Real Tflux, amrex::Real Qvflux, bool cons_qflux)
Definition: ERF_EBMOSTStress.H:62
bool spec_qflux
Definition: ERF_EBMOSTStress.H:185
most_data mdata
Definition: ERF_EBMOSTStress.H:184