ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
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  amrex::Real surf_moist_flux{0.0}; ///< Moisture flux
21 
22  amrex::Real Cnk_a{0.0185}; ///< Standard Charnock constant https://doi.org/10.1175/JAMC-D-17-0137.1
23  amrex::Real Cnk_b1{1.0/30.0}; ///< Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1
24  amrex::Real Cnk_b2{1260.0}; ///< Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1
25  amrex::Real Cnk_d{30.0}; ///< Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1
27  bool visc{false}; ///< Use viscous Charnock formulation
28 
29  const amrex::Real Bjr_beta = 1.2; // Empirical parameter from Beljaars 1995 QJRMS
30 };
31 
32 
33 /**
34  * Structure of similarity functions for Moeng formulation
35  */
37 {
38  //
39  // Similarity functions from Jimenez (2012)
40  //
41  AMREX_GPU_HOST_DEVICE
42  AMREX_FORCE_INLINE
44  calc_psi_m2 (amrex::Real zeta) const
45  {
46  if (zeta > 0) {
47  amrex::Real x = std::pow(1.0 + std::pow(zeta, 2.5), 1.0/2.5);
48  return ( -6.1*std::log(zeta + x) );
49  } else {
50  amrex::Real x = std::pow(1.0 - 16.0*zeta, 0.25);
51  amrex::Real psi_k_m = 2.0 * std::log(0.5 * (1.0 + x)) + log(0.5 * (1.0 + x * x)) -
52  2.0 * std::atan(x) + PIoTwo;
53  amrex::Real y = std::pow(1.0 - 10.0*zeta, 1.0/3.0);
54  amrex::Real psi_c_m = (3.0/2.0)*std::log((y*y + y + 1.0)/3.0)
55  - std::sqrt(3.0)*std::atan((2.0*y + 1.0)/std::sqrt(3.0))
56  + PI/std::sqrt(3.0);
57  return ( (psi_k_m + zeta*zeta*psi_c_m) / (1. + zeta*zeta) );
58  }
59  }
60 
61  AMREX_GPU_HOST_DEVICE
62  AMREX_FORCE_INLINE
64  calc_psi_h2 (amrex::Real zeta) const
65  {
66  if (zeta > 0) {
67  amrex::Real x = std::pow(1.0 + std::pow(zeta, 1.1), 1.0/1.1);
68  return ( -6.1*std::log(zeta + x) );
69  } else {
70  amrex::Real x = std::sqrt(1.0 - 16.0*zeta);
71  amrex::Real psi_k_h = 2.0 * std::log(0.5 * (1.0 + x));
72  amrex::Real y = std::pow(1.0 - 34.0*zeta, 1.0/3.0);
73  amrex::Real psi_c_h = (3.0/2.0)*std::log((y*y + y + 1.0)/3.0)
74  - std::sqrt(3.0)*std::atan((2.0*y + 1.0)/std::sqrt(3.0))
75  + PI/std::sqrt(3.0);
76  return ( (psi_k_h + zeta*zeta*psi_c_h) / (1. + zeta*zeta) );
77  }
78  }
79 
80 
81  //
82  // Similarity functions from Businger & Dyer (1966)
83  //
84  AMREX_GPU_HOST_DEVICE
85  AMREX_FORCE_INLINE
87  calc_psi_m (amrex::Real zeta) const
88  {
89  if (zeta > 0) {
90  return -beta_m * zeta;
91  } else {
92  amrex::Real x = std::sqrt(std::sqrt(1.0 - gamma_m * zeta));
93  return 2.0 * std::log(0.5 * (1.0 + x)) + log(0.5 * (1.0 + x * x)) -
94  2.0 * std::atan(x) + PIoTwo;
95  }
96  }
97 
98  AMREX_GPU_HOST_DEVICE
99  AMREX_FORCE_INLINE
101  calc_psi_h (amrex::Real zeta) const
102  {
103  if (zeta > 0) {
104  return -beta_h * zeta;
105  } else {
106  amrex::Real x = std::sqrt(1.0 - gamma_h * zeta);
107  return 2.0 * std::log(0.5 * (1.0 + x));
108  }
109  }
110 
111 private:
112  amrex::Real beta_m{5.0}; ///< Constants from Dyer, BLM, 1974
113  amrex::Real beta_h{5.0}; ///< https://doi.org/10.1007/BF00240838
116 };
117 
118 
119 /**
120  * Empirical kinematic viscosity [m2/s] formula from Andreas (1989) CRREL Rep.
121  * 89-11, valid between -173 and 277 deg C.
122  */
123 AMREX_GPU_DEVICE
124 AMREX_FORCE_INLINE
127 {
128  amrex::Real TC = T_degK - 273.15;
129  return 1.326e-5*(1. + 6.542e-3*TC + 8.301e-6*TC*TC - 4.84e-9*TC*TC*TC);
130 }
131 
132 
133 /**
134  * Adiabatic with constant roughness
135  */
136 struct adiabatic
137 {
139  amrex::Real Tflux,
140  amrex::Real Qvflux)
141  {
142  mdata.zref = zref;
143  mdata.surf_temp_flux = Tflux;
144  mdata.surf_moist_flux = Qvflux;
145  }
146 
147  AMREX_GPU_DEVICE
148  AMREX_FORCE_INLINE
149  void
150  iterate_flux (const int& i,
151  const int& j,
152  const int& k,
153  const int& /*max_iters*/,
154  const amrex::Array4<const amrex::Real>& z0_arr,
155  const amrex::Array4<const amrex::Real>& umm_arr,
156  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
157  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
158  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
159  const amrex::Array4<amrex::Real>& u_star_arr,
160  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
161  const amrex::Array4<amrex::Real>& t_star_arr,
162  const amrex::Array4<amrex::Real>& q_star_arr,
163  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
164  const amrex::Array4<amrex::Real>& /*q_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  olen_arr(i,j,k) = 1.0e16;
172  u_star_arr(i,j,k) = mdata.kappa * umm_arr(i,j,k) / std::log(mdata.zref / z0_arr(i,j,k));
173  t_star_arr(i,j,k) = 0.0;
174  q_star_arr(i,j,k) = 0.0;
175  }
176 
177 private:
180 };
181 
182 
183 /**
184  * Adiabatic with charnock roughness
185  */
187 {
189  amrex::Real Tflux,
190  amrex::Real Qvflux,
191  amrex::Real cnk_a,
192  bool cnk_visc)
193  {
194  mdata.zref = zref;
195  mdata.surf_temp_flux = Tflux;
196  mdata.surf_moist_flux = Qvflux;
197  mdata.Cnk_a = cnk_a;
198  mdata.visc = cnk_visc;
199  }
200 
201  AMREX_GPU_DEVICE
202  AMREX_FORCE_INLINE
203  void
204  iterate_flux (const int& i,
205  const int& j,
206  const int& k,
207  const int& max_iters,
208  const amrex::Array4<amrex::Real>& z0_arr,
209  const amrex::Array4<const amrex::Real>& umm_arr,
210  const amrex::Array4<const amrex::Real>& tm_arr,
211  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
212  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
213  const amrex::Array4<amrex::Real>& u_star_arr,
214  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
215  const amrex::Array4<amrex::Real>& t_star_arr,
216  const amrex::Array4<amrex::Real>& q_star_arr,
217  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
218  const amrex::Array4<amrex::Real>& /*q_surf_arr*/,
219  const amrex::Array4<amrex::Real>& olen_arr,
220  const amrex::Array4<amrex::Real>& /*pblh_arr*/,
221  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
222  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
223  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
224  {
225  amrex::Real ustar = 0.0;
226  amrex::Real z0 = z0_arr(i,j,k);
227  amrex::Real z0_old = z0;
228  amrex::Real psi_m = 0.0;
229  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
230  amrex::Real C = std::log(mdata.zref / z0);
231 
232  // Fixed point iteration on roughness
233  int iter_z = 0;
234  do {
235  // Transfer curr to old
236  z0_old = z0;
237 
238  // Update
239  C = std::log(mdata.zref / z0_old);
240  ustar = mdata.kappa * umm / (C - psi_m);
241  if (mdata.Cnk_a > 0) {
242  z0 = (mdata.Cnk_a / mdata.gravity) * ustar * ustar;
243  if (mdata.visc) {
244  z0 += air_viscosity(tm_arr(i,j,k)) / std::max(ustar, 0.05);
245  }
246  } else {
247  z0 = COARE3_roughness(mdata.zref, umm, ustar);
248  }
249 
250  ++iter_z;
251  } while ( (std::abs(z0 - z0_old) > tol_z) && (iter_z <= max_iters) );
252  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter_z < max_iters,
253  "Maximum number of MOST roughness iterations reached.");
254  C = std::log(mdata.zref / z0);
255 
256  // Populate the stored MOST arrays
257  z0_arr(i,j,k) = z0;
258  olen_arr(i,j,k) = 1.0e16;
259  u_star_arr(i,j,k) = mdata.kappa * umm / (C - psi_m);
260  t_star_arr(i,j,k) = 0.0;
261  q_star_arr(i,j,k) = 0.0;
262  }
263 
264 private:
267  const amrex::Real tol_z = 1.0e-10;
268  const amrex::Real WSMIN = 0.1; // minimum wind speed
269 };
270 
271 
272 /**
273  * Adiabatic with modified charnock roughness
274  */
276 {
278  amrex::Real Tflux,
279  amrex::Real Qvflux,
280  amrex::Real depth)
281  {
282  mdata.zref = zref;
283  mdata.surf_temp_flux = Tflux;
284  mdata.surf_moist_flux = Qvflux;
285  mdata.Cnk_d = depth;
286  mdata.Cnk_b = mdata.Cnk_b1 * std::log(mdata.Cnk_b2 / mdata.Cnk_d);
287  }
288 
289  AMREX_GPU_DEVICE
290  AMREX_FORCE_INLINE
291  void
292  iterate_flux (const int& i,
293  const int& j,
294  const int& k,
295  const int& max_iters,
296  const amrex::Array4<amrex::Real>& z0_arr,
297  const amrex::Array4<const amrex::Real>& umm_arr,
298  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
299  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
300  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
301  const amrex::Array4<amrex::Real>& u_star_arr,
302  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
303  const amrex::Array4<amrex::Real>& t_star_arr,
304  const amrex::Array4<amrex::Real>& q_star_arr,
305  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
306  const amrex::Array4<amrex::Real>& /*q_surf_arr*/,
307  const amrex::Array4<amrex::Real>& olen_arr,
308  const amrex::Array4<amrex::Real>& /*pblh_arr*/,
309  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
310  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
311  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
312  {
313  amrex::Real ustar = 0.0;
314  amrex::Real z0 = z0_arr(i,j,k);
315  amrex::Real z0_old = z0;
316  amrex::Real psi_m = 0.0;
317  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
318  amrex::Real C = std::log(mdata.zref / z0);
319 
320  // Fixed point iteration on roughness
321  int iter_z = 0;
322  do {
323  // Transfer curr to old
324  z0_old = z0;
325 
326  // Update
327  C = std::log(mdata.zref / z0_old);
328  ustar = mdata.kappa * umm / (C - psi_m);
329  z0 = std::exp( (2.7*ustar - 1.8/mdata.Cnk_b) / (ustar + 0.17/mdata.Cnk_b) );
330 
331  ++iter_z;
332  } while ( (std::abs(z0 - z0_old) > tol_z) && (iter_z <= max_iters) );
333  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter_z < max_iters,
334  "Maximum number of MOST roughness iterations reached.");
335  C = std::log(mdata.zref / z0);
336 
337  // Populate the stored MOST arrays
338  z0_arr(i,j,k) = z0;
339  olen_arr(i,j,k) = 1.0e16;
340  u_star_arr(i,j,k) = mdata.kappa * umm / (C - psi_m);
341  t_star_arr(i,j,k) = 0.0;
342  q_star_arr(i,j,k) = 0.0;
343  }
344 
345 private:
348  const amrex::Real tol_z = 1.0e-10;
349  const amrex::Real WSMIN = 0.1; // minimum wind speed
350 };
351 
352 
353 /**
354  * Adiabatic with Donelan roughness
355  */
357 {
359  amrex::Real Tflux,
360  amrex::Real Qvflux)
361  {
362  mdata.zref = zref;
363  mdata.surf_temp_flux = Tflux;
364  mdata.surf_moist_flux = Qvflux;
365  }
366 
367  AMREX_GPU_DEVICE
368  AMREX_FORCE_INLINE
369  void
370  iterate_flux (const int& i,
371  const int& j,
372  const int& k,
373  const int& max_iters,
374  const amrex::Array4<amrex::Real>& z0_arr,
375  const amrex::Array4<const amrex::Real>& umm_arr,
376  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
377  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
378  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
379  const amrex::Array4<amrex::Real>& u_star_arr,
380  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
381  const amrex::Array4<amrex::Real>& t_star_arr,
382  const amrex::Array4<amrex::Real>& q_star_arr,
383  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
384  const amrex::Array4<amrex::Real>& /*q_surf_arr*/,
385  const amrex::Array4<amrex::Real>& olen_arr,
386  const amrex::Array4<amrex::Real>& /*pblh_arr*/,
387  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
388  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
389  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
390  {
391  int iter = 0;
392  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
393  amrex::Real ustar = 0.0;
394  amrex::Real z0 = 0.0;
395  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
396  if (u_star_arr(i,j,k) == 1.E34) {
397  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
398  }
399  do {
400  ustar = u_star_arr(i,j,k);
401  z0 = Donelan_roughness(ustar);
402  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0);
403  ++iter;
404  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
405  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
406  "Maximum number of MOST iterations reached.");
407 
408  // Populate the stored MOST arrays
409  z0_arr(i,j,k) = z0;
410  olen_arr(i,j,k) = 1.0e16;
411  t_star_arr(i,j,k) = 0.0;
412  q_star_arr(i,j,k) = 0.0;
413  }
414 
415 private:
418  const amrex::Real tol = 1.0e-5;
419  const amrex::Real WSMIN = 0.1; // minimum wind speed
420 };
421 
422 
423 /**
424  * Adiabatic with wave-coupled roughness
425  */
427 {
429  amrex::Real Tflux,
430  amrex::Real Qvflux)
431  {
432  mdata.zref = zref;
433  mdata.surf_temp_flux = Tflux;
434  mdata.surf_moist_flux = Qvflux;
435  }
436 
437  AMREX_GPU_DEVICE
438  AMREX_FORCE_INLINE
439  void
440  iterate_flux (const int& i,
441  const int& j,
442  const int& k,
443  const int& max_iters,
444  const amrex::Array4<amrex::Real>& z0_arr,
445  const amrex::Array4<const amrex::Real>& umm_arr,
446  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
447  const amrex::Array4<const amrex::Real>& /*tvm_arr*/,
448  const amrex::Array4<const amrex::Real>& /*qvm_arr*/,
449  const amrex::Array4<amrex::Real>& u_star_arr,
450  const amrex::Array4<amrex::Real>& /*w_star_arr*/,
451  const amrex::Array4<amrex::Real>& t_star_arr,
452  const amrex::Array4<amrex::Real>& q_star_arr,
453  const amrex::Array4<amrex::Real>& /*t_surf_arr*/,
454  const amrex::Array4<amrex::Real>& /*q_surf_arr*/,
455  const amrex::Array4<amrex::Real>& olen_arr,
456  const amrex::Array4<amrex::Real>& /*pblh_arr*/,
457  const amrex::Array4<amrex::Real>& Hwave_arr,
458  const amrex::Array4<amrex::Real>& Lwave_arr,
459  const amrex::Array4<amrex::Real>& eta_arr) const
460  {
461  int iter = 0;
462  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
463  amrex::Real ustar = 0.0;
464  amrex::Real z0 = 0.0;
465  int ie, je;
466  ie = i < lbound(eta_arr).x ? lbound(eta_arr).x : i;
467  je = j < lbound(eta_arr).y ? lbound(eta_arr).y : j;
468  ie = ie > ubound(eta_arr).x ? ubound(eta_arr).x : ie;
469  je = je > ubound(eta_arr).y ? ubound(eta_arr).y : je;
470  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
471  do {
472  ustar = u_star_arr(i,j,k);
473  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 )
474  + 0.11 * eta_arr(ie,je,k,EddyDiff::Mom_v) / ustar, z0_eps), z0_max );
475  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0);
476  ++iter;
477  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
478  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
479  "Maximum number of MOST iterations reached.");
480 
481  // Populate the stored MOST arrays
482  z0_arr(i,j,k) = z0;
483  olen_arr(i,j,k) = 1.0e16;
484  t_star_arr(i,j,k) = 0.0;
485  q_star_arr(i,j,k) = 0.0;
486  }
487 
488 private:
491  const amrex::Real tol = 1.0e-5;
492  const amrex::Real eps = 1e-15;
493  const amrex::Real z0_eps = 1.0e-6;
494  const amrex::Real z0_max = 0.1;
495  const amrex::Real WSMIN = 0.1; // minimum wind speed
496 };
497 
498 
499 /**
500  * Surface flux with constant roughness
501  */
503 {
505  amrex::Real Tflux,
506  amrex::Real Qvflux,
507  bool cons_qflux)
508  {
509  mdata.zref = zref;
510  mdata.surf_temp_flux = Tflux;
511  mdata.surf_moist_flux = Qvflux;
512  spec_qflux = cons_qflux;
513  }
514 
515  AMREX_GPU_DEVICE
516  AMREX_FORCE_INLINE
517  void
518  iterate_flux (const int& i,
519  const int& j,
520  const int& k,
521  const int& max_iters,
522  const amrex::Array4<const amrex::Real>& z0_arr,
523  const amrex::Array4<const amrex::Real>& umm_arr,
524  const amrex::Array4<const amrex::Real>& tm_arr,
525  const amrex::Array4<const amrex::Real>& tvm_arr,
526  const amrex::Array4<const amrex::Real>& qvm_arr,
527  const amrex::Array4<amrex::Real>& u_star_arr,
528  const amrex::Array4<amrex::Real>& w_star_arr,
529  const amrex::Array4<amrex::Real>& t_star_arr,
530  const amrex::Array4<amrex::Real>& q_star_arr,
531  const amrex::Array4<amrex::Real>& t_surf_arr,
532  const amrex::Array4<amrex::Real>& q_surf_arr,
533  const amrex::Array4<amrex::Real>& olen_arr,
534  const amrex::Array4<amrex::Real>& pblh_arr,
535  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
536  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
537  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
538  {
539  int iter = 0;
540  amrex::Real ustar = 0.0;
541  amrex::Real wstar = 0.0;
542  amrex::Real tflux = 0.0;
543  amrex::Real qflux = 0.0;
544  amrex::Real zeta = 0.0;
545  amrex::Real psi_m = 0.0;
546  amrex::Real psi_h = 0.0;
547  amrex::Real Olen = 0.0;
548  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
549  if (u_star_arr(i,j,k) == 1.E34) {
550  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
551  } else {
552  Olen = olen_arr(i,j,k);
553  zeta = mdata.zref / Olen;
554  psi_m = sfuns.calc_psi_m(zeta);
555  psi_h = sfuns.calc_psi_h(zeta);
556  }
557  do {
558  ustar = u_star_arr(i,j,k);
559  qflux = (spec_qflux) ? mdata.surf_moist_flux :
560  -(qvm_arr(i,j,k) - q_surf_arr(i,j,k)) * ustar * mdata.kappa /
561  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h); // <w'Qv'>
562  tflux = mdata.surf_temp_flux*(1 + 0.61*qvm_arr(i,j,k)) + qflux*0.61*tm_arr(i,j,k);
563  if (w_star_arr) {
564  // update w* and Umagmean
565  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
566  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
567  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
568  umm = std::max(umm, WSMIN);
569  }
570  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
571  zeta = mdata.zref / Olen;
572  psi_m = sfuns.calc_psi_m(zeta);
573  psi_h = sfuns.calc_psi_h(zeta);
574  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0_arr(i,j,k)) - psi_m);
575  ++iter;
576  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
577  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
578  "Maximum number of MOST iterations reached.");
579 
580  // Populate the stored MOST arrays
581  olen_arr(i,j,k) = Olen;
582  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h) /
583  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,k);
584  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
585  if (spec_qflux) {
586  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h) /
587  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,k);
588  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
589  } else {
590  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,k) - q_surf_arr(i,j,k)) /
591  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h);
592  }
593  }
594 
595 private:
599  const amrex::Real tol = 1.0e-5;
600  const amrex::Real WSMIN = 0.1; // minimum wind speed
601 };
602 
603 
604 /**
605  * Surface flux with charnock roughness
606  */
608 {
610  amrex::Real Tflux,
611  amrex::Real Qvflux,
612  amrex::Real cnk_a,
613  bool cnk_visc,
614  bool cons_qflux)
615  {
616  mdata.zref = zref;
617  mdata.surf_temp_flux = Tflux;
618  mdata.surf_moist_flux = Qvflux;
619  mdata.Cnk_a = cnk_a;
620  mdata.visc = cnk_visc;
621  spec_qflux = cons_qflux;
622  }
623 
624  AMREX_GPU_DEVICE
625  AMREX_FORCE_INLINE
626  void
627  iterate_flux (const int& i,
628  const int& j,
629  const int& k,
630  const int& max_iters,
631  const amrex::Array4<amrex::Real>& z0_arr,
632  const amrex::Array4<const amrex::Real>& umm_arr,
633  const amrex::Array4<const amrex::Real>& tm_arr,
634  const amrex::Array4<const amrex::Real>& tvm_arr,
635  const amrex::Array4<const amrex::Real>& qvm_arr,
636  const amrex::Array4<amrex::Real>& u_star_arr,
637  const amrex::Array4<amrex::Real>& w_star_arr,
638  const amrex::Array4<amrex::Real>& t_star_arr,
639  const amrex::Array4<amrex::Real>& q_star_arr,
640  const amrex::Array4<amrex::Real>& t_surf_arr,
641  const amrex::Array4<amrex::Real>& q_surf_arr,
642  const amrex::Array4<amrex::Real>& olen_arr,
643  const amrex::Array4<amrex::Real>& pblh_arr,
644  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
645  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
646  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
647  {
648  int iter = 0;
649  amrex::Real ustar = 0.0;
650  amrex::Real wstar = 0.0;
651  amrex::Real tflux = 0.0;
652  amrex::Real qflux = 0.0;
653  amrex::Real z0 = 0.0;
654  amrex::Real zeta = 0.0;
655  amrex::Real psi_m = 0.0;
656  amrex::Real psi_h = 0.0;
657  amrex::Real Olen = 0.0;
658  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
659  if (u_star_arr(i,j,k) == 1.E34) {
660  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
661  } else {
662  Olen = olen_arr(i,j,k);
663  zeta = mdata.zref / Olen;
664  psi_m = sfuns.calc_psi_m(zeta);
665  psi_h = sfuns.calc_psi_h(zeta);
666  }
667  do {
668  ustar = u_star_arr(i,j,k);
669  if (mdata.Cnk_a > 0) {
670  z0 = (mdata.Cnk_a / mdata.gravity) * ustar * ustar;
671  if (mdata.visc) {
672  z0 += air_viscosity(tm_arr(i,j,k)) / std::max(ustar, 0.05);
673  }
674  } else {
675  z0 = COARE3_roughness(mdata.zref, umm, ustar);
676  }
677  qflux = (spec_qflux) ? mdata.surf_moist_flux :
678  -(qvm_arr(i,j,k) - q_surf_arr(i,j,k)) * ustar * mdata.kappa /
679  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h); // <w'Qv'>
680  tflux = mdata.surf_temp_flux*(1 + 0.61*qvm_arr(i,j,k)) + qflux*0.61*tm_arr(i,j,k);
681  if (w_star_arr) {
682  // update w* and Umagmean
683  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
684  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
685  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
686  umm = std::max(umm, WSMIN);
687  }
688  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
689  zeta = mdata.zref / Olen;
690  psi_m = sfuns.calc_psi_m(zeta);
691  psi_h = sfuns.calc_psi_h(zeta);
692  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
693  ++iter;
694  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
695  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
696  "Maximum number of MOST iterations reached.");
697 
698  // Populate the stored MOST arrays
699  z0_arr(i,j,k) = z0;
700  olen_arr(i,j,k) = Olen;
701  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(mdata.zref / z0) - psi_h) /
702  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,k);
703  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
704  if (spec_qflux) {
705  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h) /
706  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,k);
707  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
708  } else {
709  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,k) - q_surf_arr(i,j,k)) /
710  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h);
711  }
712  }
713 
714 private:
718  const amrex::Real tol = 1.0e-5;
719  const amrex::Real WSMIN = 0.1; // minimum wind speed
720 };
721 
722 
723 /**
724  * Surface flux with modified charnock roughness
725  */
727 {
729  amrex::Real Tflux,
730  amrex::Real Qvflux,
731  amrex::Real depth,
732  bool cons_qflux)
733  {
734  mdata.zref = zref;
735  mdata.surf_temp_flux = Tflux;
736  mdata.surf_moist_flux = Qvflux;
737  mdata.Cnk_d = depth;
738  mdata.Cnk_b = mdata.Cnk_b1 * std::log(mdata.Cnk_b2 / mdata.Cnk_d);
739  spec_qflux = cons_qflux;
740  }
741 
742  AMREX_GPU_DEVICE
743  AMREX_FORCE_INLINE
744  void
745  iterate_flux (const int& i,
746  const int& j,
747  const int& k,
748  const int& max_iters,
749  const amrex::Array4<amrex::Real>& z0_arr,
750  const amrex::Array4<const amrex::Real>& umm_arr,
751  const amrex::Array4<const amrex::Real>& tm_arr,
752  const amrex::Array4<const amrex::Real>& tvm_arr,
753  const amrex::Array4<const amrex::Real>& qvm_arr,
754  const amrex::Array4<amrex::Real>& u_star_arr,
755  const amrex::Array4<amrex::Real>& w_star_arr,
756  const amrex::Array4<amrex::Real>& t_star_arr,
757  const amrex::Array4<amrex::Real>& q_star_arr,
758  const amrex::Array4<amrex::Real>& t_surf_arr,
759  const amrex::Array4<amrex::Real>& q_surf_arr,
760  const amrex::Array4<amrex::Real>& olen_arr,
761  const amrex::Array4<amrex::Real>& pblh_arr,
762  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
763  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
764  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
765  {
766  int iter = 0;
767  amrex::Real ustar = 0.0;
768  amrex::Real wstar = 0.0;
769  amrex::Real tflux = 0.0;
770  amrex::Real qflux = 0.0;
771  amrex::Real z0 = 0.0;
772  amrex::Real zeta = 0.0;
773  amrex::Real psi_m = 0.0;
774  amrex::Real psi_h = 0.0;
775  amrex::Real Olen = 0.0;
776  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
777  if (u_star_arr(i,j,k) == 1.E34) {
778  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
779  } else {
780  Olen = olen_arr(i,j,k);
781  zeta = mdata.zref / Olen;
782  psi_m = sfuns.calc_psi_m(zeta);
783  psi_h = sfuns.calc_psi_h(zeta);
784  }
785  do {
786  ustar = u_star_arr(i,j,k);
787  z0 = std::exp( (2.7*ustar - 1.8/mdata.Cnk_b) / (ustar + 0.17/mdata.Cnk_b) );
788  qflux = (spec_qflux) ? mdata.surf_moist_flux :
789  -(qvm_arr(i,j,k) - q_surf_arr(i,j,k)) * ustar * mdata.kappa /
790  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h); // <w'Qv'>
791  tflux = mdata.surf_temp_flux*(1 + 0.61*qvm_arr(i,j,k)) + qflux*0.61*tm_arr(i,j,k);
792  if (w_star_arr) {
793  // update w* and Umagmean
794  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
795  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
796  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
797  umm = std::max(umm, WSMIN);
798  }
799  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
800  zeta = mdata.zref / Olen;
801  psi_m = sfuns.calc_psi_m(zeta);
802  psi_h = sfuns.calc_psi_h(zeta);
803  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
804  ++iter;
805  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
806  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
807  "Maximum number of MOST iterations reached.");
808 
809  // Populate the stored MOST arrays
810  z0_arr(i,j,k) = z0;
811  olen_arr(i,j,k) = Olen;
812  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(mdata.zref / z0) - psi_h) /
813  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,k);
814  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
815  if (spec_qflux) {
816  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h) /
817  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,k);
818  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
819  } else {
820  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,k) - q_surf_arr(i,j,k)) /
821  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h);
822  }
823  }
824 
825 private:
829  const amrex::Real tol = 1.0e-5;
830  const amrex::Real WSMIN = 0.1; // minimum wind speed
831 };
832 
833 
834 /**
835  * Surface flux with donelan roughness
836  */
838 {
840  amrex::Real Tflux,
841  amrex::Real Qvflux,
842  bool cons_qflux)
843  {
844  mdata.zref = zref;
845  mdata.surf_temp_flux = Tflux;
846  mdata.surf_moist_flux = Qvflux;
847  spec_qflux = cons_qflux;
848  }
849 
850  AMREX_GPU_DEVICE
851  AMREX_FORCE_INLINE
852  void
853  iterate_flux (const int& i,
854  const int& j,
855  const int& k,
856  const int& max_iters,
857  const amrex::Array4<amrex::Real>& z0_arr,
858  const amrex::Array4<const amrex::Real>& umm_arr,
859  const amrex::Array4<const amrex::Real>& tm_arr,
860  const amrex::Array4<const amrex::Real>& tvm_arr,
861  const amrex::Array4<const amrex::Real>& qvm_arr,
862  const amrex::Array4<amrex::Real>& u_star_arr,
863  const amrex::Array4<amrex::Real>& w_star_arr,
864  const amrex::Array4<amrex::Real>& t_star_arr,
865  const amrex::Array4<amrex::Real>& q_star_arr,
866  const amrex::Array4<amrex::Real>& t_surf_arr,
867  const amrex::Array4<amrex::Real>& q_surf_arr,
868  const amrex::Array4<amrex::Real>& olen_arr,
869  const amrex::Array4<amrex::Real>& pblh_arr,
870  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
871  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
872  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
873  {
874  int iter = 0;
875  amrex::Real ustar = 0.0;
876  amrex::Real wstar = 0.0;
877  amrex::Real tflux = 0.0;
878  amrex::Real qflux = 0.0;
879  amrex::Real z0 = 0.0;
880  amrex::Real zeta = 0.0;
881  amrex::Real psi_m = 0.0;
882  amrex::Real psi_h = 0.0;
883  amrex::Real Olen = 0.0;
884  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
885  if (u_star_arr(i,j,k) == 1.E34) {
886  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
887  } else {
888  Olen = olen_arr(i,j,k);
889  zeta = mdata.zref / Olen;
890  psi_m = sfuns.calc_psi_m(zeta);
891  psi_h = sfuns.calc_psi_h(zeta);
892  }
893  do {
894  ustar = u_star_arr(i,j,k);
895  z0 = Donelan_roughness(ustar);
896  qflux = (spec_qflux) ? mdata.surf_moist_flux :
897  -(qvm_arr(i,j,k) - q_surf_arr(i,j,k)) * ustar * mdata.kappa /
898  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h); // <w'Qv'>
899  tflux = mdata.surf_temp_flux*(1 + 0.61*qvm_arr(i,j,k)) + qflux*0.61*tm_arr(i,j,k);
900  if (w_star_arr) {
901  // update w* and Umagmean
902  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
903  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
904  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
905  umm = std::max(umm, WSMIN);
906  }
907  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
908  zeta = mdata.zref / Olen;
909  psi_m = sfuns.calc_psi_m(zeta);
910  psi_h = sfuns.calc_psi_h(zeta);
911  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
912  ++iter;
913  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
914  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
915  "Maximum number of MOST iterations reached.");
916 
917  // Populate the stored MOST arrays
918  z0_arr(i,j,k) = z0;
919  olen_arr(i,j,k) = Olen;
920  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(mdata.zref / z0) - psi_h) /
921  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,k);
922  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
923  if (spec_qflux) {
924  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h) /
925  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,k);
926  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
927  } else {
928  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,k) - q_surf_arr(i,j,k)) /
929  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h);
930  }
931  }
932 
933 private:
937  const amrex::Real tol = 1.0e-5;
938  const amrex::Real WSMIN = 0.1; // minimum wind speed
939 };
940 
941 
942 /**
943  * Surface flux with wave-coupled roughness
944  */
946 {
948  amrex::Real Tflux,
949  amrex::Real Qvflux,
950  bool cons_qflux)
951  {
952  mdata.zref = zref;
953  mdata.surf_temp_flux = Tflux;
954  mdata.surf_moist_flux = Qvflux;
955  spec_qflux = cons_qflux;
956  }
957 
958  AMREX_GPU_DEVICE
959  AMREX_FORCE_INLINE
960  void
961  iterate_flux (const int& i,
962  const int& j,
963  const int& k,
964  const int& max_iters,
965  const amrex::Array4<amrex::Real>& z0_arr,
966  const amrex::Array4<const amrex::Real>& umm_arr,
967  const amrex::Array4<const amrex::Real>& tm_arr,
968  const amrex::Array4<const amrex::Real>& tvm_arr,
969  const amrex::Array4<const amrex::Real>& qvm_arr,
970  const amrex::Array4<amrex::Real>& u_star_arr,
971  const amrex::Array4<amrex::Real>& w_star_arr,
972  const amrex::Array4<amrex::Real>& t_star_arr,
973  const amrex::Array4<amrex::Real>& q_star_arr,
974  const amrex::Array4<amrex::Real>& t_surf_arr,
975  const amrex::Array4<amrex::Real>& q_surf_arr,
976  const amrex::Array4<amrex::Real>& olen_arr,
977  const amrex::Array4<amrex::Real>& pblh_arr,
978  const amrex::Array4<amrex::Real>& Hwave_arr,
979  const amrex::Array4<amrex::Real>& Lwave_arr,
980  const amrex::Array4<amrex::Real>& eta_arr) const
981  {
982  int iter = 0;
983  amrex::Real ustar = 0.0;
984  amrex::Real wstar = 0.0;
985  amrex::Real tflux = 0.0;
986  amrex::Real qflux = 0.0;
987  amrex::Real z0 = 0.0;
988  amrex::Real zeta = 0.0;
989  amrex::Real psi_m = 0.0;
990  amrex::Real psi_h = 0.0;
991  amrex::Real Olen = 0.0;
992  int ie, je;
993  ie = i < lbound(eta_arr).x ? lbound(eta_arr).x : i;
994  je = j < lbound(eta_arr).y ? lbound(eta_arr).y : j;
995  ie = ie > ubound(eta_arr).x ? ubound(eta_arr).x : ie;
996  je = je > ubound(eta_arr).y ? ubound(eta_arr).y : je;
997  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
998  if (u_star_arr(i,j,k) == 1.E34) {
999  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
1000  } else {
1001  Olen = olen_arr(i,j,k);
1002  zeta = mdata.zref / Olen;
1003  psi_m = sfuns.calc_psi_m(zeta);
1004  psi_h = sfuns.calc_psi_h(zeta);
1005  }
1006  do {
1007  ustar = u_star_arr(i,j,k);
1008  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 )
1009  + 0.11 * eta_arr(ie,je,k,EddyDiff::Mom_v) / ustar, z0_eps), z0_max );
1010  qflux = (spec_qflux) ? mdata.surf_moist_flux :
1011  -(qvm_arr(i,j,k) - q_surf_arr(i,j,k)) * ustar * mdata.kappa /
1012  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h); // <w'Qv'>
1013  tflux = mdata.surf_temp_flux*(1 + 0.61*qvm_arr(i,j,k)) + qflux*0.61*tm_arr(i,j,k);
1014  if (w_star_arr) {
1015  // update w* and Umagmean
1016  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
1017  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
1018  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
1019  umm = std::max(umm, WSMIN);
1020  }
1021  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
1022  zeta = mdata.zref / Olen;
1023  psi_m = sfuns.calc_psi_m(zeta);
1024  psi_h = sfuns.calc_psi_h(zeta);
1025  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
1026  ++iter;
1027  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
1028  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
1029  "Maximum number of MOST iterations reached.");
1030 
1031  // Populate the stored MOST arrays
1032  z0_arr(i,j,k) = z0;
1033  olen_arr(i,j,k) = Olen;
1034  t_surf_arr(i,j,k) = mdata.surf_temp_flux * (std::log(mdata.zref / z0) - psi_h) /
1035  (u_star_arr(i,j,k) * mdata.kappa) + tm_arr(i,j,k);
1036  t_star_arr(i,j,k) = -mdata.surf_temp_flux / u_star_arr(i,j,k);
1037  if (spec_qflux) {
1038  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h) /
1039  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,k);
1040  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
1041  } else {
1042  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,k) - q_surf_arr(i,j,k)) /
1043  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h);
1044  }
1045  }
1046 
1047 private:
1051  const amrex::Real tol = 1.0e-5;
1052  const amrex::Real eps = 1e-15;
1053  const amrex::Real z0_eps = 1.0e-6;
1054  const amrex::Real z0_max = 0.1;
1055  const amrex::Real WSMIN = 0.1; // minimum wind speed
1056 };
1057 
1058 
1059 /**
1060  * Surface temperature with constant roughness
1061  */
1063 {
1065  amrex::Real Tflux,
1066  amrex::Real Qvflux,
1067  bool cons_qflux)
1068  {
1069  mdata.zref = zref;
1070  mdata.surf_temp_flux = Tflux;
1071  mdata.surf_moist_flux = Qvflux;
1072  spec_qflux = cons_qflux;
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<const 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>& q_surf_arr,
1093  const amrex::Array4<amrex::Real>& olen_arr,
1094  const amrex::Array4<amrex::Real>& pblh_arr,
1095  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
1096  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
1097  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
1098  {
1099  amrex::Real Rib = 0.0;
1100  amrex::Real zeta = 0.0;
1101  amrex::Real zeta_old = 0.0;
1102  amrex::Real psi_m = 0.0;
1103  amrex::Real psi_h = 0.0;
1104  amrex::Real num = 0.0;
1105  amrex::Real den = 0.0;
1106  amrex::Real z0 = z0_arr(i,j,k);
1107  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
1108  amrex::Real C = std::log(mdata.zref / z0);
1109 
1110  // First iteration we assume neutral (L -> inf)
1111  if (u_star_arr(i,j,k) == 1.E34) { olen_arr(i,j,k) = 1.0e3; }
1112  zeta = mdata.zref / olen_arr(i,j,k);
1113 
1114  // Water vapor in atmos and surface
1115  amrex::Real qv_s, qv_a;
1116  if (q_surf_arr(i,j,k) > 0.) {
1117  qv_s = q_surf_arr(i,j,k);
1118  } else {
1119  // First iteration and no qv_surf was specified
1120  // Use mean since there will be no flux
1121  qv_s = qvm_arr(i,j,k);
1122  }
1123  qv_a = qvm_arr(i,j,k);
1124 
1125  // update w* and Umagmean from Beljaars (1995)
1126  if (w_star_arr) {
1127  // NOTE: Thv flux is lagged, similar to WRF
1128  psi_m = sfuns.calc_psi_m2(zeta);
1129  psi_h = sfuns.calc_psi_h2(zeta);
1130  amrex::Real ustar = mdata.kappa * umm / (C - psi_m);
1131  amrex::Real tstar = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) / (C - psi_h);
1133  -ustar * mdata.kappa * (qv_a - qv_s) / (C - psi_h);
1134  amrex::Real tflux = -ustar*tstar*(1 + 0.61*qvm_arr(i,j,k)) + 0.61*tm_arr(i,j,k)*qflux;
1135  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
1136  amrex::Real wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
1137  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
1138  umm = std::max(umm, WSMIN);
1139  }
1140 
1141  // Bulk Richardson number w/ moisture
1142  amrex::Real thv_s = t_surf_arr(i,j,k) * (1.0 + 0.61*qv_s);
1143  amrex::Real thv_a = tm_arr(i,j,k) * (1.0 + 0.61*qv_a);
1144  Rib = ( (mdata.gravity * mdata.zref) / tm_arr(i,j,k) ) *
1145  ( (thv_a - thv_s) / (umm * umm) );
1146  Rib = std::min(std::max(Rib,-4.0),4.0);
1147 
1148  // Fixed point iteration on zeta
1149  int iter = 0;
1150  do {
1151  // Transfer curr to old
1152  zeta_old = zeta;
1153 
1154  // Stability functions
1155  psi_m = sfuns.calc_psi_m2(zeta_old);
1156  psi_h = sfuns.calc_psi_h2(zeta_old);
1157 
1158  // Limiting
1159  num = std::max(C - psi_m, 1.0);
1160  den = std::max(C - psi_h, 1.0);
1161 
1162  // Update with under relaxation
1163  zeta = (1.0 - alpha) * zeta_old + alpha * Rib * num * num / den;
1164 
1165  ++iter;
1166  } while ( (std::abs(zeta - zeta_old) > tol) && (iter <= max_iters) );
1167  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
1168  "Maximum number of MOST iterations reached.");
1169 
1170  // Populate the stored MOST arrays
1171  olen_arr(i,j,k) = mdata.zref / zeta;
1172  u_star_arr(i,j,k) = mdata.kappa * umm / (C - psi_m);
1173  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) / (C - psi_h);
1174  if (spec_qflux) {
1175  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (C - psi_h) /
1176  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,k);
1177  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
1178  } else {
1179  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,k) - q_surf_arr(i,j,k)) / (C - psi_h);
1180  }
1181  }
1182 
1183 private:
1187  const amrex::Real tol = 1.0e-3;
1188  const amrex::Real alpha = 0.5;
1189  const amrex::Real WSMIN = 0.1; // minimum wind speed
1190 };
1191 
1192 
1193 /**
1194  * Surface temperature with charnock roughness
1195  */
1197 {
1199  amrex::Real Tflux,
1200  amrex::Real Qvflux,
1201  amrex::Real cnk_a,
1202  bool cnk_visc,
1203  bool cons_qflux)
1204  {
1205  mdata.zref = zref;
1206  mdata.surf_temp_flux = Tflux;
1207  mdata.surf_moist_flux = Qvflux;
1208  mdata.Cnk_a = cnk_a;
1209  mdata.visc = cnk_visc;
1210  spec_qflux = cons_qflux;
1211  }
1212 
1213  AMREX_GPU_DEVICE
1214  AMREX_FORCE_INLINE
1215  void
1216  iterate_flux (const int& i,
1217  const int& j,
1218  const int& k,
1219  const int& max_iters,
1220  const amrex::Array4<amrex::Real>& z0_arr,
1221  const amrex::Array4<const amrex::Real>& umm_arr,
1222  const amrex::Array4<const amrex::Real>& tm_arr,
1223  const amrex::Array4<const amrex::Real>& tvm_arr,
1224  const amrex::Array4<const amrex::Real>& qvm_arr,
1225  const amrex::Array4<amrex::Real>& u_star_arr,
1226  const amrex::Array4<amrex::Real>& w_star_arr,
1227  const amrex::Array4<amrex::Real>& t_star_arr,
1228  const amrex::Array4<amrex::Real>& q_star_arr,
1229  const amrex::Array4<amrex::Real>& t_surf_arr,
1230  const amrex::Array4<amrex::Real>& q_surf_arr,
1231  const amrex::Array4<amrex::Real>& olen_arr,
1232  const amrex::Array4<amrex::Real>& pblh_arr,
1233  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
1234  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
1235  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
1236  {
1237  amrex::Real Rib = 0.0;
1238  amrex::Real zeta = 0.0;
1239  amrex::Real zeta_old = 0.0;
1240  amrex::Real psi_m = 0.0;
1241  amrex::Real psi_h = 0.0;
1242  amrex::Real num = 0.0;
1243  amrex::Real den = 0.0;
1244  amrex::Real z0 = z0_arr(i,j,k);
1245  amrex::Real z0_old = z0;
1246  amrex::Real ustar = 0.0;
1247  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
1248  amrex::Real C = std::log(mdata.zref / z0);
1249 
1250  // First iteration we assume neutral (L -> inf)
1251  if (u_star_arr(i,j,k) == 1.E34) { olen_arr(i,j,k) = 1.0e3; }
1252  zeta = mdata.zref / olen_arr(i,j,k);
1253 
1254  // Water vapor in atmos and surface
1255  amrex::Real qv_s, qv_a;
1256  if (q_surf_arr(i,j,k) > 0.) {
1257  qv_s = q_surf_arr(i,j,k);
1258  } else {
1259  // First iteration and no qv_surf was specified
1260  // Use mean since there will be no flux
1261  qv_s = qvm_arr(i,j,k);
1262  }
1263  qv_a = qvm_arr(i,j,k);
1264 
1265  // update w* and Umagmean from Beljaars (1995)
1266  if (w_star_arr) {
1267  // NOTE: Thv flux is lagged, similar to WRF
1268  psi_m = sfuns.calc_psi_m2(zeta);
1269  psi_h = sfuns.calc_psi_h2(zeta);
1270  ustar = mdata.kappa * umm / (C - psi_m);
1271  amrex::Real tstar = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) / (C - psi_h);
1273  -ustar * mdata.kappa * (qv_a - qv_s) / (C - psi_h);
1274  amrex::Real tflux = -ustar*tstar*(1 + 0.61*qvm_arr(i,j,k)) + 0.61*tm_arr(i,j,k)*qflux;
1275  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
1276  amrex::Real wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
1277  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
1278  umm = std::max(umm, WSMIN);
1279  }
1280 
1281  // Bulk Richardson number w/ moisture
1282  amrex::Real thv_s = t_surf_arr(i,j,k) * (1.0 + 0.61*qv_s);
1283  amrex::Real thv_a = tm_arr(i,j,k) * (1.0 + 0.61*qv_a);
1284  Rib = ( (mdata.gravity * mdata.zref) / tm_arr(i,j,k) ) *
1285  ( (thv_a - thv_s) / (umm * umm) );
1286  Rib = std::min(std::max(Rib,-4.0),4.0);
1287 
1288  // Fixed point iteration on zeta
1289  int iter = 0;
1290  do {
1291  // Transfer curr to old
1292  zeta_old = zeta;
1293 
1294  // Stability functions
1295  psi_m = sfuns.calc_psi_m2(zeta_old);
1296  psi_h = sfuns.calc_psi_h2(zeta_old);
1297 
1298  // Fixed point iteration on roughness
1299  int iter_z = 0;
1300  do {
1301  // Transfer curr to old
1302  z0_old = z0;
1303 
1304  // Update
1305  C = std::log(mdata.zref / z0_old);
1306  ustar = mdata.kappa * umm / (C - psi_m);
1307  if (mdata.Cnk_a > 0) {
1308  z0 = (mdata.Cnk_a / mdata.gravity) * ustar * ustar;
1309  if (mdata.visc) {
1310  z0 += air_viscosity(tm_arr(i,j,k)) / std::max(ustar, 0.05);
1311  }
1312  } else {
1313  z0 = COARE3_roughness(mdata.zref, umm, ustar);
1314  }
1315 
1316  ++iter_z;
1317  } while ( (std::abs(z0 - z0_old) > tol_z) && (iter_z <= max_iters) );
1318  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter_z < max_iters,
1319  "Maximum number of MOST roughness iterations reached.");
1320  C = std::log(mdata.zref / z0);
1321 
1322  // Limiting
1323  num = std::max(C - psi_m, 1.0);
1324  den = std::max(C - psi_h, 1.0);
1325 
1326  // Update with under relaxation
1327  zeta = (1.0 - alpha) * zeta_old + alpha * Rib * num * num / den;
1328 
1329  ++iter;
1330  } while ( (std::abs(zeta - zeta_old) > tol) && (iter <= max_iters) );
1331  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
1332  "Maximum number of MOST iterations reached.");
1333 
1334  // Populate the stored MOST arrays
1335  z0_arr(i,j,k) = z0;
1336  olen_arr(i,j,k) = mdata.zref / zeta;
1337  u_star_arr(i,j,k) = mdata.kappa * umm / (C - psi_m);
1338  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) / (C - psi_h);
1339  if (spec_qflux) {
1340  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (C - psi_h) /
1341  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,k);
1342  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
1343  } else {
1344  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,k) - q_surf_arr(i,j,k)) / (C - psi_h);
1345  }
1346  }
1347 
1348 private:
1352  const amrex::Real tol = 1.0e-3;
1353  const amrex::Real tol_z = 1.0e-10;
1354  const amrex::Real alpha = 0.5;
1355  const amrex::Real WSMIN = 0.1; // minimum wind speed
1356 };
1357 
1358 
1359 /**
1360  * Surface temperature with modified charnock roughness
1361  */
1363 {
1365  amrex::Real Tflux,
1366  amrex::Real Qvflux,
1367  amrex::Real depth,
1368  bool cons_qflux)
1369  {
1370  mdata.zref = zref;
1371  mdata.surf_temp_flux = Tflux;
1372  mdata.surf_moist_flux = Qvflux;
1373  mdata.Cnk_d = depth;
1374  mdata.Cnk_b = mdata.Cnk_b1 * std::log(mdata.Cnk_b2 / mdata.Cnk_d);
1375  spec_qflux = cons_qflux;
1376  }
1377 
1378  AMREX_GPU_DEVICE
1379  AMREX_FORCE_INLINE
1380  void
1381  iterate_flux (const int& i,
1382  const int& j,
1383  const int& k,
1384  const int& max_iters,
1385  const amrex::Array4<amrex::Real>& z0_arr,
1386  const amrex::Array4<const amrex::Real>& umm_arr,
1387  const amrex::Array4<const amrex::Real>& tm_arr,
1388  const amrex::Array4<const amrex::Real>& tvm_arr,
1389  const amrex::Array4<const amrex::Real>& qvm_arr,
1390  const amrex::Array4<amrex::Real>& u_star_arr,
1391  const amrex::Array4<amrex::Real>& w_star_arr,
1392  const amrex::Array4<amrex::Real>& t_star_arr,
1393  const amrex::Array4<amrex::Real>& q_star_arr,
1394  const amrex::Array4<amrex::Real>& t_surf_arr,
1395  const amrex::Array4<amrex::Real>& q_surf_arr,
1396  const amrex::Array4<amrex::Real>& olen_arr,
1397  const amrex::Array4<amrex::Real>& pblh_arr,
1398  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
1399  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
1400  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
1401  {
1402  amrex::Real Rib = 0.0;
1403  amrex::Real zeta = 0.0;
1404  amrex::Real zeta_old = 0.0;
1405  amrex::Real psi_m = 0.0;
1406  amrex::Real psi_h = 0.0;
1407  amrex::Real num = 0.0;
1408  amrex::Real den = 0.0;
1409  amrex::Real z0 = z0_arr(i,j,k);
1410  amrex::Real z0_old = z0;
1411  amrex::Real ustar = 0.0;
1412  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
1413  amrex::Real C = std::log(mdata.zref / z0);
1414 
1415  // First iteration we assume neutral (L -> inf)
1416  if (u_star_arr(i,j,k) == 1.E34) { olen_arr(i,j,k) = 1.0e3; }
1417  zeta = mdata.zref / olen_arr(i,j,k);
1418 
1419  // Water vapor in atmos and surface
1420  amrex::Real qv_s, qv_a;
1421  if (q_surf_arr(i,j,k) > 0.) {
1422  qv_s = q_surf_arr(i,j,k);
1423  } else {
1424  // First iteration and no qv_surf was specified
1425  // Use mean since there will be no flux
1426  qv_s = qvm_arr(i,j,k);
1427  }
1428  qv_a = qvm_arr(i,j,k);
1429 
1430  // update w* and Umagmean from Beljaars (1995)
1431  if (w_star_arr) {
1432  // NOTE: Thv flux is lagged, similar to WRF
1433  psi_m = sfuns.calc_psi_m2(zeta);
1434  psi_h = sfuns.calc_psi_h2(zeta);
1435  ustar = mdata.kappa * umm / (C - psi_m);
1436  amrex::Real tstar = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) / (C - psi_h);
1438  -ustar * mdata.kappa * (qv_a - qv_s) / (C - psi_h);
1439  amrex::Real tflux = -ustar*tstar*(1 + 0.61*qvm_arr(i,j,k)) + 0.61*tm_arr(i,j,k)*qflux;
1440  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
1441  amrex::Real wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
1442  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
1443  umm = std::max(umm, WSMIN);
1444  }
1445 
1446  // Bulk Richardson number w/ moisture
1447  amrex::Real thv_s = t_surf_arr(i,j,k) * (1.0 + 0.61*qv_s);
1448  amrex::Real thv_a = tm_arr(i,j,k) * (1.0 + 0.61*qv_a);
1449  Rib = ( (mdata.gravity * mdata.zref) / tm_arr(i,j,k) ) *
1450  ( (thv_a - thv_s) / (umm * umm) );
1451  Rib = std::min(std::max(Rib,-4.0),4.0);
1452 
1453  // Fixed point iteration on zeta
1454  int iter = 0;
1455  do {
1456  // Transfer curr to old
1457  zeta_old = zeta;
1458 
1459  // Stability functions
1460  psi_m = sfuns.calc_psi_m2(zeta_old);
1461  psi_h = sfuns.calc_psi_h2(zeta_old);
1462 
1463  // Fixed point iteration on roughness
1464  int iter_z = 0;
1465  do {
1466  // Transfer curr to old
1467  z0_old = z0;
1468 
1469  // Update
1470  C = std::log(mdata.zref / z0_old);
1471  ustar = mdata.kappa * umm / (C - psi_m);
1472  z0 = std::exp( (2.7*ustar - 1.8/mdata.Cnk_b) / (ustar + 0.17/mdata.Cnk_b) );
1473 
1474  ++iter_z;
1475  } while ( (std::abs(z0 - z0_old) > tol_z) && (iter_z <= max_iters) );
1476  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter_z < max_iters,
1477  "Maximum number of MOST roughness iterations reached.");
1478  C = std::log(mdata.zref / z0);
1479 
1480  // Limiting
1481  num = std::max(C - psi_m, 1.0);
1482  den = std::max(C - psi_h, 1.0);
1483 
1484  // Update with under relaxation
1485  zeta = (1.0 - alpha) * zeta_old + alpha * Rib * num * num / den;
1486 
1487  ++iter;
1488  } while ( (std::abs(zeta - zeta_old) > tol) && (iter <= max_iters) );
1489  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
1490  "Maximum number of MOST iterations reached.");
1491 
1492  // Populate the stored MOST arrays
1493  z0_arr(i,j,k) = z0;
1494  olen_arr(i,j,k) = mdata.zref / zeta;
1495  u_star_arr(i,j,k) = mdata.kappa * umm / (C - psi_m);
1496  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) / (C - psi_h);
1497  if (spec_qflux) {
1498  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (C - psi_h) /
1499  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,k);
1500  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
1501  } else {
1502  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,k) - q_surf_arr(i,j,k)) / (C - psi_h);
1503  }
1504  }
1505 
1506 private:
1510  const amrex::Real tol = 1.0e-3;
1511  const amrex::Real tol_z = 1.0e-10;
1512  const amrex::Real alpha = 0.5;
1513  const amrex::Real WSMIN = 0.1; // minimum wind speed
1514 };
1515 
1516 
1517 /**
1518  * Surface temperature with donelan roughness
1519  */
1521 {
1523  amrex::Real Tflux,
1524  amrex::Real Qvflux,
1525  bool cons_qflux)
1526  {
1527  mdata.zref = zref;
1528  mdata.surf_temp_flux = Tflux;
1529  mdata.surf_moist_flux = Qvflux;
1530  spec_qflux = cons_qflux;
1531  }
1532 
1533  AMREX_GPU_DEVICE
1534  AMREX_FORCE_INLINE
1535  void
1536  iterate_flux (const int& i,
1537  const int& j,
1538  const int& k,
1539  const int& max_iters,
1540  const amrex::Array4<amrex::Real>& z0_arr,
1541  const amrex::Array4<const amrex::Real>& umm_arr,
1542  const amrex::Array4<const amrex::Real>& tm_arr,
1543  const amrex::Array4<const amrex::Real>& tvm_arr,
1544  const amrex::Array4<const amrex::Real>& qvm_arr,
1545  const amrex::Array4<amrex::Real>& u_star_arr,
1546  const amrex::Array4<amrex::Real>& w_star_arr,
1547  const amrex::Array4<amrex::Real>& t_star_arr,
1548  const amrex::Array4<amrex::Real>& q_star_arr,
1549  const amrex::Array4<amrex::Real>& t_surf_arr,
1550  const amrex::Array4<amrex::Real>& q_surf_arr,
1551  const amrex::Array4<amrex::Real>& olen_arr,
1552  const amrex::Array4<amrex::Real>& pblh_arr,
1553  const amrex::Array4<amrex::Real>& /*Hwave_arr*/,
1554  const amrex::Array4<amrex::Real>& /*Lwave_arr*/,
1555  const amrex::Array4<amrex::Real>& /*eta_arr*/) const
1556  {
1557  int iter = 0;
1558  amrex::Real ustar = 0.0;
1559  amrex::Real wstar = 0.0;
1560  amrex::Real z0 = 0.0;
1561  amrex::Real tflux = 0.0;
1562  amrex::Real qflux = 0.0;
1563  amrex::Real zeta = 0.0;
1564  amrex::Real psi_m = 0.0;
1565  amrex::Real psi_h = 0.0;
1566  amrex::Real Olen = 0.0;
1567  amrex::Real Oleno = 0.0;
1568  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
1569  if (u_star_arr(i,j,k) == 1.E34) {
1570  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
1571  } else {
1572  Olen = olen_arr(i,j,k);
1573  Oleno = Olen;
1574  zeta = mdata.zref / Olen;
1575  psi_m = sfuns.calc_psi_m(zeta);
1576  psi_h = sfuns.calc_psi_h(zeta);
1577  }
1578  do {
1579  ustar = u_star_arr(i,j,k);
1580  z0 = Donelan_roughness(ustar);
1581  tflux = -(tm_arr(i,j,k) - t_surf_arr(i,j,k)) * ustar * mdata.kappa /
1582  (std::log(mdata.zref / z0) - psi_h); // <w'T'>
1583  tflux *= (1 + 0.61*qvm_arr(i,j,k));
1584  qflux = (spec_qflux) ? mdata.surf_moist_flux :
1585  -(qvm_arr(i,j,k) - q_surf_arr(i,j,k)) * ustar * mdata.kappa /
1586  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h); // <w'Qv'>
1587  tflux += 0.61*tm_arr(i,j,k) * qflux; // ~= <w'Tv'>
1588  if (w_star_arr) {
1589  // update w* and Umagmean
1590  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
1591  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
1592  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
1593  umm = std::max(umm, WSMIN);
1594  }
1595  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
1596  if ( (((Olen >= 0.0) && (Oleno <= 0.0)) ||
1597  ((Olen <= 0.0) && (Oleno >= 0.0))) &&
1598  std::fabs(Olen) + std::fabs(Oleno) < 1.) {
1599  Olen = 0.5 * (Olen + Oleno);
1600  }
1601  Oleno = Olen;
1602  zeta = mdata.zref / Olen;
1603  psi_m = sfuns.calc_psi_m(zeta);
1604  psi_h = sfuns.calc_psi_h(zeta);
1605  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
1606  ++iter;
1607  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
1608  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
1609  "Maximum number of MOST iterations reached.");
1610 
1611  // Populate the stored MOST arrays
1612  z0_arr(i,j,k) = z0;
1613  olen_arr(i,j,k) = Olen;
1614  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) /
1615  (std::log(mdata.zref / z0) - psi_h);
1616  if (spec_qflux) {
1617  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h) /
1618  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,k);
1619  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
1620  } else {
1621  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,k) - q_surf_arr(i,j,k)) /
1622  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h);
1623  }
1624  }
1625 
1626 private:
1630  const amrex::Real tol = 1.0e-5;
1631  const amrex::Real WSMIN = 0.1; // minimum wind speed
1632 };
1633 
1634 
1635 /**
1636  * Surface temperature with wave-coupled roughness
1637  */
1639 {
1641  amrex::Real Tflux,
1642  amrex::Real Qvflux,
1643  bool cons_qflux)
1644  {
1645  mdata.zref = zref;
1646  mdata.surf_temp_flux = Tflux;
1647  mdata.surf_moist_flux = Qvflux;
1648  spec_qflux = cons_qflux;
1649  }
1650 
1651  AMREX_GPU_DEVICE
1652  AMREX_FORCE_INLINE
1653  void
1654  iterate_flux (const int& i,
1655  const int& j,
1656  const int& k,
1657  const int& max_iters,
1658  const amrex::Array4<amrex::Real>& z0_arr,
1659  const amrex::Array4<const amrex::Real>& umm_arr,
1660  const amrex::Array4<const amrex::Real>& tm_arr,
1661  const amrex::Array4<const amrex::Real>& tvm_arr,
1662  const amrex::Array4<const amrex::Real>& qvm_arr,
1663  const amrex::Array4<amrex::Real>& u_star_arr,
1664  const amrex::Array4<amrex::Real>& w_star_arr,
1665  const amrex::Array4<amrex::Real>& t_star_arr,
1666  const amrex::Array4<amrex::Real>& q_star_arr,
1667  const amrex::Array4<amrex::Real>& t_surf_arr,
1668  const amrex::Array4<amrex::Real>& q_surf_arr,
1669  const amrex::Array4<amrex::Real>& olen_arr,
1670  const amrex::Array4<amrex::Real>& pblh_arr,
1671  const amrex::Array4<amrex::Real>& Hwave_arr,
1672  const amrex::Array4<amrex::Real>& Lwave_arr,
1673  const amrex::Array4<amrex::Real>& eta_arr) const
1674  {
1675  int iter = 0;
1676  amrex::Real ustar = 0.0;
1677  amrex::Real wstar = 0.0;
1678  amrex::Real z0 = 0.0;
1679  amrex::Real tflux = 0.0;
1680  amrex::Real qflux = 0.0;
1681  amrex::Real zeta = 0.0;
1682  amrex::Real psi_m = 0.0;
1683  amrex::Real psi_h = 0.0;
1684  amrex::Real Olen = 0.0;
1685  amrex::Real Oleno = 0.0;
1686  int ie, je;
1687  ie = i < lbound(eta_arr).x ? lbound(eta_arr).x : i;
1688  je = j < lbound(eta_arr).y ? lbound(eta_arr).y : j;
1689  ie = ie > ubound(eta_arr).x ? ubound(eta_arr).x : ie;
1690  je = je > ubound(eta_arr).y ? ubound(eta_arr).y : je;
1691  amrex::Real umm = std::max(umm_arr(i,j,k), WSMIN);
1692  if (u_star_arr(i,j,k) == 1.E34) {
1693  u_star_arr(i,j,k) = mdata.kappa * umm / std::log(mdata.zref / z0_arr(i,j,k));
1694  } else {
1695  Olen = olen_arr(i,j,k);
1696  Oleno = Olen;
1697  zeta = mdata.zref / Olen;
1698  psi_m = sfuns.calc_psi_m(zeta);
1699  psi_h = sfuns.calc_psi_h(zeta);
1700  }
1701  do {
1702  ustar = u_star_arr(i,j,k);
1703  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 )
1704  + 0.11 * eta_arr(ie,je,k,EddyDiff::Mom_v) / ustar, z0_eps), z0_max );
1705  tflux = -(tm_arr(i,j,k) - t_surf_arr(i,j,k)) * ustar * mdata.kappa /
1706  (std::log(mdata.zref / z0) - psi_h); // <w'T'>
1707  tflux *= (1 + 0.61*qvm_arr(i,j,k));
1708  qflux = (spec_qflux) ? mdata.surf_moist_flux :
1709  -(qvm_arr(i,j,k) - q_surf_arr(i,j,k)) * ustar * mdata.kappa /
1710  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h); // <w'Qv'>
1711  tflux += 0.61*tm_arr(i,j,k) * qflux; // ~= <w'Tv'>
1712  if (w_star_arr) {
1713  // update w* and Umagmean
1714  w_star_arr(i,j,k) = calc_wstar(tflux, pblh_arr(i,j,k), tvm_arr(i,j,k));
1715  wstar = mdata.Bjr_beta * w_star_arr(i,j,k);
1716  umm = std::sqrt(umm_arr(i,j,k)*umm_arr(i,j,k) + wstar*wstar);
1717  umm = std::max(umm, WSMIN);
1718  }
1719  Olen = -ustar * ustar * ustar * tvm_arr(i,j,k) / (mdata.kappa * mdata.gravity * tflux);
1720  if ( (((Olen >= 0.0) && (Oleno <= 0.0)) ||
1721  ((Olen <= 0.0) && (Oleno >= 0.0))) &&
1722  std::fabs(Olen) + std::fabs(Oleno) < 1.) {
1723  Olen = 0.5 * (Olen + Oleno);
1724  }
1725  Oleno = Olen;
1726  zeta = mdata.zref / Olen;
1727  psi_m = sfuns.calc_psi_m(zeta);
1728  psi_h = sfuns.calc_psi_h(zeta);
1729  u_star_arr(i,j,k) = mdata.kappa * umm / (std::log(mdata.zref / z0) - psi_m);
1730  ++iter;
1731  } while ((std::abs(u_star_arr(i,j,k) - ustar) > tol) && iter <= max_iters);
1732  AMREX_ALWAYS_ASSERT_WITH_MESSAGE(iter < max_iters,
1733  "Maximum number of MOST iterations reached.");
1734 
1735  // Populate the stored MOST arrays
1736  z0_arr(i,j,k) = z0;
1737  olen_arr(i,j,k) = Olen;
1738  t_star_arr(i,j,k) = mdata.kappa * (tm_arr(i,j,k) - t_surf_arr(i,j,k)) /
1739  (std::log(mdata.zref / z0) - psi_h);
1740  if (spec_qflux) {
1741  q_surf_arr(i,j,k) = mdata.surf_moist_flux * (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h) /
1742  (u_star_arr(i,j,k) * mdata.kappa) + qvm_arr(i,j,k);
1743  q_star_arr(i,j,k) = -mdata.surf_moist_flux / u_star_arr(i,j,k);
1744  } else {
1745  q_star_arr(i,j,k) = mdata.kappa * (qvm_arr(i,j,k) - q_surf_arr(i,j,k)) /
1746  (std::log(mdata.zref / z0_arr(i,j,k)) - psi_h);
1747  }
1748  }
1749 
1750 private:
1754  const amrex::Real tol = 1.0e-5;
1755  const amrex::Real eps = 1e-15;
1756  const amrex::Real z0_eps = 1.0e-6;
1757  const amrex::Real z0_max = 0.1;
1758  const amrex::Real WSMIN = 0.1; // minimum wind speed
1759 };
1760 
1761 
1762 /**
1763  * Moeng flux formulation
1764  */
1766 {
1768 
1769  AMREX_GPU_DEVICE
1770  AMREX_FORCE_INLINE
1771  amrex::Real
1772  compute_q_flux (const int& i,
1773  const int& j,
1774  const int& k,
1775  const amrex::Array4<const amrex::Real>& cons_arr,
1776  const amrex::Array4<const amrex::Real>& velx_arr,
1777  const amrex::Array4<const amrex::Real>& vely_arr,
1778  const amrex::Array4<const amrex::Real>& umm_arr,
1779  const amrex::Array4<const amrex::Real>& qvm_arr,
1780  const amrex::Array4<const amrex::Real>& u_star_arr,
1781  const amrex::Array4<const amrex::Real>& q_star_arr,
1782  const amrex::Array4<const amrex::Real>& q_surf_arr) const
1783  {
1784  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
1785  amrex::Real qv = cons_arr(i,j,k,RhoQ1_comp) / rho;
1786  amrex::Real velx = 0.5 * ( velx_arr(i,j,k) + velx_arr(i+1,j ,k) );
1787  amrex::Real vely = 0.5 * ( vely_arr(i,j,k) + vely_arr(i ,j+1,k) );
1788 
1789  amrex::Real qv_mean = qvm_arr(i,j,k);
1790  amrex::Real ustar = u_star_arr(i,j,k);
1791  amrex::Real qstar = q_star_arr(i,j,k);
1792  amrex::Real qv_surf = q_surf_arr(i,j,k);
1793  amrex::Real wsp_mean = umm_arr(i,j,k);
1794  wsp_mean = std::max(wsp_mean, WSMIN);
1795 
1796  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1797  amrex::Real num1 = wsp * (qv_mean-qv_surf);
1798  amrex::Real num2 = wsp_mean * (qv-qv_mean);
1799 
1800  // NOTE: this is rho*<Qv'w'> = -K dQvdz
1801  amrex::Real moflux = (std::abs(qstar) > eps) ?
1802  -rho*qstar*ustar*(num1+num2)/((qv_mean-qv_surf)*wsp_mean) : 0.0;
1803 
1804  return moflux;
1805  }
1806 
1807  AMREX_GPU_DEVICE
1808  AMREX_FORCE_INLINE
1809  amrex::Real
1810  compute_t_flux (const int& i,
1811  const int& j,
1812  const int& k,
1813  const amrex::Array4<const amrex::Real>& cons_arr,
1814  const amrex::Array4<const amrex::Real>& velx_arr,
1815  const amrex::Array4<const amrex::Real>& vely_arr,
1816  const amrex::Array4<const amrex::Real>& umm_arr,
1817  const amrex::Array4<const amrex::Real>& tm_arr,
1818  const amrex::Array4<const amrex::Real>& u_star_arr,
1819  const amrex::Array4<const amrex::Real>& t_star_arr,
1820  const amrex::Array4<const amrex::Real>& t_surf_arr) const
1821  {
1822  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
1823  amrex::Real theta = cons_arr(i,j,k,RhoTheta_comp) / rho;
1824  amrex::Real velx = 0.5 * ( velx_arr(i,j,k) + velx_arr(i+1,j ,k) );
1825  amrex::Real vely = 0.5 * ( vely_arr(i,j,k) + vely_arr(i ,j+1,k) );
1826 
1827  amrex::Real theta_mean = tm_arr(i,j,k);
1828  amrex::Real ustar = u_star_arr(i,j,k);
1829  amrex::Real tstar = t_star_arr(i,j,k);
1830  amrex::Real theta_surf = t_surf_arr(i,j,k);
1831  amrex::Real wsp_mean = umm_arr(i,j,k);
1832  wsp_mean = std::max(wsp_mean, WSMIN);
1833 
1834  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1835  amrex::Real num1 = wsp * (theta_mean-theta_surf);
1836  amrex::Real num2 = wsp_mean * (theta-theta_mean);
1837 
1838  // NOTE: this is rho*<T'w'> = -K dTdz
1839  amrex::Real moflux = (std::abs(tstar) > eps) ?
1840  -rho*tstar*ustar*(num1+num2)/((theta_mean-theta_surf)*wsp_mean) : 0.0;
1841 
1842  return moflux;
1843  }
1844 
1845  AMREX_GPU_DEVICE
1846  AMREX_FORCE_INLINE
1847  amrex::Real
1848  compute_u_flux (const int& i,
1849  const int& j,
1850  const int& k,
1851  const amrex::Array4<const amrex::Real>& cons_arr,
1852  const amrex::Array4<const amrex::Real>& velx_arr,
1853  const amrex::Array4<const amrex::Real>& vely_arr,
1854  const amrex::Array4<const amrex::Real>& umm_arr,
1855  const amrex::Array4<const amrex::Real>& um_arr,
1856  const amrex::Array4<const amrex::Real>& u_star_arr) const
1857  {
1858  amrex::Real velx = velx_arr(i,j,k);
1859  amrex::Real vely = 0.25 * ( vely_arr(i ,j,k) + vely_arr(i ,j+1,k)
1860  + vely_arr(i-1,j,k) + vely_arr(i-1,j+1,k) );
1861  amrex::Real rho = 0.5 * ( cons_arr(i-1,j,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
1862 
1863  amrex::Real umean = um_arr(i,j,k);
1864  amrex::Real ustar = 0.5 * ( u_star_arr(i-1,j,k) + u_star_arr(i,j,k) );
1865  amrex::Real wsp_mean = 0.5 * ( umm_arr(i-1,j,k) + umm_arr(i,j,k) );
1866  wsp_mean = std::max(wsp_mean, WSMIN);
1867 
1868  // Note: The surface mean shear stress is decomposed into tau_xz by
1869  // multiplying the modeled shear stress (rho*ustar^2) with
1870  // a factor of umean/wsp_mean for directionality; this factor
1871  // modifies the denominator from what is in Moeng 1984.
1872  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1873  amrex::Real num1 = wsp * umean;
1874  amrex::Real num2 = wsp_mean * (velx-umean);
1875 
1876  // NOTE: this is rho*<u'w'> = -K dudz
1877  amrex::Real stressx = -rho*ustar*ustar * (num1+num2)/(wsp_mean*wsp_mean);
1878 
1879  return stressx;
1880  }
1881 
1882  AMREX_GPU_DEVICE
1883  AMREX_FORCE_INLINE
1884  amrex::Real
1885  compute_v_flux (const int& i,
1886  const int& j,
1887  const int& k,
1888  const amrex::Array4<const amrex::Real>& cons_arr,
1889  const amrex::Array4<const amrex::Real>& velx_arr,
1890  const amrex::Array4<const amrex::Real>& vely_arr,
1891  const amrex::Array4<const amrex::Real>& umm_arr,
1892  const amrex::Array4<const amrex::Real>& vm_arr,
1893  const amrex::Array4<const amrex::Real>& u_star_arr) const
1894  {
1895  amrex::Real velx = 0.25 * ( velx_arr(i,j ,k) + velx_arr(i+1,j ,k)
1896  + velx_arr(i,j-1,k) + velx_arr(i+1,j-1,k) );
1897  amrex::Real vely = vely_arr(i,j,k);
1898  amrex::Real rho = 0.5 * ( cons_arr(i,j-1,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
1899 
1900  amrex::Real vmean = vm_arr(i,j,k);
1901  amrex::Real ustar = 0.5 * ( u_star_arr(i,j-1,k) + u_star_arr(i,j,k) );
1902  amrex::Real wsp_mean = 0.5 * ( umm_arr(i,j-1,k) + umm_arr(i,j,k) );
1903  wsp_mean = std::max(wsp_mean, WSMIN);
1904 
1905  // Note: The surface mean shear stress is decomposed into tau_yz by
1906  // multiplying the modeled shear stress (rho*ustar^2) with
1907  // a factor of vmean/wsp_mean for directionality; this factor
1908  // modifies the denominator from what is in Moeng 1984.
1909  amrex::Real wsp = sqrt(velx*velx+vely*vely);
1910  amrex::Real num1 = wsp * vmean;
1911  amrex::Real num2 = wsp_mean * (vely-vmean);
1912 
1913  // NOTE: this is rho*<v'w'> = -K dvdz
1914  amrex::Real stressy = -rho*ustar*ustar * (num1+num2)/(wsp_mean*wsp_mean);
1915 
1916  return stressy;
1917  }
1918 
1919 private:
1920  const amrex::Real eps = 1e-15;
1921  const amrex::Real WSMIN = 0.1; // minimum wind speed
1922 };
1923 
1924 
1925 /**
1926  * Donelan flux formulation
1927  */
1929 {
1931 
1932  AMREX_GPU_DEVICE
1933  AMREX_FORCE_INLINE
1934  amrex::Real
1935  compute_q_flux (const int& /*i*/,
1936  const int& /*j*/,
1937  const int& /*k*/,
1938  const amrex::Array4<const amrex::Real>& /*cons_arr*/,
1939  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1940  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1941  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
1942  const amrex::Array4<const amrex::Real>& /*qm_arr*/,
1943  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
1944  const amrex::Array4<const amrex::Real>& /*q_star_arr*/,
1945  const amrex::Array4<const amrex::Real>& /*q_surf_arr*/) const
1946  {
1947  // NOTE: this is rho*<q'w'> = -K dqdz
1948  amrex::Real moflux = 0.0;
1949 
1950  return moflux;
1951  }
1952 
1953  AMREX_GPU_DEVICE
1954  AMREX_FORCE_INLINE
1955  amrex::Real
1956  compute_t_flux (const int& i,
1957  const int& j,
1958  const int& k,
1959  const amrex::Array4<const amrex::Real>& cons_arr,
1960  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
1961  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
1962  const amrex::Array4<const amrex::Real>& umm_arr,
1963  const amrex::Array4<const amrex::Real>& tm_arr,
1964  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
1965  const amrex::Array4<const amrex::Real>& /*t_star_arr*/,
1966  const amrex::Array4<const amrex::Real>& t_surf_arr) const
1967  {
1968  amrex::Real Ch = 0.0012;
1969  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
1970  amrex::Real theta_surf = t_surf_arr(i,j,k);
1971  amrex::Real theta_mean = tm_arr(i,j,k);
1972  amrex::Real wsp_mean = umm_arr(i,j,k);
1973 
1974  // NOTE: this is rho*<T'w'> = -K dTdz
1975  amrex::Real moflux = -rho * Ch * wsp_mean * (theta_mean - theta_surf);
1976 
1977  return moflux;
1978  }
1979 
1980  AMREX_GPU_DEVICE
1981  AMREX_FORCE_INLINE
1982  amrex::Real
1983  compute_u_flux (const int& i,
1984  const int& j,
1985  const int& k,
1986  const amrex::Array4<const amrex::Real>& cons_arr,
1987  const amrex::Array4<const amrex::Real>& velx_arr,
1988  const amrex::Array4<const amrex::Real>& vely_arr,
1989  const amrex::Array4<const amrex::Real>& umm_arr,
1990  const amrex::Array4<const amrex::Real>& /*um_arr*/,
1991  const amrex::Array4<const amrex::Real>& /*u_star_arr*/) const
1992  {
1993  amrex::Real velx = velx_arr(i,j,k);
1994  amrex::Real vely = 0.25 * ( vely_arr(i ,j,k) + vely_arr(i ,j+1,k)
1995  + vely_arr(i-1,j,k) + vely_arr(i-1,j+1,k) );
1996  amrex::Real rho = 0.5 * ( cons_arr(i-1,j,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
1997 
1998  amrex::Real Cd = 0.001;
1999  const amrex::Real c = 7e-5;
2000  amrex::Real wsp = sqrt(velx*velx+vely*vely);
2001  amrex::Real wsp_mean = 0.5 * ( umm_arr(i-1,j,k) + umm_arr(i,j,k) );
2002  if (wsp_mean <= 5.0) {
2003  Cd = 0.001;
2004  } else if (wsp_mean < 25.0 && wsp_mean > 5.0) {
2005  Cd = 0.001 + c * (wsp_mean - 5.0);
2006  } else {
2007  Cd = 0.0024;
2008  }
2009 
2010  // NOTE: this is rho*<u'w'> = -K dudz
2011  amrex::Real stressx = -rho * Cd * velx * wsp;
2012 
2013  return stressx;
2014  }
2015 
2016  AMREX_GPU_DEVICE
2017  AMREX_FORCE_INLINE
2018  amrex::Real
2019  compute_v_flux (const int& i,
2020  const int& j,
2021  const int& k,
2022  const amrex::Array4<const amrex::Real>& cons_arr,
2023  const amrex::Array4<const amrex::Real>& velx_arr,
2024  const amrex::Array4<const amrex::Real>& vely_arr,
2025  const amrex::Array4<const amrex::Real>& umm_arr,
2026  const amrex::Array4<const amrex::Real>& /*vm_arr*/,
2027  const amrex::Array4<const amrex::Real>& /*u_star_arr*/) const
2028  {
2029  amrex::Real velx = 0.25 * ( velx_arr(i,j ,k) + velx_arr(i+1,j ,k)
2030  + velx_arr(i,j-1,k) + velx_arr(i+1,j-1,k) );
2031  amrex::Real vely = vely_arr(i,j,k);
2032  amrex::Real rho = 0.5 * ( cons_arr(i,j-1,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
2033 
2034  amrex::Real Cd = 0.001;
2035  const amrex::Real c = 7e-5;
2036  amrex::Real wsp = sqrt(velx*velx+vely*vely);
2037  amrex::Real wsp_mean = 0.5 * ( umm_arr(i,j-1,k) + umm_arr(i,j,k) );
2038  if (wsp_mean <= 5.0) {
2039  Cd = 0.001;
2040  } else if (wsp_mean < 25.0 && wsp_mean > 5.0) {
2041  Cd = 0.001 + c * (wsp_mean - 5.0);
2042  } else {
2043  Cd = 0.0024;
2044  }
2045 
2046  // NOTE: this is rho*<v'w'> = -K dvdz
2047  amrex::Real stressy = -rho * Cd * vely * wsp;
2048 
2049  return stressy;
2050  }
2051 
2052 private:
2053 
2054 };
2055 
2056 
2057 /**
2058  * Custom flux formulation
2059  */
2061 {
2062  custom_flux (bool specified_rho_surf)
2063  : fluxes_include_rho(specified_rho_surf)
2064  {}
2065 
2066  AMREX_GPU_DEVICE
2067  AMREX_FORCE_INLINE
2068  amrex::Real
2069  compute_q_flux (const int& i,
2070  const int& j,
2071  const int& k,
2072  const amrex::Array4<const amrex::Real>& cons_arr,
2073  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
2074  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
2075  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
2076  const amrex::Array4<const amrex::Real>& /*qm_arr*/,
2077  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
2078  const amrex::Array4<const amrex::Real>& q_star_arr,
2079  const amrex::Array4<const amrex::Real>& /*q_surf_arr*/) const
2080  {
2081  amrex::Real rho = (fluxes_include_rho) ? 1.0 : cons_arr(i,j,k,Rho_comp);
2082  amrex::Real qstar = q_star_arr(i,j,k);
2083 
2084  // NOTE: this is rho*<q'w'> = -K dqdz
2085  amrex::Real moflux = (std::abs(qstar) > eps) ? rho * qstar : 0.0;
2086 
2087  return moflux;
2088  }
2089 
2090  AMREX_GPU_DEVICE
2091  AMREX_FORCE_INLINE
2092  amrex::Real
2093  compute_t_flux (const int& i,
2094  const int& j,
2095  const int& k,
2096  const amrex::Array4<const amrex::Real>& cons_arr,
2097  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
2098  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
2099  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
2100  const amrex::Array4<const amrex::Real>& /*tm_arr*/,
2101  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
2102  const amrex::Array4<const amrex::Real>& t_star_arr,
2103  const amrex::Array4<const amrex::Real>& /*t_surf_arr*/) const
2104  {
2105  amrex::Real rho = (fluxes_include_rho) ? 1.0 : cons_arr(i,j,k,Rho_comp);
2106  amrex::Real tstar = t_star_arr(i,j,k);
2107 
2108  // NOTE: this is rho*<T'w'> = -K dTdz
2109  amrex::Real moflux = (std::abs(tstar) > eps) ? rho * tstar : 0.0;
2110 
2111  return moflux;
2112  }
2113 
2114  AMREX_GPU_DEVICE
2115  AMREX_FORCE_INLINE
2116  amrex::Real
2117  compute_u_flux (const int& i,
2118  const int& j,
2119  const int& k,
2120  const amrex::Array4<const amrex::Real>& cons_arr,
2121  const amrex::Array4<const amrex::Real>& velx_arr,
2122  const amrex::Array4<const amrex::Real>& vely_arr,
2123  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
2124  const amrex::Array4<const amrex::Real>& /*um_arr*/,
2125  const amrex::Array4<const amrex::Real>& u_star_arr) const
2126  {
2127  amrex::Real velx = velx_arr(i,j,k);
2128  amrex::Real vely = 0.25 * ( vely_arr(i ,j,k) + vely_arr(i ,j+1,k)
2129  + vely_arr(i-1,j,k) + vely_arr(i-1,j+1,k) );
2130  amrex::Real rho = (fluxes_include_rho) ? 1.0 : 0.5 * ( cons_arr(i-1,j,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
2131 
2132  amrex::Real ustar = 0.5 * ( u_star_arr(i-1,j,k) + u_star_arr(i,j,k) );
2133  amrex::Real wsp = sqrt(velx*velx+vely*vely);
2134 
2135  // NOTE: this is rho*<u'w'> = -K dudz
2136  amrex::Real stressx = -rho * ustar * ustar * velx / wsp;
2137 
2138  return stressx;
2139  }
2140 
2141  AMREX_GPU_DEVICE
2142  AMREX_FORCE_INLINE
2143  amrex::Real
2144  compute_v_flux (const int& i,
2145  const int& j,
2146  const int& k,
2147  const amrex::Array4<const amrex::Real>& cons_arr,
2148  const amrex::Array4<const amrex::Real>& velx_arr,
2149  const amrex::Array4<const amrex::Real>& vely_arr,
2150  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
2151  const amrex::Array4<const amrex::Real>& /*vm_arr*/,
2152  const amrex::Array4<const amrex::Real>& u_star_arr) const
2153  {
2154  amrex::Real velx = 0.25 * ( velx_arr(i,j ,k) + velx_arr(i+1,j ,k)
2155  + velx_arr(i,j-1,k) + velx_arr(i+1,j-1,k) );
2156  amrex::Real vely = vely_arr(i,j,k);
2157  amrex::Real rho = (fluxes_include_rho) ? 1.0 : 0.5 * ( cons_arr(i,j-1,k,Rho_comp) + cons_arr(i,j,k,Rho_comp) );
2158 
2159  amrex::Real ustar = 0.5 * ( u_star_arr(i,j-1,k) + u_star_arr(i,j,k) );
2160  amrex::Real wsp = sqrt(velx*velx+vely*vely);
2161 
2162  // NOTE: this is rho*<v'w'> = -K dvdz
2163  amrex::Real stressy = -rho * ustar * ustar * vely / wsp;
2164 
2165  return stressy;
2166  }
2167 
2168 private:
2169  const amrex::Real eps = 1e-15;
2170  const bool fluxes_include_rho{false};
2171 
2172 };
2173 
2174 /**
2175  * Rotate flux formulation
2176  */
2178 {
2180 
2181  AMREX_GPU_DEVICE
2182  AMREX_FORCE_INLINE
2183  amrex::Real
2184  compute_q_flux (const int& i,
2185  const int& j,
2186  const int& k,
2187  const amrex::Array4<const amrex::Real>& cons_arr,
2188  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
2189  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
2190  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
2191  const amrex::Array4<const amrex::Real>& qvm_arr,
2192  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
2193  const amrex::Array4<const amrex::Real>& q_star_arr,
2194  const amrex::Array4<const amrex::Real>& q_surf_arr) const
2195  {
2196  // NOTE: this is the total stress
2197  amrex::Real qv_mean = qvm_arr(i,j,k);
2198  amrex::Real qv_surf = q_surf_arr(i,j,k);
2199  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
2200  amrex::Real qstar = q_star_arr(i,j,k);
2201  amrex::Real moflux = (std::abs(qstar) > eps) ? -rho * qstar * (qv_mean-qv_surf): 0.0;
2202 
2203  return moflux;
2204  }
2205 
2206  AMREX_GPU_DEVICE
2207  AMREX_FORCE_INLINE
2208  amrex::Real
2209  compute_t_flux (const int& i,
2210  const int& j,
2211  const int& k,
2212  const amrex::Array4<const amrex::Real>& cons_arr,
2213  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
2214  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
2215  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
2216  const amrex::Array4<const amrex::Real>& tm_arr,
2217  const amrex::Array4<const amrex::Real>& /*u_star_arr*/,
2218  const amrex::Array4<const amrex::Real>& t_star_arr,
2219  const amrex::Array4<const amrex::Real>& t_surf_arr) const
2220  {
2221  // NOTE: this is the total stress
2222  amrex::Real theta_mean = tm_arr(i,j,k);
2223  amrex::Real theta_surf = t_surf_arr(i,j,k);
2224  amrex::Real rho = cons_arr(i,j,k,Rho_comp);
2225  amrex::Real tstar = t_star_arr(i,j,k);
2226  amrex::Real moflux = (std::abs(tstar) > eps) ? -rho * tstar * (theta_mean-theta_surf) : 0.0;
2227 
2228  return moflux;
2229  }
2230 
2231  AMREX_GPU_DEVICE
2232  AMREX_FORCE_INLINE
2233  amrex::Real
2234  compute_u_flux (const int& i,
2235  const int& j,
2236  const int& k,
2237  const amrex::Array4<const amrex::Real>& cons_arr,
2238  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
2239  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
2240  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
2241  const amrex::Array4<const amrex::Real>& /*um_arr*/,
2242  const amrex::Array4<const amrex::Real>& u_star_arr) const
2243  {
2244  // NOTE: this is the total stress
2245  amrex::Real rho = 0.5 *( cons_arr(i-1,j,k,Rho_comp) + cons_arr(i ,j,k,Rho_comp) );
2246  amrex::Real ustar = 0.5 * ( u_star_arr(i-1,j,k) + u_star_arr(i,j,k) );
2247  amrex::Real stressx = -rho * ustar * ustar;
2248 
2249  return stressx;
2250  }
2251 
2252  AMREX_GPU_DEVICE
2253  AMREX_FORCE_INLINE
2254  amrex::Real
2255  compute_v_flux (const int& /*i*/,
2256  const int& /*j*/,
2257  const int& /*k*/,
2258  const amrex::Array4<const amrex::Real>& /*cons_arr*/,
2259  const amrex::Array4<const amrex::Real>& /*velx_arr*/,
2260  const amrex::Array4<const amrex::Real>& /*vely_arr*/,
2261  const amrex::Array4<const amrex::Real>& /*umm_arr*/,
2262  const amrex::Array4<const amrex::Real>& /*vm_arr*/,
2263  const amrex::Array4<const amrex::Real>& /*u_star_arr*/) const
2264  {
2265  // NOTE: this is the total stress
2266  amrex::Real stressy = 0.0;
2267 
2268  return stressy;
2269  }
2270 
2271 private:
2272  const amrex::Real eps = 1e-15;
2273 };
2274 #endif
constexpr amrex::Real KAPPA
Definition: ERF_Constants.H:20
constexpr amrex::Real PI
Definition: ERF_Constants.H:6
constexpr amrex::Real CONST_GRAV
Definition: ERF_Constants.H:21
constexpr amrex::Real PIoTwo
Definition: ERF_Constants.H:7
@ num
Definition: ERF_DataStruct.H:22
#define Rho_comp
Definition: ERF_IndexDefines.H:36
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:37
#define RhoQ1_comp
Definition: ERF_IndexDefines.H:42
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:126
amrex::Real Real
Definition: ERF_ShocInterface.H:16
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
@ qv
Definition: ERF_Kessler.H:28
Definition: ERF_MOSTStress.H:187
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:268
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 > &q_star_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:204
most_data mdata
Definition: ERF_MOSTStress.H:265
adiabatic_charnock(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, amrex::Real cnk_a, bool cnk_visc)
Definition: ERF_MOSTStress.H:188
const amrex::Real tol_z
Definition: ERF_MOSTStress.H:267
similarity_funs sfuns
Definition: ERF_MOSTStress.H:266
Definition: ERF_MOSTStress.H:357
adiabatic_donelan(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux)
Definition: ERF_MOSTStress.H:358
similarity_funs sfuns
Definition: ERF_MOSTStress.H:417
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:419
most_data mdata
Definition: ERF_MOSTStress.H:416
const amrex::Real tol
Definition: ERF_MOSTStress.H:418
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 > &q_star_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:370
Definition: ERF_MOSTStress.H:276
adiabatic_mod_charnock(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, amrex::Real depth)
Definition: ERF_MOSTStress.H:277
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 > &q_star_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:292
similarity_funs sfuns
Definition: ERF_MOSTStress.H:347
const amrex::Real tol_z
Definition: ERF_MOSTStress.H:348
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:349
most_data mdata
Definition: ERF_MOSTStress.H:346
Definition: ERF_MOSTStress.H:427
const amrex::Real eps
Definition: ERF_MOSTStress.H:492
const amrex::Real z0_eps
Definition: ERF_MOSTStress.H:493
most_data mdata
Definition: ERF_MOSTStress.H:489
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 > &q_star_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &Hwave_arr, const amrex::Array4< amrex::Real > &Lwave_arr, const amrex::Array4< amrex::Real > &eta_arr) const
Definition: ERF_MOSTStress.H:440
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:495
similarity_funs sfuns
Definition: ERF_MOSTStress.H:490
const amrex::Real tol
Definition: ERF_MOSTStress.H:491
adiabatic_wave_coupled(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux)
Definition: ERF_MOSTStress.H:428
const amrex::Real z0_max
Definition: ERF_MOSTStress.H:494
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 &, const amrex::Array4< const amrex::Real > &z0_arr, const amrex::Array4< const amrex::Real > &umm_arr, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &, const amrex::Array4< amrex::Real > &u_star_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &t_star_arr, const amrex::Array4< amrex::Real > &q_star_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:150
similarity_funs sfuns
Definition: ERF_MOSTStress.H:179
adiabatic(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux)
Definition: ERF_MOSTStress.H:138
most_data mdata
Definition: ERF_MOSTStress.H:178
Definition: ERF_MOSTStress.H:2061
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:2117
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:2069
custom_flux(bool specified_rho_surf)
Definition: ERF_MOSTStress.H:2062
const bool fluxes_include_rho
Definition: ERF_MOSTStress.H:2170
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:2093
const amrex::Real eps
Definition: ERF_MOSTStress.H:2169
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:2144
Definition: ERF_MOSTStress.H:1929
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:1983
donelan_flux()
Definition: ERF_MOSTStress.H:1930
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:2019
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:1956
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:1935
Definition: ERF_MOSTStress.H:1766
moeng_flux()
Definition: ERF_MOSTStress.H:1767
const amrex::Real eps
Definition: ERF_MOSTStress.H:1920
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:1848
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1921
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:1810
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real compute_q_flux(const int &i, const int &j, const int &k, const amrex::Array4< const amrex::Real > &cons_arr, const amrex::Array4< const amrex::Real > &velx_arr, const amrex::Array4< const amrex::Real > &vely_arr, const amrex::Array4< const amrex::Real > &umm_arr, const amrex::Array4< const amrex::Real > &qvm_arr, const amrex::Array4< const amrex::Real > &u_star_arr, const amrex::Array4< const amrex::Real > &q_star_arr, const amrex::Array4< const amrex::Real > &q_surf_arr) const
Definition: ERF_MOSTStress.H:1772
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:1885
Definition: ERF_MOSTStress.H:13
amrex::Real surf_moist_flux
Moisture flux.
Definition: ERF_MOSTStress.H:20
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:24
amrex::Real Cnk_b
Definition: ERF_MOSTStress.H:26
amrex::Real Cnk_d
Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1.
Definition: ERF_MOSTStress.H:25
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:22
const amrex::Real Bjr_beta
Definition: ERF_MOSTStress.H:29
amrex::Real Cnk_b1
Modified Charnock Eq (4) https://doi.org/10.1175/JAMC-D-17-0137.1.
Definition: ERF_MOSTStress.H:23
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:27
amrex::Real surf_temp_flux
Heat flux.
Definition: ERF_MOSTStress.H:19
Definition: ERF_MOSTStress.H:2178
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:2255
rotate_flux()
Definition: ERF_MOSTStress.H:2179
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:2234
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:2209
const amrex::Real eps
Definition: ERF_MOSTStress.H:2272
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 > &qvm_arr, const amrex::Array4< const amrex::Real > &, const amrex::Array4< const amrex::Real > &q_star_arr, const amrex::Array4< const amrex::Real > &q_surf_arr) const
Definition: ERF_MOSTStress.H:2184
Definition: ERF_MOSTStress.H:37
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real calc_psi_m(amrex::Real zeta) const
Definition: ERF_MOSTStress.H:87
amrex::Real beta_m
Constants from Dyer, BLM, 1974.
Definition: ERF_MOSTStress.H:112
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real calc_psi_h2(amrex::Real zeta) const
Definition: ERF_MOSTStress.H:64
amrex::Real beta_h
https://doi.org/10.1007/BF00240838
Definition: ERF_MOSTStress.H:113
amrex::Real gamma_h
Definition: ERF_MOSTStress.H:115
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real calc_psi_h(amrex::Real zeta) const
Definition: ERF_MOSTStress.H:101
amrex::Real gamma_m
Definition: ERF_MOSTStress.H:114
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real calc_psi_m2(amrex::Real zeta) const
Definition: ERF_MOSTStress.H:44
Definition: ERF_MOSTStress.H:608
const amrex::Real tol
Definition: ERF_MOSTStress.H:718
similarity_funs sfuns
Definition: ERF_MOSTStress.H:717
most_data mdata
Definition: ERF_MOSTStress.H:715
bool spec_qflux
Definition: ERF_MOSTStress.H:716
surface_flux_charnock(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, amrex::Real cnk_a, bool cnk_visc, bool cons_qflux)
Definition: ERF_MOSTStress.H:609
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:719
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 > &q_surf_arr, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:627
Definition: ERF_MOSTStress.H:838
bool spec_qflux
Definition: ERF_MOSTStress.H:935
surface_flux_donelan(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, bool cons_qflux)
Definition: ERF_MOSTStress.H:839
const amrex::Real tol
Definition: ERF_MOSTStress.H:937
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:938
most_data mdata
Definition: ERF_MOSTStress.H:934
similarity_funs sfuns
Definition: ERF_MOSTStress.H:936
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 > &q_surf_arr, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:853
Definition: ERF_MOSTStress.H:727
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 > &q_surf_arr, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:745
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:830
most_data mdata
Definition: ERF_MOSTStress.H:826
bool spec_qflux
Definition: ERF_MOSTStress.H:827
surface_flux_mod_charnock(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, amrex::Real depth, bool cons_qflux)
Definition: ERF_MOSTStress.H:728
const amrex::Real tol
Definition: ERF_MOSTStress.H:829
similarity_funs sfuns
Definition: ERF_MOSTStress.H:828
Definition: ERF_MOSTStress.H:946
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1055
surface_flux_wave_coupled(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, bool cons_qflux)
Definition: ERF_MOSTStress.H:947
const amrex::Real z0_max
Definition: ERF_MOSTStress.H:1054
const amrex::Real z0_eps
Definition: ERF_MOSTStress.H:1053
most_data mdata
Definition: ERF_MOSTStress.H:1048
similarity_funs sfuns
Definition: ERF_MOSTStress.H:1050
const amrex::Real tol
Definition: ERF_MOSTStress.H:1051
bool spec_qflux
Definition: ERF_MOSTStress.H:1049
const amrex::Real eps
Definition: ERF_MOSTStress.H:1052
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 > &q_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:961
Definition: ERF_MOSTStress.H:503
similarity_funs sfuns
Definition: ERF_MOSTStress.H:598
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:600
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 > &q_surf_arr, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:518
most_data mdata
Definition: ERF_MOSTStress.H:596
const amrex::Real tol
Definition: ERF_MOSTStress.H:599
bool spec_qflux
Definition: ERF_MOSTStress.H:597
surface_flux(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, bool cons_qflux)
Definition: ERF_MOSTStress.H:504
Definition: ERF_MOSTStress.H:1197
most_data mdata
Definition: ERF_MOSTStress.H:1349
surface_temp_charnock(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, amrex::Real cnk_a, bool cnk_visc, bool cons_qflux)
Definition: ERF_MOSTStress.H:1198
bool spec_qflux
Definition: ERF_MOSTStress.H:1350
const amrex::Real tol_z
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 > &q_surf_arr, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:1216
const amrex::Real alpha
Definition: ERF_MOSTStress.H:1354
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1355
const amrex::Real tol
Definition: ERF_MOSTStress.H:1352
similarity_funs sfuns
Definition: ERF_MOSTStress.H:1351
Definition: ERF_MOSTStress.H:1521
surface_temp_donelan(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, bool cons_qflux)
Definition: ERF_MOSTStress.H:1522
similarity_funs sfuns
Definition: ERF_MOSTStress.H:1629
most_data mdata
Definition: ERF_MOSTStress.H:1627
const amrex::Real tol
Definition: ERF_MOSTStress.H:1630
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1631
bool spec_qflux
Definition: ERF_MOSTStress.H:1628
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 > &q_surf_arr, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:1536
Definition: ERF_MOSTStress.H:1363
similarity_funs sfuns
Definition: ERF_MOSTStress.H:1509
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 > &q_surf_arr, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:1381
const amrex::Real tol
Definition: ERF_MOSTStress.H:1510
const amrex::Real tol_z
Definition: ERF_MOSTStress.H:1511
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1513
surface_temp_mod_charnock(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, amrex::Real depth, bool cons_qflux)
Definition: ERF_MOSTStress.H:1364
most_data mdata
Definition: ERF_MOSTStress.H:1507
bool spec_qflux
Definition: ERF_MOSTStress.H:1508
const amrex::Real alpha
Definition: ERF_MOSTStress.H:1512
Definition: ERF_MOSTStress.H:1639
const amrex::Real eps
Definition: ERF_MOSTStress.H:1755
const amrex::Real tol
Definition: ERF_MOSTStress.H:1754
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 > &q_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:1654
most_data mdata
Definition: ERF_MOSTStress.H:1751
const amrex::Real z0_eps
Definition: ERF_MOSTStress.H:1756
similarity_funs sfuns
Definition: ERF_MOSTStress.H:1753
bool spec_qflux
Definition: ERF_MOSTStress.H:1752
const amrex::Real z0_max
Definition: ERF_MOSTStress.H:1757
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1758
surface_temp_wave_coupled(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, bool cons_qflux)
Definition: ERF_MOSTStress.H:1640
Definition: ERF_MOSTStress.H:1063
similarity_funs sfuns
Definition: ERF_MOSTStress.H:1186
surface_temp(amrex::Real zref, amrex::Real Tflux, amrex::Real Qvflux, bool cons_qflux)
Definition: ERF_MOSTStress.H:1064
const amrex::Real tol
Definition: ERF_MOSTStress.H:1187
bool spec_qflux
Definition: ERF_MOSTStress.H:1185
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 > &q_surf_arr, const amrex::Array4< amrex::Real > &olen_arr, const amrex::Array4< amrex::Real > &pblh_arr, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &, const amrex::Array4< amrex::Real > &) const
Definition: ERF_MOSTStress.H:1078
const amrex::Real alpha
Definition: ERF_MOSTStress.H:1188
const amrex::Real WSMIN
Definition: ERF_MOSTStress.H:1189
most_data mdata
Definition: ERF_MOSTStress.H:1184