ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
ERF_MOSTStress.H
Go to the documentation of this file.
1 #ifndef ERF_MOSTStress_H
2 #define ERF_MOSTStress_H
3 
4 #include <ERF_Constants.H>
5 #include <ERF_IndexDefines.H>
6 #include <ERF_MOSTRoughness.H>
7 #include <ERF_Wstar.H>
8 
9 /**
10  * Structure of plain old data relevant to MOST BCs
11  */
12 struct most_data
13 {
14 public:
15  amrex::Real zref{10.0}; ///< Reference height (m)
16  amrex::Real z0_const{0.1}; ///< Roughness height -- default constant value(m)
17  amrex::Real kappa{KAPPA}; ///< von Karman constant
18  amrex::Real gravity{CONST_GRAV}; ///< Acceleration due to gravity (m/s^2)
19  amrex::Real surf_temp_flux{0.0}; ///< Heat flux
20 
21  amrex::Real Cnk_a{0.0185}; ///< Standard Charnock constant https://doi.org/10.1175/JAMC-D-17-0137.1
22  amrex::Real Cnk_b1{1.0/30.0}; ///< Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1
23  amrex::Real Cnk_b2{1260.0}; ///< Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1
24  amrex::Real Cnk_d{30.0}; ///< Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1
25  amrex::Real Cnk_b;
26  bool visc{false}; ///< Use viscous Charnock formulation
27 
28  const amrex::Real Bjr_beta = 1.2; // Empirical parameter from Beljaars 1995 QJRMS
29 };
30 
31 
32 /**
33  * Structure of similarity functions for Moeng formulation
34  */
36 {
37  AMREX_GPU_HOST_DEVICE
38  AMREX_FORCE_INLINE
39  amrex::Real
40  calc_psi_m (amrex::Real zeta) const
41  {
42  if (zeta > 0) {
43  return -beta_m * zeta;
44  } else {
45  amrex::Real x = std::sqrt(std::sqrt(1.0 - gamma_m * zeta));
46  return 2.0 * std::log(0.5 * (1.0 + x)) + log(0.5 * (1.0 + x * x)) -
47  2.0 * std::atan(x) + PIoTwo;
48  }
49  }
50 
51  AMREX_GPU_HOST_DEVICE
52  AMREX_FORCE_INLINE
53  amrex::Real
54  calc_psi_h (amrex::Real zeta) const
55  {
56  if (zeta > 0) {
57  return -beta_h * zeta;
58  } else {
59  amrex::Real x = std::sqrt(1.0 - gamma_h * zeta);
60  return 2.0 * std::log(0.5 * (1.0 + x));
61  }
62  }
63 
64 private:
65  amrex::Real beta_m{5.0}; ///< Constants from Dyer, BLM, 1974
66  amrex::Real beta_h{5.0}; ///< https://doi.org/10.1007/BF00240838
67  amrex::Real gamma_m{16.0};
68  amrex::Real gamma_h{16.0};
69 };
70 
71 
72 /**
73  * Empirical kinematic viscosity [m2/s] formula from Andreas (1989) CRREL Rep.
74  * 89-11, valid between -173 and 277 deg C.
75  */
76 AMREX_GPU_DEVICE
77 AMREX_FORCE_INLINE
78 amrex::Real
79 air_viscosity (amrex::Real T_degK)
80 {
81  amrex::Real TC = T_degK - 273.15;
82  return 1.326e-5*(1. + 6.542e-3*TC + 8.301e-6*TC*TC - 4.84e-9*TC*TC*TC);
83 }
84 
85 
86 /**
87  * Adiabatic with constant roughness
88  */
89 struct adiabatic
90 {
91  adiabatic (amrex::Real zref,
92  amrex::Real flux)
93  {
94  mdata.zref = zref;
95  mdata.surf_temp_flux = flux;
96  }
97 
98  AMREX_GPU_DEVICE
99  AMREX_FORCE_INLINE
100  void
101  iterate_flux (const int& i,
102  const int& j,
103  const int& k,
104  const int& /*max_iters*/,
105  const amrex::Array4<const amrex::Real>& z0_arr,
106  const amrex::Array4<const amrex::Real>& umm_arr,
107  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
108  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
109  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
110  const amrex::Array4<amrex::Real>& u_star_arr,
111  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
112  const amrex::Array4<amrex::Real>& t_star_arr,
113  const amrex::Array4<amrex::Real>& /*q_star_arr*/,
114  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
115  const amrex::Array4<amrex::Real>& olen_arr,
116  const amrex::Array4<amrex::Real>& /*pblh_arr*/,
117  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
118  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
119  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
120  {
121  u_star_arr(i,j,k) = mdata.kappa * umm_arr(i,j,k) / std::log(mdata.zref / z0_arr(i,j,k));
122  t_star_arr(i,j,k) = 0.0;
123  olen_arr(i,j,k) = 1.0e16;
124  }
125 
126 private:
129 };
130 
131 
132 /**
133  * Adiabatic with charnock roughness
134  */
136 {
137  adiabatic_charnock (amrex::Real zref,
138  amrex::Real flux,
139  amrex::Real cnk_a,
140  bool cnk_visc)
141  {
142  mdata.zref = zref;
143  mdata.surf_temp_flux = flux;
144  mdata.Cnk_a = cnk_a;
145  mdata.visc = cnk_visc;
146  }
147 
148  AMREX_GPU_DEVICE
149  AMREX_FORCE_INLINE
150  void
151  iterate_flux (const int& i,
152  const int& j,
153  const int& k,
154  const int& max_iters,
155  const amrex::Array4<amrex::Real>& z0_arr,
156  const amrex::Array4<const amrex::Real>& umm_arr,
157  const amrex::Array4<const amrex::Real>& tm_arr,
158  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
159  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
160  const amrex::Array4<amrex::Real>& u_star_arr,
161  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
162  const amrex::Array4<amrex::Real>& t_star_arr,
163  const amrex::Array4<amrex::Real>& /*q_star_arr*/,
164  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
165  const amrex::Array4<amrex::Real>& olen_arr,
166  const amrex::Array4<amrex::Real>& /*pblh_arr*/,
167  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
168  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
169  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
170  {
171  int iter = 0;
172  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
173  amrex::Real ustar = 0.0;
174  amrex::Real z0 = 0.0;
175  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
176  do {
177  ustar = u_star_arr(i,j,k);
178  if (mdata.Cnk_a > 0) {
179  z0 = (mdata.Cnk_a / mdata.gravity) * ustar * ustar;
180  if (mdata.visc) {
181  z0 += air_viscosity(tm_arr(i,j,k)) / std::max(ustar, 0.05);
182  }
183  } else {
184  z0 = COARE3_roughness(mdata.zref, umm, ustar);
185  }
186  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0);
187  ++iter;
188  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
189 
190  t_star_arr(i,j,k) = 0.0;
191  olen_arr(i,j,k) = 1.0e16;
192  z0_arr(i,j,k) = z0;
193  }
194 
195 private:
198  const amrex::Real tol = 1.0e-5;
199  const amrex::Real WSMIN = 0.1; // minimum wind speed
200 };
201 
202 
203 /**
204  * Adiabatic with modified charnock roughness
205  */
207 {
208  adiabatic_mod_charnock (amrex::Real zref,
209  amrex::Real flux,
210  amrex::Real depth)
211  {
212  mdata.zref = zref;
213  mdata.surf_temp_flux = flux;
214  mdata.Cnk_d = depth;
215  mdata.Cnk_b = mdata.Cnk_b1 * std::log(mdata.Cnk_b2 / mdata.Cnk_d);
216  }
217 
218  AMREX_GPU_DEVICE
219  AMREX_FORCE_INLINE
220  void
221  iterate_flux (const int& i,
222  const int& j,
223  const int& k,
224  const int& max_iters,
225  const amrex::Array4<amrex::Real>& z0_arr,
226  const amrex::Array4<const amrex::Real>& umm_arr,
227  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
228  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
229  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
230  const amrex::Array4<amrex::Real>& u_star_arr,
231  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
232  const amrex::Array4<amrex::Real>& t_star_arr,
233  const amrex::Array4<amrex::Real>& /*q_star_arr*/,
234  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
235  const amrex::Array4<amrex::Real>& olen_arr,
236  const amrex::Array4<amrex::Real>& /*pblh_arr*/,
237  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
238  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
239  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
240  {
241  int iter = 0;
242  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
243  amrex::Real ustar = 0.0;
244  amrex::Real z0 = 0.0;
245  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
246  do {
247  ustar = u_star_arr(i,j,k);
248  z0 = std::exp( (2.7*ustar - 1.8/mdata.Cnk_b) / (ustar + 0.17/mdata.Cnk_b) );
249  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0);
250  ++iter;
251  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
252 
253  t_star_arr(i,j,k) = 0.0;
254  olen_arr(i,j,k) = 1.0e16;
255  z0_arr(i,j,k) = z0;
256  }
257 
258 private:
261  const amrex::Real tol = 1.0e-5;
262  const amrex::Real WSMIN = 0.1; // minimum wind speed
263 };
264 
265 
266 /**
267  * Adiabatic with Donelan roughness
268  */
270 {
271  adiabatic_donelan (amrex::Real zref,
272  amrex::Real flux)
273  {
274  mdata.zref = zref;
275  mdata.surf_temp_flux = flux;
276  }
277 
278  AMREX_GPU_DEVICE
279  AMREX_FORCE_INLINE
280  void
281  iterate_flux (const int& i,
282  const int& j,
283  const int& k,
284  const int& max_iters,
285  const amrex::Array4<amrex::Real>& z0_arr,
286  const amrex::Array4<const amrex::Real>& umm_arr,
287  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
288  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
289  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
290  const amrex::Array4<amrex::Real>& u_star_arr,
291  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
292  const amrex::Array4<amrex::Real>& t_star_arr,
293  const amrex::Array4<amrex::Real>& /*q_star_arr*/,
294  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
295  const amrex::Array4<amrex::Real>& olen_arr,
296  const amrex::Array4<amrex::Real>& /*pblh_arr*/,
297  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
298  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
299  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
300  {
301  int iter = 0;
302  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
303  amrex::Real ustar = 0.0;
304  amrex::Real z0 = 0.0;
305  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
306  if (u_star_arr(i,j,k) == 1.E34) {
307  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
308  }
309  do {
310  ustar = u_star_arr(i,j,k);
311  z0 = Donelan_roughness(ustar);
312  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0);
313  ++iter;
314  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
315 
316  t_star_arr(i,j,k) = 0.0;
317  olen_arr(i,j,k) = 1.0e16;
318  z0_arr(i,j,k) = z0;
319  }
320 
321 private:
324  const amrex::Real tol = 1.0e-5;
325  const amrex::Real WSMIN = 0.1; // minimum wind speed
326 };
327 
328 
329 /**
330  * Adiabatic with wave-coupled roughness
331  */
333 {
334  adiabatic_wave_coupled (amrex::Real zref,
335  amrex::Real flux)
336  {
337  mdata.zref = zref;
338  mdata.surf_temp_flux = flux;
339  }
340 
341  AMREX_GPU_DEVICE
342  AMREX_FORCE_INLINE
343  void
344  iterate_flux (const int& i,
345  const int& j,
346  const int& k,
347  const int& max_iters,
348  const amrex::Array4<amrex::Real>& z0_arr,
349  const amrex::Array4<const amrex::Real>& umm_arr,
350  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
351  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
352  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
353  const amrex::Array4<amrex::Real>& u_star_arr,
354  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
355  const amrex::Array4<amrex::Real>& t_star_arr,
356  const amrex::Array4<amrex::Real>& /*q_star_arr*/,
357  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
358  const amrex::Array4<amrex::Real>& olen_arr,
359  const amrex::Array4<amrex::Real>& /*pblh_arr*/,
360  const amrex::Array4<amrex::Real>& Hwave_arr,
361  const amrex::Array4<amrex::Real>& Lwave_arr,
362  const amrex::Array4<amrex::Real>& eta_arr) const
363  {
364  int iter = 0;
365  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
366  amrex::Real ustar = 0.0;
367  amrex::Real z0 = 0.0;
368  int ie, je;
369  ie = i < lbound(eta_arr).x ? lbound(eta_arr).x : i;
370  je = j < lbound(eta_arr).y ? lbound(eta_arr).y : j;
371  ie = ie > ubound(eta_arr).x ? ubound(eta_arr).x : ie;
372  je = je > ubound(eta_arr).y ? ubound(eta_arr).y : je;
373  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
374  do {
375  ustar = u_star_arr(i,j,k);
376  z0 = std::min( std::max(1200.0 * Hwave_arr(i,j,k) * std::pow( Hwave_arr(i,j,k)/(Lwave_arr(i,j,k)+eps), 4.5 )
377  + 0.11 * eta_arr(ie,je,k,EddyDiff::Mom_v) / ustar, z0_eps), z0_max );
378  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0);
379  ++iter;
380  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
381 
382  t_star_arr(i,j,k) = 0.0;
383  olen_arr(i,j,k) = 1.0e16;
384  z0_arr(i,j,k) = z0;
385  }
386 
387 private:
390  const amrex::Real tol = 1.0e-5;
391  const amrex::Real eps = 1e-15;
392  const amrex::Real z0_eps = 1.0e-6;
393  const amrex::Real z0_max = 0.1;
394  const amrex::Real WSMIN = 0.1; // minimum wind speed
395 };
396 
397 
398 /**
399  * Surface flux with constant roughness
400  */
402 {
403  surface_flux (amrex::Real zref,
404  amrex::Real flux)
405  {
406  mdata.zref = zref;
407  mdata.surf_temp_flux = flux;
408  }
409 
410  AMREX_GPU_DEVICE
411  AMREX_FORCE_INLINE
412  void
413  iterate_flux (const int& i,
414  const int& j,
415  const int& k,
416  const int& max_iters,
417  const amrex::Array4<const amrex::Real>& z0_arr,
418  const amrex::Array4<const amrex::Real>& umm_arr,
419  const amrex::Array4<const amrex::Real>& tm_arr,
420  const amrex::Array4<const amrex::Real>& tvm_arr,
421  const amrex::Array4<const amrex::Real>& qvm_arr,
422  const amrex::Array4<amrex::Real>& u_star_arr,
423  const amrex::Array4<amrex::Real>& w_star_arr,
424  const amrex::Array4<amrex::Real>& t_star_arr,
425  const amrex::Array4<amrex::Real>& q_star_arr,
426  const amrex::Array4<amrex::Real>& t_surf_arr,
427  const amrex::Array4<amrex::Real>& olen_arr,
428  const amrex::Array4<amrex::Real>& pblh_arr,
429  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
430  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
431  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
432  {
433  int iter = 0;
434  amrex::Real ustar = 0.0;
435  amrex::Real wstar = 0.0;
436  amrex::Real tflux = 0.0;
437  amrex::Real zeta = 0.0;
438  amrex::Real psi_m = 0.0;
439  amrex::Real psi_h = 0.0;
440  amrex::Real Olen = 0.0;
441  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
442  if (u_star_arr(i,j,k) == 1.E34) {
443  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
444  } else {
445  Olen = olen_arr(i,j,k);
446  zeta = mdata.zref / Olen;
447  psi_m = sfuns.calc_psi_m(zeta);
448  psi_h = sfuns.calc_psi_h(zeta);
449  }
450  do {
451  ustar = u_star_arr(i,j,k);
452  tflux = mdata.surf_temp_flux*(1 + 0.61*qvm_arr(i,j,k)) - 0.61*tm_arr(i,j,k)*ustar*q_star_arr(i,j,k);
453  if (w_star_arr) {
454  // update w* and Umagmean
455  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
456  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
457  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
458  umm = std::max(umm, WSMIN);
459  }
460  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
461  zeta = mdata.zref / Olen;
462  psi_m = sfuns.calc_psi_m(zeta);
463  psi_h = sfuns.calc_psi_h(zeta);
464  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0_arr(i,j,k)) - psi_m);
465  ++iter;
466  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
467  AMREX_ASSERT_WITH_MESSAGE(iter < max_iters, "Maximum number of MOST iterations reached.");
468 
469  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h) /
470  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,k);
471  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
472  olen_arr(i,j,k) = Olen;
473  }
474 
475 private:
478  const amrex::Real tol = 1.0e-5;
479  const amrex::Real WSMIN = 0.1; // minimum wind speed
480 };
481 
482 
483 /**
484  * Surface flux with charnock roughness
485  */
487 {
488  surface_flux_charnock (amrex::Real zref,
489  amrex::Real flux,
490  amrex::Real cnk_a,
491  bool cnk_visc)
492  {
493  mdata.zref = zref;
494  mdata.surf_temp_flux = flux;
495  mdata.Cnk_a = cnk_a;
496  mdata.visc = cnk_visc;
497  }
498 
499  AMREX_GPU_DEVICE
500  AMREX_FORCE_INLINE
501  void
502  iterate_flux (const int& i,
503  const int& j,
504  const int& k,
505  const int& max_iters,
506  const amrex::Array4<amrex::Real>& z0_arr,
507  const amrex::Array4<const amrex::Real>& umm_arr,
508  const amrex::Array4<const amrex::Real>& tm_arr,
509  const amrex::Array4<const amrex::Real>& tvm_arr,
510  const amrex::Array4<const amrex::Real>& qvm_arr,
511  const amrex::Array4<amrex::Real>& u_star_arr,
512  const amrex::Array4<amrex::Real>& w_star_arr,
513  const amrex::Array4<amrex::Real>& t_star_arr,
514  const amrex::Array4<amrex::Real>& q_star_arr,
515  const amrex::Array4<amrex::Real>& t_surf_arr,
516  const amrex::Array4<amrex::Real>& olen_arr,
517  const amrex::Array4<amrex::Real>& pblh_arr,
518  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
519  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
520  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
521  {
522  int iter = 0;
523  amrex::Real ustar = 0.0;
524  amrex::Real wstar = 0.0;
525  amrex::Real tflux = 0.0;
526  amrex::Real z0 = 0.0;
527  amrex::Real zeta = 0.0;
528  amrex::Real psi_m = 0.0;
529  amrex::Real psi_h = 0.0;
530  amrex::Real Olen = 0.0;
531  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
532  if (u_star_arr(i,j,k) == 1.E34) {
533  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
534  } else {
535  Olen = olen_arr(i,j,k);
536  zeta = mdata.zref / Olen;
537  psi_m = sfuns.calc_psi_m(zeta);
538  psi_h = sfuns.calc_psi_h(zeta);
539  }
540  do {
541  ustar = u_star_arr(i,j,k);
542  if (mdata.Cnk_a > 0) {
543  z0 = (mdata.Cnk_a / mdata.gravity) * ustar * ustar;
544  if (mdata.visc) {
545  z0 += air_viscosity(tm_arr(i,j,k)) / std::max(ustar, 0.05);
546  }
547  } else {
548  z0 = COARE3_roughness(mdata.zref, umm, ustar);
549  }
550  tflux = mdata.surf_temp_flux*(1 + 0.61*qvm_arr(i,j,k)) - 0.61*tm_arr(i,j,k)*ustar*q_star_arr(i,j,k);
551  if (w_star_arr) {
552  // update w* and Umagmean
553  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
554  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
555  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
556  umm = std::max(umm, WSMIN);
557  }
558  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
559  zeta = mdata.zref / Olen;
560  psi_m = sfuns.calc_psi_m(zeta);
561  psi_h = sfuns.calc_psi_h(zeta);
562  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
563  ++iter;
564  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
565  AMREX_ASSERT_WITH_MESSAGE(iter < max_iters, "Maximum number of MOST iterations reached.");
566 
567  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(mdata.zref / z0) - psi_h) /
568  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,k);
569  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
570  olen_arr(i,j,k) = Olen;
571  z0_arr(i,j,k) = z0;
572  }
573 
574 private:
577  const amrex::Real tol = 1.0e-5;
578  const amrex::Real WSMIN = 0.1; // minimum wind speed
579 };
580 
581 
582 /**
583  * Surface flux with modified charnock roughness
584  */
586 {
587  surface_flux_mod_charnock (amrex::Real zref,
588  amrex::Real flux,
589  amrex::Real depth)
590  {
591  mdata.zref = zref;
592  mdata.surf_temp_flux = flux;
593  mdata.Cnk_d = depth;
594  mdata.Cnk_b = mdata.Cnk_b1 * std::log(mdata.Cnk_b2 / mdata.Cnk_d);
595  }
596 
597  AMREX_GPU_DEVICE
598  AMREX_FORCE_INLINE
599  void
600  iterate_flux (const int& i,
601  const int& j,
602  const int& k,
603  const int& max_iters,
604  const amrex::Array4<amrex::Real>& z0_arr,
605  const amrex::Array4<const amrex::Real>& umm_arr,
606  const amrex::Array4<const amrex::Real>& tm_arr,
607  const amrex::Array4<const amrex::Real>& tvm_arr,
608  const amrex::Array4<const amrex::Real>& qvm_arr,
609  const amrex::Array4<amrex::Real>& u_star_arr,
610  const amrex::Array4<amrex::Real>& w_star_arr,
611  const amrex::Array4<amrex::Real>& t_star_arr,
612  const amrex::Array4<amrex::Real>& q_star_arr,
613  const amrex::Array4<amrex::Real>& t_surf_arr,
614  const amrex::Array4<amrex::Real>& olen_arr,
615  const amrex::Array4<amrex::Real>& pblh_arr,
616  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
617  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
618  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
619  {
620  int iter = 0;
621  amrex::Real ustar = 0.0;
622  amrex::Real wstar = 0.0;
623  amrex::Real tflux = 0.0;
624  amrex::Real z0 = 0.0;
625  amrex::Real zeta = 0.0;
626  amrex::Real psi_m = 0.0;
627  amrex::Real psi_h = 0.0;
628  amrex::Real Olen = 0.0;
629  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
630  if (u_star_arr(i,j,k) == 1.E34) {
631  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
632  } else {
633  Olen = olen_arr(i,j,k);
634  zeta = mdata.zref / Olen;
635  psi_m = sfuns.calc_psi_m(zeta);
636  psi_h = sfuns.calc_psi_h(zeta);
637  }
638  do {
639  ustar = u_star_arr(i,j,k);
640  z0 = std::exp( (2.7*ustar - 1.8/mdata.Cnk_b) / (ustar + 0.17/mdata.Cnk_b) );
641  tflux = mdata.surf_temp_flux*(1 + 0.61*qvm_arr(i,j,k)) - 0.61*tm_arr(i,j,k)*ustar*q_star_arr(i,j,k);
642  if (w_star_arr) {
643  // update w* and Umagmean
644  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
645  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
646  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
647  umm = std::max(umm, WSMIN);
648  }
649  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
650  zeta = mdata.zref / Olen;
651  psi_m = sfuns.calc_psi_m(zeta);
652  psi_h = sfuns.calc_psi_h(zeta);
653  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
654  ++iter;
655  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
656  AMREX_ASSERT_WITH_MESSAGE(iter < max_iters, "Maximum number of MOST iterations reached.");
657 
658  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(mdata.zref / z0) - psi_h) /
659  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,k);
660  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
661  olen_arr(i,j,k) = Olen;
662  z0_arr(i,j,k) = z0;
663  }
664 
665 private:
668  const amrex::Real tol = 1.0e-5;
669  const amrex::Real WSMIN = 0.1; // minimum wind speed
670 };
671 
672 
673 /**
674  * Surface flux with donelan roughness
675  */
677 {
678  surface_flux_donelan (amrex::Real zref,
679  amrex::Real flux)
680  {
681  mdata.zref = zref;
682  mdata.surf_temp_flux = flux;
683  }
684 
685  AMREX_GPU_DEVICE
686  AMREX_FORCE_INLINE
687  void
688  iterate_flux (const int& i,
689  const int& j,
690  const int& k,
691  const int& max_iters,
692  const amrex::Array4<amrex::Real>& z0_arr,
693  const amrex::Array4<const amrex::Real>& umm_arr,
694  const amrex::Array4<const amrex::Real>& tm_arr,
695  const amrex::Array4<const amrex::Real>& tvm_arr,
696  const amrex::Array4<const amrex::Real>& qvm_arr,
697  const amrex::Array4<amrex::Real>& u_star_arr,
698  const amrex::Array4<amrex::Real>& w_star_arr,
699  const amrex::Array4<amrex::Real>& t_star_arr,
700  const amrex::Array4<amrex::Real>& q_star_arr,
701  const amrex::Array4<amrex::Real>& t_surf_arr,
702  const amrex::Array4<amrex::Real>& olen_arr,
703  const amrex::Array4<amrex::Real>& pblh_arr,
704  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
705  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
706  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
707  {
708  int iter = 0;
709  amrex::Real ustar = 0.0;
710  amrex::Real wstar = 0.0;
711  amrex::Real tflux = 0.0;
712  amrex::Real z0 = 0.0;
713  amrex::Real zeta = 0.0;
714  amrex::Real psi_m = 0.0;
715  amrex::Real psi_h = 0.0;
716  amrex::Real Olen = 0.0;
717  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
718  if (u_star_arr(i,j,k) == 1.E34) {
719  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
720  } else {
721  Olen = olen_arr(i,j,k);
722  zeta = mdata.zref / Olen;
723  psi_m = sfuns.calc_psi_m(zeta);
724  psi_h = sfuns.calc_psi_h(zeta);
725  }
726  do {
727  ustar = u_star_arr(i,j,k);
728  z0 = Donelan_roughness(ustar);
729  tflux = mdata.surf_temp_flux*(1 + 0.61*qvm_arr(i,j,k)) - 0.61*tm_arr(i,j,k)*ustar*q_star_arr(i,j,k);
730  if (w_star_arr) {
731  // update w* and Umagmean
732  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
733  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
734  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
735  umm = std::max(umm, WSMIN);
736  }
737  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
738  zeta = mdata.zref / Olen;
739  psi_m = sfuns.calc_psi_m(zeta);
740  psi_h = sfuns.calc_psi_h(zeta);
741  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
742  ++iter;
743  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
744  AMREX_ASSERT_WITH_MESSAGE(iter < max_iters, "Maximum number of MOST iterations reached.");
745 
746  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(mdata.zref / z0) - psi_h) /
747  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,k);
748  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
749  olen_arr(i,j,k) = Olen;
750  z0_arr(i,j,k) = z0;
751  }
752 
753 private:
756  const amrex::Real tol = 1.0e-5;
757  const amrex::Real WSMIN = 0.1; // minimum wind speed
758 };
759 
760 
761 /**
762  * Surface flux with wave-coupled roughness
763  */
765 {
766  surface_flux_wave_coupled (amrex::Real zref,
767  amrex::Real flux)
768  {
769  mdata.zref = zref;
770  mdata.surf_temp_flux = flux;
771  }
772 
773  AMREX_GPU_DEVICE
774  AMREX_FORCE_INLINE
775  void
776  iterate_flux (const int& i,
777  const int& j,
778  const int& k,
779  const int& max_iters,
780  const amrex::Array4<amrex::Real>& z0_arr,
781  const amrex::Array4<const amrex::Real>& umm_arr,
782  const amrex::Array4<const amrex::Real>& tm_arr,
783  const amrex::Array4<const amrex::Real>& tvm_arr,
784  const amrex::Array4<const amrex::Real>& qvm_arr,
785  const amrex::Array4<amrex::Real>& u_star_arr,
786  const amrex::Array4<amrex::Real>& w_star_arr,
787  const amrex::Array4<amrex::Real>& t_star_arr,
788  const amrex::Array4<amrex::Real>& q_star_arr,
789  const amrex::Array4<amrex::Real>& t_surf_arr,
790  const amrex::Array4<amrex::Real>& olen_arr,
791  const amrex::Array4<amrex::Real>& pblh_arr,
792  const amrex::Array4<amrex::Real>& Hwave_arr,
793  const amrex::Array4<amrex::Real>& Lwave_arr,
794  const amrex::Array4<amrex::Real>& eta_arr) const
795  {
796  int iter = 0;
797  amrex::Real ustar = 0.0;
798  amrex::Real wstar = 0.0;
799  amrex::Real tflux = 0.0;
800  amrex::Real z0 = 0.0;
801  amrex::Real zeta = 0.0;
802  amrex::Real psi_m = 0.0;
803  amrex::Real psi_h = 0.0;
804  amrex::Real Olen = 0.0;
805  int ie, je;
806  ie = i < lbound(eta_arr).x ? lbound(eta_arr).x : i;
807  je = j < lbound(eta_arr).y ? lbound(eta_arr).y : j;
808  ie = ie > ubound(eta_arr).x ? ubound(eta_arr).x : ie;
809  je = je > ubound(eta_arr).y ? ubound(eta_arr).y : je;
810  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
811  if (u_star_arr(i,j,k) == 1.E34) {
812  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
813  } else {
814  Olen = olen_arr(i,j,k);
815  zeta = mdata.zref / Olen;
816  psi_m = sfuns.calc_psi_m(zeta);
817  psi_h = sfuns.calc_psi_h(zeta);
818  }
819  do {
820  ustar = u_star_arr(i,j,k);
821  z0 = std::min( std::max(1200.0 * Hwave_arr(i,j,k) * std::pow( Hwave_arr(i,j,k)/(Lwave_arr(i,j,k)+eps), 4.5 )
822  + 0.11 * eta_arr(ie,je,k,EddyDiff::Mom_v) / ustar, z0_eps), z0_max );
823  tflux = mdata.surf_temp_flux*(1 + 0.61*qvm_arr(i,j,k)) - 0.61*tm_arr(i,j,k)*ustar*q_star_arr(i,j,k);
824  if (w_star_arr) {
825  // update w* and Umagmean
826  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
827  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
828  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
829  umm = std::max(umm, WSMIN);
830  }
831  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
832  zeta = mdata.zref / Olen;
833  psi_m = sfuns.calc_psi_m(zeta);
834  psi_h = sfuns.calc_psi_h(zeta);
835  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
836  ++iter;
837  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
838  AMREX_ASSERT_WITH_MESSAGE(iter < max_iters, "Maximum number of MOST iterations reached.");
839 
840  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(mdata.zref / z0) - psi_h) /
841  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,k);
842  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
843  olen_arr(i,j,k) = Olen;
844  z0_arr(i,j,k) = z0;
845  }
846 
847 private:
850  const amrex::Real tol = 1.0e-5;
851  const amrex::Real eps = 1e-15;
852  const amrex::Real z0_eps = 1.0e-6;
853  const amrex::Real z0_max = 0.1;
854  const amrex::Real WSMIN = 0.1; // minimum wind speed
855 };
856 
857 
858 /**
859  * Surface temperature with constant roughness
860  */
862 {
863  surface_temp (amrex::Real zref,
864  amrex::Real flux)
865  {
866  mdata.zref = zref;
867  mdata.surf_temp_flux = flux;
868  }
869 
870  AMREX_GPU_DEVICE
871  AMREX_FORCE_INLINE
872  void
873  iterate_flux (const int& i,
874  const int& j,
875  const int& k,
876  const int& max_iters,
877  const amrex::Array4<const amrex::Real>& z0_arr,
878  const amrex::Array4<const amrex::Real>& umm_arr,
879  const amrex::Array4<const amrex::Real>& tm_arr,
880  const amrex::Array4<const amrex::Real>& tvm_arr,
881  const amrex::Array4<const amrex::Real>& qvm_arr,
882  const amrex::Array4<amrex::Real>& u_star_arr,
883  const amrex::Array4<amrex::Real>& w_star_arr,
884  const amrex::Array4<amrex::Real>& t_star_arr,
885  const amrex::Array4<amrex::Real>& q_star_arr,
886  const amrex::Array4<amrex::Real>& t_surf_arr,
887  const amrex::Array4<amrex::Real>& olen_arr,
888  const amrex::Array4<amrex::Real>& pblh_arr,
889  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
890  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
891  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
892  {
893  int iter = 0;
894  amrex::Real ustar = 0.0;
895  amrex::Real wstar = 0.0;
896  amrex::Real tflux = 0.0;
897  amrex::Real zeta = 0.0;
898  amrex::Real psi_m = 0.0;
899  amrex::Real psi_h = 0.0;
900  amrex::Real Olen = 0.0;
901  amrex::Real Oleno = 0.0;
902  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
903  if (u_star_arr(i,j,k) == 1.E34) {
904  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
905  } else {
906  Olen = olen_arr(i,j,k);
907  Oleno = Olen;
908  zeta = mdata.zref / Olen;
909  psi_m = sfuns.calc_psi_m(zeta);
910  psi_h = sfuns.calc_psi_h(zeta);
911  }
912  do {
913  ustar = u_star_arr(i,j,k);
914  tflux = -(tm_arr(i,j,k) - t_surf_arr(i,j,k)) * ustar * mdata.kappa /
915  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h); // <w'T'>
916  tflux *= (1 + 0.61*qvm_arr(i,j,k));
917  tflux += 0.61*tm_arr(i,j,k) * -ustar*q_star_arr(i,j,k); // ~= <w'Tv'>
918  if (w_star_arr) {
919  // update w* and Umagmean
920  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
921  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
922  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
923  umm = std::max(umm, WSMIN);
924  }
925  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
926  if ( ((Olen >= 0.0) && (Oleno <= 0.0)) ||
927  ((Olen <= 0.0) && (Oleno >= 0.0)) ) {
928  Olen = 0.5 * (Olen + Oleno);
929  }
930  Oleno = Olen;
931  zeta = mdata.zref / Olen;
932  psi_m = sfuns.calc_psi_m(zeta);
933  psi_h = sfuns.calc_psi_h(zeta);
934  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0_arr(i,j,k)) - psi_m);
935  ++iter;
936  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
937  AMREX_ASSERT_WITH_MESSAGE(iter < max_iters, "Maximum number of MOST iterations reached.");
938 
939  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) /
940  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h);
941  olen_arr(i,j,k) = Olen;
942  }
943 
944 private:
947  const amrex::Real tol = 1.0e-5;
948  const amrex::Real WSMIN = 0.1; // minimum wind speed
949 };
950 
951 
952 /**
953  * Surface temperature with charnock roughness
954  */
956 {
957  surface_temp_charnock (amrex::Real zref,
958  amrex::Real flux,
959  amrex::Real cnk_a,
960  bool cnk_visc)
961  {
962  mdata.zref = zref;
963  mdata.surf_temp_flux = flux;
964  mdata.Cnk_a = cnk_a;
965  mdata.visc = cnk_visc;
966  }
967 
968  AMREX_GPU_DEVICE
969  AMREX_FORCE_INLINE
970  void
971  iterate_flux (const int& i,
972  const int& j,
973  const int& k,
974  const int& max_iters,
975  const amrex::Array4<amrex::Real>& z0_arr,
976  const amrex::Array4<const amrex::Real>& umm_arr,
977  const amrex::Array4<const amrex::Real>& tm_arr,
978  const amrex::Array4<const amrex::Real>& tvm_arr,
979  const amrex::Array4<const amrex::Real>& qvm_arr,
980  const amrex::Array4<amrex::Real>& u_star_arr,
981  const amrex::Array4<amrex::Real>& w_star_arr,
982  const amrex::Array4<amrex::Real>& t_star_arr,
983  const amrex::Array4<amrex::Real>& q_star_arr,
984  const amrex::Array4<amrex::Real>& t_surf_arr,
985  const amrex::Array4<amrex::Real>& olen_arr,
986  const amrex::Array4<amrex::Real>& pblh_arr,
987  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
988  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
989  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
990  {
991  int iter = 0;
992  amrex::Real ustar = 0.0;
993  amrex::Real wstar = 0.0;
994  amrex::Real z0 = 0.0;
995  amrex::Real tflux = 0.0;
996  amrex::Real zeta = 0.0;
997  amrex::Real psi_m = 0.0;
998  amrex::Real psi_h = 0.0;
999  amrex::Real Olen = 0.0;
1000  amrex::Real Oleno = 0.0;
1001  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
1002  if (u_star_arr(i,j,k) == 1.E34) {
1003  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
1004  } else {
1005  Olen = olen_arr(i,j,k);
1006  Oleno = Olen;
1007  zeta = mdata.zref / Olen;
1008  psi_m = sfuns.calc_psi_m(zeta);
1009  psi_h = sfuns.calc_psi_h(zeta);
1010  }
1011  do {
1012  ustar = u_star_arr(i,j,k);
1013  if (mdata.Cnk_a > 0) {
1014  z0 = (mdata.Cnk_a / mdata.gravity) * ustar * ustar;
1015  if (mdata.visc) {
1016  z0 += air_viscosity(tm_arr(i,j,k)) / std::max(ustar, 0.05);
1017  }
1018  } else {
1019  z0 = COARE3_roughness(mdata.zref, umm, ustar);
1020  }
1021  tflux = -(tm_arr(i,j,k) - t_surf_arr(i,j,k)) * ustar * mdata.kappa /
1022  (std::log(mdata.zref / z0) - psi_h); // <w'T'>
1023  tflux *= (1 + 0.61*qvm_arr(i,j,k));
1024  tflux += 0.61*tm_arr(i,j,k) * -ustar*q_star_arr(i,j,k); // ~= <w'Tv'>
1025  if (w_star_arr) {
1026  // update w* and Umagmean
1027  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
1028  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
1029  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
1030  umm = std::max(umm, WSMIN);
1031  }
1032  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
1033  if ( ((Olen >= 0.0) && (Oleno <= 0.0)) ||
1034  ((Olen <= 0.0) && (Oleno >= 0.0)) ) {
1035  Olen = 0.5 * (Olen + Oleno);
1036  }
1037  Oleno = Olen;
1038  zeta = mdata.zref / Olen;
1039  psi_m = sfuns.calc_psi_m(zeta);
1040  psi_h = sfuns.calc_psi_h(zeta);
1041  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
1042  ++iter;
1043  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
1044  AMREX_ASSERT_WITH_MESSAGE(iter < max_iters, "Maximum number of MOST iterations reached.");
1045 
1046  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) /
1047  (std::log(mdata.zref / z0) - psi_h);
1048  olen_arr(i,j,k) = Olen;
1049  z0_arr(i,j,k) = z0;
1050  }
1051 
1052 private:
1055  const amrex::Real tol = 1.0e-5;
1056  const amrex::Real WSMIN = 0.1; // minimum wind speed
1057 };
1058 
1059 
1060 /**
1061  * Surface temperature with modified charnock roughness
1062  */
1064 {
1065  surface_temp_mod_charnock (amrex::Real zref,
1066  amrex::Real flux,
1067  amrex::Real depth)
1068  {
1069  mdata.zref = zref;
1070  mdata.surf_temp_flux = flux;
1071  mdata.Cnk_d = depth;
1072  mdata.Cnk_b = mdata.Cnk_b1 * std::log(mdata.Cnk_b2 / mdata.Cnk_d);
1073  }
1074 
1075  AMREX_GPU_DEVICE
1076  AMREX_FORCE_INLINE
1077  void
1078  iterate_flux (const int& i,
1079  const int& j,
1080  const int& k,
1081  const int& max_iters,
1082  const amrex::Array4<amrex::Real>& z0_arr,
1083  const amrex::Array4<const amrex::Real>& umm_arr,
1084  const amrex::Array4<const amrex::Real>& tm_arr,
1085  const amrex::Array4<const amrex::Real>& tvm_arr,
1086  const amrex::Array4<const amrex::Real>& qvm_arr,
1087  const amrex::Array4<amrex::Real>& u_star_arr,
1088  const amrex::Array4<amrex::Real>& w_star_arr,
1089  const amrex::Array4<amrex::Real>& t_star_arr,
1090  const amrex::Array4<amrex::Real>& q_star_arr,
1091  const amrex::Array4<amrex::Real>& t_surf_arr,
1092  const amrex::Array4<amrex::Real>& olen_arr,
1093  const amrex::Array4<amrex::Real>& pblh_arr,
1094  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
1095  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
1096  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
1097  {
1098  int iter = 0;
1099  amrex::Real ustar = 0.0;
1100  amrex::Real wstar = 0.0;
1101  amrex::Real z0 = 0.0;
1102  amrex::Real tflux = 0.0;
1103  amrex::Real zeta = 0.0;
1104  amrex::Real psi_m = 0.0;
1105  amrex::Real psi_h = 0.0;
1106  amrex::Real Olen = 0.0;
1107  amrex::Real Oleno = 0.0;
1108  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
1109  if (u_star_arr(i,j,k) == 1.E34) {
1110  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
1111  } else {
1112  Olen = olen_arr(i,j,k);
1113  Oleno = Olen;
1114  zeta = mdata.zref / Olen;
1115  psi_m = sfuns.calc_psi_m(zeta);
1116  psi_h = sfuns.calc_psi_h(zeta);
1117  }
1118  do {
1119  ustar = u_star_arr(i,j,k);
1120  z0 = std::exp( (2.7*ustar - 1.8/mdata.Cnk_b) / (ustar + 0.17/mdata.Cnk_b) );
1121  tflux = -(tm_arr(i,j,k) - t_surf_arr(i,j,k)) * ustar * mdata.kappa /
1122  (std::log(mdata.zref / z0) - psi_h); // <w'T'>
1123  tflux *= (1 + 0.61*qvm_arr(i,j,k));
1124  tflux += 0.61*tm_arr(i,j,k) * -ustar*q_star_arr(i,j,k); // ~= <w'Tv'>
1125  if (w_star_arr) {
1126  // update w* and Umagmean
1127  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
1128  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
1129  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
1130  umm = std::max(umm, WSMIN);
1131  }
1132  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
1133  if ( ((Olen >= 0.0) && (Oleno <= 0.0)) ||
1134  ((Olen <= 0.0) && (Oleno >= 0.0)) ) {
1135  Olen = 0.5 * (Olen + Oleno);
1136  }
1137  Oleno = Olen;
1138  zeta = mdata.zref / Olen;
1139  psi_m = sfuns.calc_psi_m(zeta);
1140  psi_h = sfuns.calc_psi_h(zeta);
1141  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
1142  ++iter;
1143  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
1144  AMREX_ASSERT_WITH_MESSAGE(iter < max_iters, "Maximum number of MOST iterations reached.");
1145 
1146  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) /
1147  (std::log(mdata.zref / z0) - psi_h);
1148  olen_arr(i,j,k) = Olen;
1149  z0_arr(i,j,k) = z0;
1150  }
1151 
1152 private:
1155  const amrex::Real tol = 1.0e-5;
1156  const amrex::Real WSMIN = 0.1; // minimum wind speed
1157 };
1158 
1159 
1160 /**
1161  * Surface temperature with donelan roughness
1162  */
1164 {
1165  surface_temp_donelan (amrex::Real zref,
1166  amrex::Real flux)
1167  {
1168  mdata.zref = zref;
1169  mdata.surf_temp_flux = flux;
1170  }
1171 
1172  AMREX_GPU_DEVICE
1173  AMREX_FORCE_INLINE
1174  void
1175  iterate_flux (const int& i,
1176  const int& j,
1177  const int& k,
1178  const int& max_iters,
1179  const amrex::Array4<amrex::Real>& z0_arr,
1180  const amrex::Array4<const amrex::Real>& umm_arr,
1181  const amrex::Array4<const amrex::Real>& tm_arr,
1182  const amrex::Array4<const amrex::Real>& tvm_arr,
1183  const amrex::Array4<const amrex::Real>& qvm_arr,
1184  const amrex::Array4<amrex::Real>& u_star_arr,
1185  const amrex::Array4<amrex::Real>& w_star_arr,
1186  const amrex::Array4<amrex::Real>& t_star_arr,
1187  const amrex::Array4<amrex::Real>& q_star_arr,
1188  const amrex::Array4<amrex::Real>& t_surf_arr,
1189  const amrex::Array4<amrex::Real>& olen_arr,
1190  const amrex::Array4<amrex::Real>& pblh_arr,
1191  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
1192  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
1193  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
1194  {
1195  int iter = 0;
1196  amrex::Real ustar = 0.0;
1197  amrex::Real wstar = 0.0;
1198  amrex::Real z0 = 0.0;
1199  amrex::Real tflux = 0.0;
1200  amrex::Real zeta = 0.0;
1201  amrex::Real psi_m = 0.0;
1202  amrex::Real psi_h = 0.0;
1203  amrex::Real Olen = 0.0;
1204  amrex::Real Oleno = 0.0;
1205  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
1206  if (u_star_arr(i,j,k) == 1.E34) {
1207  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
1208  } else {
1209  Olen = olen_arr(i,j,k);
1210  Oleno = Olen;
1211  zeta = mdata.zref / Olen;
1212  psi_m = sfuns.calc_psi_m(zeta);
1213  psi_h = sfuns.calc_psi_h(zeta);
1214  }
1215  do {
1216  ustar = u_star_arr(i,j,k);
1217  z0 = Donelan_roughness(ustar);
1218  tflux = -(tm_arr(i,j,k) - t_surf_arr(i,j,k)) * ustar * mdata.kappa /
1219  (std::log(mdata.zref / z0) - psi_h); // <w'T'>
1220  tflux *= (1 + 0.61*qvm_arr(i,j,k));
1221  tflux += 0.61*tm_arr(i,j,k) * -ustar*q_star_arr(i,j,k); // ~= <w'Tv'>
1222  if (w_star_arr) {
1223  // update w* and Umagmean
1224  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
1225  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
1226  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
1227  umm = std::max(umm, WSMIN);
1228  }
1229  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
1230  if ( ((Olen >= 0.0) && (Oleno <= 0.0)) ||
1231  ((Olen <= 0.0) && (Oleno >= 0.0)) ) {
1232  Olen = 0.5 * (Olen + Oleno);
1233  }
1234  Oleno = Olen;
1235  zeta = mdata.zref / Olen;
1236  psi_m = sfuns.calc_psi_m(zeta);
1237  psi_h = sfuns.calc_psi_h(zeta);
1238  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
1239  ++iter;
1240  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
1241  AMREX_ASSERT_WITH_MESSAGE(iter < max_iters, "Maximum number of MOST iterations reached.");
1242 
1243  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) /
1244  (std::log(mdata.zref / z0) - psi_h);
1245  olen_arr(i,j,k) = Olen;
1246  z0_arr(i,j,k) = z0;
1247  }
1248 
1249 private:
1252  const amrex::Real tol = 1.0e-5;
1253  const amrex::Real WSMIN = 0.1; // minimum wind speed
1254 };
1255 
1256 
1257 /**
1258  * Surface temperature with wave-coupled roughness
1259  */
1261 {
1262  surface_temp_wave_coupled (amrex::Real zref,
1263  amrex::Real flux)
1264  {
1265  mdata.zref = zref;
1266  mdata.surf_temp_flux = flux;
1267  }
1268 
1269  AMREX_GPU_DEVICE
1270  AMREX_FORCE_INLINE
1271  void
1272  iterate_flux (const int& i,
1273  const int& j,
1274  const int& k,
1275  const int& max_iters,
1276  const amrex::Array4<amrex::Real>& z0_arr,
1277  const amrex::Array4<const amrex::Real>& umm_arr,
1278  const amrex::Array4<const amrex::Real>& tm_arr,
1279  const amrex::Array4<const amrex::Real>& tvm_arr,
1280  const amrex::Array4<const amrex::Real>& qvm_arr,
1281  const amrex::Array4<amrex::Real>& u_star_arr,
1282  const amrex::Array4<amrex::Real>& w_star_arr,
1283  const amrex::Array4<amrex::Real>& t_star_arr,
1284  const amrex::Array4<amrex::Real>& q_star_arr,
1285  const amrex::Array4<amrex::Real>& t_surf_arr,
1286  const amrex::Array4<amrex::Real>& olen_arr,
1287  const amrex::Array4<amrex::Real>& pblh_arr,
1288  const amrex::Array4<amrex::Real>& Hwave_arr,
1289  const amrex::Array4<amrex::Real>& Lwave_arr,
1290  const amrex::Array4<amrex::Real>& eta_arr) const
1291  {
1292  int iter = 0;
1293  amrex::Real ustar = 0.0;
1294  amrex::Real wstar = 0.0;
1295  amrex::Real z0 = 0.0;
1296  amrex::Real tflux = 0.0;
1297  amrex::Real zeta = 0.0;
1298  amrex::Real psi_m = 0.0;
1299  amrex::Real psi_h = 0.0;
1300  amrex::Real Olen = 0.0;
1301  amrex::Real Oleno = 0.0;
1302  int ie, je;
1303  ie = i < lbound(eta_arr).x ? lbound(eta_arr).x : i;
1304  je = j < lbound(eta_arr).y ? lbound(eta_arr).y : j;
1305  ie = ie > ubound(eta_arr).x ? ubound(eta_arr).x : ie;
1306  je = je > ubound(eta_arr).y ? ubound(eta_arr).y : je;
1307  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
1308  if (u_star_arr(i,j,k) == 1.E34) {
1309  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
1310  } else {
1311  Olen = olen_arr(i,j,k);
1312  Oleno = Olen;
1313  zeta = mdata.zref / Olen;
1314  psi_m = sfuns.calc_psi_m(zeta);
1315  psi_h = sfuns.calc_psi_h(zeta);
1316  }
1317  do {
1318  ustar = u_star_arr(i,j,k);
1319  z0 = std::min( std::max(1200.0 * Hwave_arr(i,j,k) * std::pow( Hwave_arr(i,j,k)/(Lwave_arr(i,j,k)+eps), 4.5 )
1320  + 0.11 * eta_arr(ie,je,k,EddyDiff::Mom_v) / ustar, z0_eps), z0_max );
1321  tflux = -(tm_arr(i,j,k) - t_surf_arr(i,j,k)) * ustar * mdata.kappa /
1322  (std::log(mdata.zref / z0) - psi_h); // <w'T'>
1323  tflux *= (1 + 0.61*qvm_arr(i,j,k));
1324  tflux += 0.61*tm_arr(i,j,k) * -ustar*q_star_arr(i,j,k); // ~= <w'Tv'>
1325  if (w_star_arr) {
1326  // update w* and Umagmean
1327  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
1328  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
1329  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
1330  umm = std::max(umm, WSMIN);
1331  }
1332  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
1333  if ( ((Olen >= 0.0) && (Oleno <= 0.0)) ||
1334  ((Olen <= 0.0) && (Oleno >= 0.0)) ) {
1335  Olen = 0.5 * (Olen + Oleno);
1336  }
1337  Oleno = Olen;
1338  zeta = mdata.zref / Olen;
1339  psi_m = sfuns.calc_psi_m(zeta);
1340  psi_h = sfuns.calc_psi_h(zeta);
1341  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
1342  ++iter;
1343  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
1344  AMREX_ASSERT_WITH_MESSAGE(iter < max_iters, "Maximum number of MOST iterations reached.");
1345 
1346  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) /
1347  (std::log(mdata.zref / z0) - psi_h);
1348  olen_arr(i,j,k) = Olen;
1349  z0_arr(i,j,k) = z0;
1350  }
1351 
1352 private:
1355  const amrex::Real tol = 1.0e-5;
1356  const amrex::Real eps = 1e-15;
1357  const amrex::Real z0_eps = 1.0e-6;
1358  const amrex::Real z0_max = 0.1;
1359  const amrex::Real WSMIN = 0.1; // minimum wind speed
1360 };
1361 
1362 
1363 /**
1364  * Moeng flux formulation
1365  */
1367 {
1369 
1370  AMREX_GPU_DEVICE
1371  AMREX_FORCE_INLINE
1372  amrex::Real
1373  compute_q_flux (const int& /*i*/,
1374  const int& /*j*/,
1375  const int& /*k*/,
1376  const amrex::Array4<const amrex::Real>& /*cons_arr*/,
1377  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1378  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1379  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1380  const amrex::Array4<const amrex::Real>& /*qm_arr*/,
1381  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
1382  const amrex::Array4<const amrex::Real>& /*q_star_arr*/,
1383  const amrex::Array4<const amrex::Real>& /*t_surf_arr*/) const
1384  {
1385  // NOTE: this is rho*<q'w'> = -K dqdz
1386  amrex::Real moflux = 0.0;
1387 
1388  return moflux;
1389  }
1390 
1391  AMREX_GPU_DEVICE
1392  AMREX_FORCE_INLINE
1393  amrex::Real
1394  compute_t_flux (const int& i,
1395  const int& j,
1396  const int& k,
1397  const amrex::Array4<const amrex::Real>& cons_arr,
1398  const amrex::Array4<const amrex::Real>& velx_arr,
1399  const amrex::Array4<const amrex::Real>& vely_arr,
1400  const amrex::Array4<const amrex::Real>& umm_arr,
1401  const amrex::Array4<const amrex::Real>& tm_arr,
1402  const amrex::Array4<const amrex::Real>& u_star_arr,
1403  const amrex::Array4<const amrex::Real>& t_star_arr,
1404  const amrex::Array4<const amrex::Real>& t_surf_arr) const
1405  {
1406  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
1407  amrex::Real theta = cons_arr(i,j,k,RhoTheta_comp) / rho;
1408  amrex::Real velx = 0.5 * ( velx_arr(i,j,k) + velx_arr(i+1,j ,k) );
1409  amrex::Real vely = 0.5 * ( vely_arr(i,j,k) + vely_arr(i ,j+1,k) );
1410 
1411  amrex::Real theta_mean = tm_arr(i,j,k);
1412  amrex::Real ustar = u_star_arr(i,j,k);
1413  amrex::Real tstar = t_star_arr(i,j,k);
1414  amrex::Real theta_surf = t_surf_arr(i,j,k);
1415  amrex::Real wsp_mean = umm_arr(i,j,k);
1416  wsp_mean = std::max(wsp_mean, WSMIN);
1417 
1418  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1419  amrex::Real num1 = wsp * (theta_mean-theta_surf);
1420  amrex::Real num2 = wsp_mean * (theta-theta_mean);
1421 
1422  // NOTE: this is rho*<T'w'> = -K dTdz
1423  amrex::Real moflux = (std::abs(tstar) > eps) ?
1424  -rho*tstar*ustar*(num1+num2)/((theta_mean-theta_surf)*wsp_mean) : 0.0;
1425 
1426  return moflux;
1427  }
1428 
1429  AMREX_GPU_DEVICE
1430  AMREX_FORCE_INLINE
1431  amrex::Real
1432  compute_u_flux (const int& i,
1433  const int& j,
1434  const int& k,
1435  const amrex::Array4<const amrex::Real>& cons_arr,
1436  const amrex::Array4<const amrex::Real>& velx_arr,
1437  const amrex::Array4<const amrex::Real>& vely_arr,
1438  const amrex::Array4<const amrex::Real>& umm_arr,
1439  const amrex::Array4<const amrex::Real>& um_arr,
1440  const amrex::Array4<const amrex::Real>& u_star_arr) const
1441  {
1442  amrex::Real velx = velx_arr(i,j,k);
1443  amrex::Real vely = 0.25 * ( vely_arr(i ,j,k) + vely_arr(i ,j+1,k)
1444  + vely_arr(i-1,j,k) + vely_arr(i-1,j+1,k) );
1445  amrex::Real rho = 0.5 * ( cons_arr(i-1,j,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
1446 
1447  amrex::Real umean = um_arr(i,j,k);
1448  amrex::Real ustar = 0.5 * ( u_star_arr(i-1,j,k) + u_star_arr(i,j,k) );
1449  amrex::Real wsp_mean = 0.5 * ( umm_arr(i-1,j,k) + umm_arr(i,j,k) );
1450  wsp_mean = std::max(wsp_mean, WSMIN);
1451 
1452  // Note: The surface mean shear stress is decomposed into tau_xz by
1453  // multiplying the modeled shear stress (rho*ustar^2) with
1454  // a factor of umean/wsp_mean for directionality; this factor
1455  // modifies the denominator from what is in Moeng 1984.
1456  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1457  amrex::Real num1 = wsp * umean;
1458  amrex::Real num2 = wsp_mean * (velx-umean);
1459 
1460  // NOTE: this is rho*<u'w'> = -K dudz
1461  amrex::Real stressx = -rho*ustar*ustar * (num1+num2)/(wsp_mean*wsp_mean);
1462 
1463  return stressx;
1464  }
1465 
1466  AMREX_GPU_DEVICE
1467  AMREX_FORCE_INLINE
1468  amrex::Real
1469  compute_v_flux (const int& i,
1470  const int& j,
1471  const int& k,
1472  const amrex::Array4<const amrex::Real>& cons_arr,
1473  const amrex::Array4<const amrex::Real>& velx_arr,
1474  const amrex::Array4<const amrex::Real>& vely_arr,
1475  const amrex::Array4<const amrex::Real>& umm_arr,
1476  const amrex::Array4<const amrex::Real>& vm_arr,
1477  const amrex::Array4<const amrex::Real>& u_star_arr) const
1478  {
1479  amrex::Real velx = 0.25 * ( velx_arr(i,j ,k) + velx_arr(i+1,j ,k)
1480  + velx_arr(i,j-1,k) + velx_arr(i+1,j-1,k) );
1481  amrex::Real vely = vely_arr(i,j,k);
1482  amrex::Real rho = 0.5 * ( cons_arr(i,j-1,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
1483 
1484  amrex::Real vmean = vm_arr(i,j,k);
1485  amrex::Real ustar = 0.5 * ( u_star_arr(i,j-1,k) + u_star_arr(i,j,k) );
1486  amrex::Real wsp_mean = 0.5 * ( umm_arr(i,j-1,k) + umm_arr(i,j,k) );
1487  wsp_mean = std::max(wsp_mean, WSMIN);
1488 
1489  // Note: The surface mean shear stress is decomposed into tau_yz by
1490  // multiplying the modeled shear stress (rho*ustar^2) with
1491  // a factor of vmean/wsp_mean for directionality; this factor
1492  // modifies the denominator from what is in Moeng 1984.
1493  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1494  amrex::Real num1 = wsp * vmean;
1495  amrex::Real num2 = wsp_mean * (vely-vmean);
1496 
1497  // NOTE: this is rho*<v'w'> = -K dvdz
1498  amrex::Real stressy = -rho*ustar*ustar * (num1+num2)/(wsp_mean*wsp_mean);
1499 
1500  return stressy;
1501  }
1502 
1503 private:
1504  const amrex::Real eps = 1e-15;
1505  const amrex::Real WSMIN = 0.1; // minimum wind speed
1506 };
1507 
1508 
1509 /**
1510  * Donelan flux formulation
1511  */
1513 {
1515 
1516  AMREX_GPU_DEVICE
1517  AMREX_FORCE_INLINE
1518  amrex::Real
1519  compute_q_flux (const int& /*i*/,
1520  const int& /*j*/,
1521  const int& /*k*/,
1522  const amrex::Array4<const amrex::Real>& /*cons_arr*/,
1523  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1524  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1525  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1526  const amrex::Array4<const amrex::Real>& /*qm_arr*/,
1527  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
1528  const amrex::Array4<const amrex::Real>& /*q_star_arr*/,
1529  const amrex::Array4<const amrex::Real>& /*t_surf_arr*/) const
1530  {
1531  // NOTE: this is rho*<q'w'> = -K dqdz
1532  amrex::Real moflux = 0.0;
1533 
1534  return moflux;
1535  }
1536 
1537  AMREX_GPU_DEVICE
1538  AMREX_FORCE_INLINE
1539  amrex::Real
1540  compute_t_flux (const int& i,
1541  const int& j,
1542  const int& k,
1543  const amrex::Array4<const amrex::Real>& cons_arr,
1544  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1545  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1546  const amrex::Array4<const amrex::Real>& umm_arr,
1547  const amrex::Array4<const amrex::Real>& tm_arr,
1548  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
1549  const amrex::Array4<const amrex::Real>& /*t_star_arr*/,
1550  const amrex::Array4<const amrex::Real>& t_surf_arr) const
1551  {
1552  amrex::Real Ch = 0.0012;
1553  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
1554  amrex::Real theta_surf = t_surf_arr(i,j,k);
1555  amrex::Real theta_mean = tm_arr(i,j,k);
1556  amrex::Real wsp_mean = umm_arr(i,j,k);
1557 
1558  // NOTE: this is rho*<T'w'> = -K dTdz
1559  amrex::Real moflux = -rho * Ch * wsp_mean * (theta_mean - theta_surf);
1560 
1561  return moflux;
1562  }
1563 
1564  AMREX_GPU_DEVICE
1565  AMREX_FORCE_INLINE
1566  amrex::Real
1567  compute_u_flux (const int& i,
1568  const int& j,
1569  const int& k,
1570  const amrex::Array4<const amrex::Real>& cons_arr,
1571  const amrex::Array4<const amrex::Real>& velx_arr,
1572  const amrex::Array4<const amrex::Real>& vely_arr,
1573  const amrex::Array4<const amrex::Real>& umm_arr,
1574  const amrex::Array4<const amrex::Real>& /*um_arr*/,
1575  const amrex::Array4<const amrex::Real>& /*u_star_arr*/) const
1576  {
1577  amrex::Real velx = velx_arr(i,j,k);
1578  amrex::Real vely = 0.25 * ( vely_arr(i ,j,k) + vely_arr(i ,j+1,k)
1579  + vely_arr(i-1,j,k) + vely_arr(i-1,j+1,k) );
1580  amrex::Real rho = 0.5 * ( cons_arr(i-1,j,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
1581 
1582  amrex::Real Cd = 0.001;
1583  const amrex::Real c = 7e-5;
1584  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1585  amrex::Real wsp_mean = 0.5 * ( umm_arr(i-1,j,k) + umm_arr(i,j,k) );
1586  if (wsp_mean <= 5.0) {
1587  Cd = 0.001;
1588  } else if (wsp_mean < 25.0 && wsp_mean > 5.0) {
1589  Cd = 0.001 + c * (wsp_mean - 5.0);
1590  } else {
1591  Cd = 0.0024;
1592  }
1593 
1594  // NOTE: this is rho*<u'w'> = -K dudz
1595  amrex::Real stressx = -rho * Cd * velx * wsp;
1596 
1597  return stressx;
1598  }
1599 
1600  AMREX_GPU_DEVICE
1601  AMREX_FORCE_INLINE
1602  amrex::Real
1603  compute_v_flux (const int& i,
1604  const int& j,
1605  const int& k,
1606  const amrex::Array4<const amrex::Real>& cons_arr,
1607  const amrex::Array4<const amrex::Real>& velx_arr,
1608  const amrex::Array4<const amrex::Real>& vely_arr,
1609  const amrex::Array4<const amrex::Real>& umm_arr,
1610  const amrex::Array4<const amrex::Real>& /*vm_arr*/,
1611  const amrex::Array4<const amrex::Real>& /*u_star_arr*/) const
1612  {
1613  amrex::Real velx = 0.25 * ( velx_arr(i,j ,k) + velx_arr(i+1,j ,k)
1614  + velx_arr(i,j-1,k) + velx_arr(i+1,j-1,k) );
1615  amrex::Real vely = vely_arr(i,j,k);
1616  amrex::Real rho = 0.5 * ( cons_arr(i,j-1,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
1617 
1618  amrex::Real Cd = 0.001;
1619  const amrex::Real c = 7e-5;
1620  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1621  amrex::Real wsp_mean = 0.5 * ( umm_arr(i,j-1,k) + umm_arr(i,j,k) );
1622  if (wsp_mean <= 5.0) {
1623  Cd = 0.001;
1624  } else if (wsp_mean < 25.0 && wsp_mean > 5.0) {
1625  Cd = 0.001 + c * (wsp_mean - 5.0);
1626  } else {
1627  Cd = 0.0024;
1628  }
1629 
1630  // NOTE: this is rho*<v'w'> = -K dvdz
1631  amrex::Real stressy = -rho * Cd * vely * wsp;
1632 
1633  return stressy;
1634  }
1635 
1636 private:
1637 
1638 };
1639 
1640 
1641 /**
1642  * Custom flux formulation
1643  */
1645 {
1647 
1648  AMREX_GPU_DEVICE
1649  AMREX_FORCE_INLINE
1650  amrex::Real
1651  compute_q_flux (const int& i,
1652  const int& j,
1653  const int& k,
1654  const amrex::Array4<const amrex::Real>& cons_arr,
1655  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1656  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1657  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1658  const amrex::Array4<const amrex::Real>& /*qm_arr*/,
1659  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
1660  const amrex::Array4<const amrex::Real>& q_star_arr,
1661  const amrex::Array4<const amrex::Real>& /*t_surf_arr*/) const
1662  {
1663  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
1664  amrex::Real qstar = q_star_arr(i,j,k);
1665 
1666  // NOTE: this is rho*<q'w'> = -K dqdz
1667  amrex::Real moflux = (std::abs(qstar) > eps) ? rho * qstar : 0.0;
1668 
1669  return moflux;
1670  }
1671 
1672  AMREX_GPU_DEVICE
1673  AMREX_FORCE_INLINE
1674  amrex::Real
1675  compute_t_flux (const int& i,
1676  const int& j,
1677  const int& k,
1678  const amrex::Array4<const amrex::Real>& cons_arr,
1679  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1680  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1681  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1682  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
1683  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
1684  const amrex::Array4<const amrex::Real>& t_star_arr,
1685  const amrex::Array4<const amrex::Real>& /*t_surf_arr*/) const
1686  {
1687  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
1688  amrex::Real tstar = t_star_arr(i,j,k);
1689 
1690  // NOTE: this is rho*<T'w'> = -K dTdz
1691  amrex::Real moflux = (std::abs(tstar) > eps) ? rho * tstar : 0.0;
1692 
1693  return moflux;
1694  }
1695 
1696  AMREX_GPU_DEVICE
1697  AMREX_FORCE_INLINE
1698  amrex::Real
1699  compute_u_flux (const int& i,
1700  const int& j,
1701  const int& k,
1702  const amrex::Array4<const amrex::Real>& cons_arr,
1703  const amrex::Array4<const amrex::Real>& velx_arr,
1704  const amrex::Array4<const amrex::Real>& vely_arr,
1705  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1706  const amrex::Array4<const amrex::Real>& /*um_arr*/,
1707  const amrex::Array4<const amrex::Real>& u_star_arr) const
1708  {
1709  amrex::Real velx = velx_arr(i,j,k);
1710  amrex::Real vely = 0.25 * ( vely_arr(i ,j,k) + vely_arr(i ,j+1,k)
1711  + vely_arr(i-1,j,k) + vely_arr(i-1,j+1,k) );
1712  amrex::Real rho = 0.5 * ( cons_arr(i-1,j,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
1713 
1714  amrex::Real ustar = 0.5 * ( u_star_arr(i-1,j,k) + u_star_arr(i,j,k) );
1715  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1716 
1717  // NOTE: this is rho*<u'w'> = -K dudz
1718  amrex::Real stressx = -rho * ustar * ustar * velx / wsp;
1719 
1720  return stressx;
1721  }
1722 
1723  AMREX_GPU_DEVICE
1724  AMREX_FORCE_INLINE
1725  amrex::Real
1726  compute_v_flux (const int& i,
1727  const int& j,
1728  const int& k,
1729  const amrex::Array4<const amrex::Real>& cons_arr,
1730  const amrex::Array4<const amrex::Real>& velx_arr,
1731  const amrex::Array4<const amrex::Real>& vely_arr,
1732  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1733  const amrex::Array4<const amrex::Real>& /*vm_arr*/,
1734  const amrex::Array4<const amrex::Real>& u_star_arr) const
1735  {
1736  amrex::Real velx = 0.25 * ( velx_arr(i,j ,k) + velx_arr(i+1,j ,k)
1737  + velx_arr(i,j-1,k) + velx_arr(i+1,j-1,k) );
1738  amrex::Real vely = vely_arr(i,j,k);
1739  amrex::Real rho = 0.5 * ( cons_arr(i,j-1,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
1740 
1741  amrex::Real ustar = 0.5 * ( u_star_arr(i,j-1,k) + u_star_arr(i,j,k) );
1742  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1743 
1744  // NOTE: this is rho*<v'w'> = -K dvdz
1745  amrex::Real stressy = -rho * ustar * ustar * vely / wsp;
1746 
1747  return stressy;
1748  }
1749 
1750 private:
1751  const amrex::Real eps = 1e-15;
1752 
1753 };
1754 
1755 /**
1756  * Rotate flux formulation
1757  */
1759 {
1761 
1762  AMREX_GPU_DEVICE
1763  AMREX_FORCE_INLINE
1764  amrex::Real
1765  compute_q_flux (const int& i,
1766  const int& j,
1767  const int& k,
1768  const amrex::Array4<const amrex::Real>& cons_arr,
1769  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1770  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1771  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1772  const amrex::Array4<const amrex::Real>& /*qm_arr*/,
1773  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
1774  const amrex::Array4<const amrex::Real>& /*q_star_arr*/,
1775  const amrex::Array4<const amrex::Real>& /*t_surf_arr*/) const
1776  {
1777  // NOTE: this is the total stress
1778  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
1779  amrex::Real qstar = 0.0; //q_star_arr(i,j,k);
1780  amrex::Real moflux = (std::abs(qstar) > eps) ? -rho * qstar : 0.0;
1781 
1782  return moflux;
1783  }
1784 
1785  AMREX_GPU_DEVICE
1786  AMREX_FORCE_INLINE
1787  amrex::Real
1788  compute_t_flux (const int& i,
1789  const int& j,
1790  const int& k,
1791  const amrex::Array4<const amrex::Real>& cons_arr,
1792  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1793  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1794  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1795  const amrex::Array4<const amrex::Real>& tm_arr,
1796  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
1797  const amrex::Array4<const amrex::Real>& t_star_arr,
1798  const amrex::Array4<const amrex::Real>& t_surf_arr) const
1799  {
1800  // NOTE: this is the total stress
1801  amrex::Real theta_mean = tm_arr(i,j,k);
1802  amrex::Real theta_surf = t_surf_arr(i,j,k);
1803  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
1804  amrex::Real tstar = t_star_arr(i,j,k);
1805  amrex::Real moflux = (std::abs(tstar) > eps) ? -rho * tstar * (theta_mean-theta_surf) : 0.0;
1806 
1807  return moflux;
1808  }
1809 
1810  AMREX_GPU_DEVICE
1811  AMREX_FORCE_INLINE
1812  amrex::Real
1813  compute_u_flux (const int& i,
1814  const int& j,
1815  const int& k,
1816  const amrex::Array4<const amrex::Real>& cons_arr,
1817  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1818  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1819  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1820  const amrex::Array4<const amrex::Real>& /*um_arr*/,
1821  const amrex::Array4<const amrex::Real>& u_star_arr) const
1822  {
1823  // NOTE: this is the total stress
1824  amrex::Real rho = 0.5 *( cons_arr(i-1,j,k,Rho_comp) + cons_arr(i ,j,k,Rho_comp) );
1825  amrex::Real ustar = 0.5 * ( u_star_arr(i-1,j,k) + u_star_arr(i,j,k) );
1826  amrex::Real stressx = -rho * ustar * ustar;
1827 
1828  return stressx;
1829  }
1830 
1831  AMREX_GPU_DEVICE
1832  AMREX_FORCE_INLINE
1833  amrex::Real
1834  compute_v_flux (const int& /*i*/,
1835  const int& /*j*/,
1836  const int& /*k*/,
1837  const amrex::Array4<const amrex::Real>& /*cons_arr*/,
1838  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1839  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1840  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1841  const amrex::Array4<const amrex::Real>& /*vm_arr*/,
1842  const amrex::Array4<const amrex::Real>& /*u_star_arr*/) const
1843  {
1844  // NOTE: this is the total stress
1845  amrex::Real stressy = 0.0;
1846 
1847  return stressy;
1848  }
1849 
1850 private:
1851  const amrex::Real eps = 1e-15;
1852 };
1853 #endif
constexpr amrex::Real KAPPA
Definition: ERF_Constants.H:20
constexpr amrex::Real CONST_GRAV
Definition: ERF_Constants.H:21
constexpr amrex::Real PIoTwo
Definition: ERF_Constants.H:7
#define Rho_comp
Definition: ERF_IndexDefines.H:36
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:37
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real Donelan_roughness(amrex::Real ustar)
Definition: ERF_MOSTRoughness.H:29
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real COARE3_roughness(amrex::Real zref, amrex::Real umm, amrex::Real ustar)
Definition: ERF_MOSTRoughness.H:9
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real air_viscosity(amrex::Real T_degK)
Definition: ERF_MOSTStress.H:79
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=0.0)
Definition: ERF_Wstar.H:13
@ Mom_v
Definition: ERF_IndexDefines.H:175
@ theta
Definition: ERF_MM5.H:20
@ rho
Definition: ERF_Kessler.H:22
Definition: ERF_MOSTStress.H:136
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:199
adiabatic_charnock(amrex::Real zref, amrex::Real flux, amrex::Real cnk_a, bool cnk_visc)
Definition: ERF_MOSTStress.H:137
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< 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 > &, 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 > &, 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_MOSTStress.H:151
most_data mdata
Definition: ERF_MOSTStress.H:196
const amrex::Real tol
Definition: ERF_MOSTStress.H:198
similarity_funs sfuns
Definition: ERF_MOSTStress.H:197
Definition: ERF_MOSTStress.H:270
similarity_funs sfuns
Definition: ERF_MOSTStress.H:323
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:325
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< 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 > &, 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_MOSTStress.H:281
adiabatic_donelan(amrex::Real zref, amrex::Real flux)
Definition: ERF_MOSTStress.H:271
most_data mdata
Definition: ERF_MOSTStress.H:322
const amrex::Real tol
Definition: ERF_MOSTStress.H:324
Definition: ERF_MOSTStress.H:207
adiabatic_mod_charnock(amrex::Real zref, amrex::Real flux, amrex::Real depth)
Definition: ERF_MOSTStress.H:208
similarity_funs sfuns
Definition: ERF_MOSTStress.H:260
const amrex::Real tol
Definition: ERF_MOSTStress.H:261
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:262
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< 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 > &, 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_MOSTStress.H:221
most_data mdata
Definition: ERF_MOSTStress.H:259
Definition: ERF_MOSTStress.H:333
const amrex::Real eps
Definition: ERF_MOSTStress.H:391
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< 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 > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &Hwave_arr, const amrex::Array4< amrex::Real > &Lwave_arr, const amrex::Array4< amrex::Real > &eta_arr) const
Definition: ERF_MOSTStress.H:344
const amrex::Real z0_eps
Definition: ERF_MOSTStress.H:392
most_data mdata
Definition: ERF_MOSTStress.H:388
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:394
adiabatic_wave_coupled(amrex::Real zref, amrex::Real flux)
Definition: ERF_MOSTStress.H:334
similarity_funs sfuns
Definition: ERF_MOSTStress.H:389
const amrex::Real tol
Definition: ERF_MOSTStress.H:390
const amrex::Real z0_max
Definition: ERF_MOSTStress.H:393
Definition: ERF_MOSTStress.H:90
adiabatic(amrex::Real zref, amrex::Real flux)
Definition: ERF_MOSTStress.H:91
similarity_funs sfuns
Definition: ERF_MOSTStress.H:128
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 > &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 > &, 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_MOSTStress.H:101
most_data mdata
Definition: ERF_MOSTStress.H:127
Definition: ERF_MOSTStress.H:1645
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_u_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 > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &u_star_arr) const
Definition: ERF_MOSTStress.H:1699
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 > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &q_star_arr, const amrex::Array4< const amrex::Real > &) const
Definition: ERF_MOSTStress.H:1651
custom_flux()
Definition: ERF_MOSTStress.H:1646
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 > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &t_star_arr, const amrex::Array4< const amrex::Real > &) const
Definition: ERF_MOSTStress.H:1675
const amrex::Real eps
Definition: ERF_MOSTStress.H:1751
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_v_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 > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &u_star_arr) const
Definition: ERF_MOSTStress.H:1726
Definition: ERF_MOSTStress.H:1513
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_u_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 > &, const amrex::Array4< const amrex::Real > &) const
Definition: ERF_MOSTStress.H:1567
donelan_flux()
Definition: ERF_MOSTStress.H:1514
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_v_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 > &, const amrex::Array4< const amrex::Real > &) const
Definition: ERF_MOSTStress.H:1603
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 > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &umm_arr, const amrex::Array4< const amrex::Real > &tm_arr, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &t_surf_arr) const
Definition: ERF_MOSTStress.H:1540
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_q_flux(const int &, const int &, const int &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &) const
Definition: ERF_MOSTStress.H:1519
Definition: ERF_MOSTStress.H:1367
moeng_flux()
Definition: ERF_MOSTStress.H:1368
const amrex::Real eps
Definition: ERF_MOSTStress.H:1504
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_u_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 > &um_arr, const amrex::Array4< const amrex::Real > &u_star_arr) const
Definition: ERF_MOSTStress.H:1432
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1505
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 > &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
Definition: ERF_MOSTStress.H:1394
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_q_flux(const int &, const int &, const int &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &) const
Definition: ERF_MOSTStress.H:1373
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_v_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 > &vm_arr, const amrex::Array4< const amrex::Real > &u_star_arr) const
Definition: ERF_MOSTStress.H:1469
Definition: ERF_MOSTStress.H:13
amrex::Real zref
Reference height (m)
Definition: ERF_MOSTStress.H:15
amrex::Real Cnk_b2
Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1.
Definition: ERF_MOSTStress.H:23
amrex::Real Cnk_b
Definition: ERF_MOSTStress.H:25
amrex::Real Cnk_d
Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1.
Definition: ERF_MOSTStress.H:24
amrex::Real kappa
von Karman constant
Definition: ERF_MOSTStress.H:17
amrex::Real gravity
Acceleration due to gravity (m/s^2)
Definition: ERF_MOSTStress.H:18
amrex::Real Cnk_a
Standard Charnock constant https://doi.org/10.1175/JAMC-D-17-0137.1.
Definition: ERF_MOSTStress.H:21
const amrex::Real Bjr_beta
Definition: ERF_MOSTStress.H:28
amrex::Real Cnk_b1
Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1.
Definition: ERF_MOSTStress.H:22
amrex::Real z0_const
Roughness height – default constant value(m)
Definition: ERF_MOSTStress.H:16
bool visc
Use viscous Charnock formulation.
Definition: ERF_MOSTStress.H:26
amrex::Real surf_temp_flux
Heat flux.
Definition: ERF_MOSTStress.H:19
Definition: ERF_MOSTStress.H:1759
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 > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &) const
Definition: ERF_MOSTStress.H:1765
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_v_flux(const int &, const int &, const int &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &) const
Definition: ERF_MOSTStress.H:1834
rotate_flux()
Definition: ERF_MOSTStress.H:1760
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_u_flux(const int &i, const int &j, const int &k, const amrex::Array4< const amrex::Real > &cons_arr, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &u_star_arr) const
Definition: ERF_MOSTStress.H:1813
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 > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &tm_arr, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &t_star_arr, const amrex::Array4< const amrex::Real > &t_surf_arr) const
Definition: ERF_MOSTStress.H:1788
const amrex::Real eps
Definition: ERF_MOSTStress.H:1851
Definition: ERF_MOSTStress.H:36
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real calc_psi_m(amrex::Real zeta) const
Definition: ERF_MOSTStress.H:40
amrex::Real beta_m
Constants from Dyer, BLM, 1974.
Definition: ERF_MOSTStress.H:65
amrex::Real beta_h
https://doi.org/10.1007/BF00240838
Definition: ERF_MOSTStress.H:66
amrex::Real gamma_h
Definition: ERF_MOSTStress.H:68
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real calc_psi_h(amrex::Real zeta) const
Definition: ERF_MOSTStress.H:54
amrex::Real gamma_m
Definition: ERF_MOSTStress.H:67
Definition: ERF_MOSTStress.H:487
const amrex::Real tol
Definition: ERF_MOSTStress.H:577
similarity_funs sfuns
Definition: ERF_MOSTStress.H:576
most_data mdata
Definition: ERF_MOSTStress.H:575
surface_flux_charnock(amrex::Real zref, amrex::Real flux, amrex::Real cnk_a, bool cnk_visc)
Definition: ERF_MOSTStress.H:488
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< 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 > &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_MOSTStress.H:502
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:578
Definition: ERF_MOSTStress.H:677
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< 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 > &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_MOSTStress.H:688
const amrex::Real tol
Definition: ERF_MOSTStress.H:756
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:757
most_data mdata
Definition: ERF_MOSTStress.H:754
similarity_funs sfuns
Definition: ERF_MOSTStress.H:755
surface_flux_donelan(amrex::Real zref, amrex::Real flux)
Definition: ERF_MOSTStress.H:678
Definition: ERF_MOSTStress.H:586
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:669
most_data mdata
Definition: ERF_MOSTStress.H:666
surface_flux_mod_charnock(amrex::Real zref, amrex::Real flux, amrex::Real depth)
Definition: ERF_MOSTStress.H:587
const amrex::Real tol
Definition: ERF_MOSTStress.H:668
similarity_funs sfuns
Definition: ERF_MOSTStress.H:667
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< 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 > &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_MOSTStress.H:600
Definition: ERF_MOSTStress.H:765
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:854
surface_flux_wave_coupled(amrex::Real zref, amrex::Real flux)
Definition: ERF_MOSTStress.H:766
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< 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 > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &Hwave_arr, const amrex::Array4< amrex::Real > &Lwave_arr, const amrex::Array4< amrex::Real > &eta_arr) const
Definition: ERF_MOSTStress.H:776
const amrex::Real z0_max
Definition: ERF_MOSTStress.H:853
const amrex::Real z0_eps
Definition: ERF_MOSTStress.H:852
most_data mdata
Definition: ERF_MOSTStress.H:848
similarity_funs sfuns
Definition: ERF_MOSTStress.H:849
const amrex::Real tol
Definition: ERF_MOSTStress.H:850
const amrex::Real eps
Definition: ERF_MOSTStress.H:851
Definition: ERF_MOSTStress.H:402
similarity_funs sfuns
Definition: ERF_MOSTStress.H:477
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 > &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 > &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_MOSTStress.H:413
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:479
most_data mdata
Definition: ERF_MOSTStress.H:476
surface_flux(amrex::Real zref, amrex::Real flux)
Definition: ERF_MOSTStress.H:403
const amrex::Real tol
Definition: ERF_MOSTStress.H:478
Definition: ERF_MOSTStress.H:956
most_data mdata
Definition: ERF_MOSTStress.H:1053
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< 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 > &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_MOSTStress.H:971
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1056
const amrex::Real tol
Definition: ERF_MOSTStress.H:1055
surface_temp_charnock(amrex::Real zref, amrex::Real flux, amrex::Real cnk_a, bool cnk_visc)
Definition: ERF_MOSTStress.H:957
similarity_funs sfuns
Definition: ERF_MOSTStress.H:1054
Definition: ERF_MOSTStress.H:1164
similarity_funs sfuns
Definition: ERF_MOSTStress.H:1251
most_data mdata
Definition: ERF_MOSTStress.H:1250
surface_temp_donelan(amrex::Real zref, amrex::Real flux)
Definition: ERF_MOSTStress.H:1165
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< 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 > &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_MOSTStress.H:1175
const amrex::Real tol
Definition: ERF_MOSTStress.H:1252
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1253
Definition: ERF_MOSTStress.H:1064
similarity_funs sfuns
Definition: ERF_MOSTStress.H:1154
surface_temp_mod_charnock(amrex::Real zref, amrex::Real flux, amrex::Real depth)
Definition: ERF_MOSTStress.H:1065
const amrex::Real tol
Definition: ERF_MOSTStress.H:1155
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1156
most_data mdata
Definition: ERF_MOSTStress.H:1153
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< 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 > &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_MOSTStress.H:1078
Definition: ERF_MOSTStress.H:1261
const amrex::Real eps
Definition: ERF_MOSTStress.H:1356
const amrex::Real tol
Definition: ERF_MOSTStress.H:1355
most_data mdata
Definition: ERF_MOSTStress.H:1353
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< 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 > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &Hwave_arr, const amrex::Array4< amrex::Real > &Lwave_arr, const amrex::Array4< amrex::Real > &eta_arr) const
Definition: ERF_MOSTStress.H:1272
const amrex::Real z0_eps
Definition: ERF_MOSTStress.H:1357
similarity_funs sfuns
Definition: ERF_MOSTStress.H:1354
const amrex::Real z0_max
Definition: ERF_MOSTStress.H:1358
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1359
surface_temp_wave_coupled(amrex::Real zref, amrex::Real flux)
Definition: ERF_MOSTStress.H:1262
Definition: ERF_MOSTStress.H:862
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 > &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 > &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_MOSTStress.H:873
similarity_funs sfuns
Definition: ERF_MOSTStress.H:946
surface_temp(amrex::Real zref, amrex::Real flux)
Definition: ERF_MOSTStress.H:863
const amrex::Real tol
Definition: ERF_MOSTStress.H:947
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:948
most_data mdata
Definition: ERF_MOSTStress.H:945