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