ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_TurbStruct.H
Go to the documentation of this file.
1 #ifndef ERF_TURB_STRUCT_H_
2 #define ERF_TURB_STRUCT_H_
3 
4 #include <ERF_MYNNStruct.H>
5 
6 AMREX_ENUM(LESType, None, Smagorinsky, Smagorinsky2D, Deardorff);
7 
8 AMREX_ENUM(RANSType, None, kEqn);
9 
10 AMREX_ENUM(PBLType, None, MYJ, MYNN25, MYNNEDMF, YSU, MRF);
11 
12 AMREX_ENUM(StratType, theta, thetav, thetal);
13 
14 template <typename T>
15 void
17  const amrex::ParmParse& pp,
18  const char* query_string,
19  T& query_var,
20  const int lev,
21  const int maxlev)
22 {
23  int count = pp.countval(query_string);
24  if (count == 0) {
25  return; // nothing to do
26  } else if (count == 1) {
27  pp.query(query_string, query_var);
28  } else if (count == maxlev + 1) {
29  pp.query(query_string, query_var, lev);
30  } else {
31  amrex::Error(
32  "For parmparse variable " + pp.prefixedName(query_string) +
33  ": if specified, specify once total or once for each level");
34  }
35 }
36 
37 /**
38  * Container holding quantities related to turbulence parametrizations
39  */
40 struct TurbChoice
41 {
42 public:
43  void init_params (int lev, int max_level, std::string pp_prefix)
44  {
45  amrex::ParmParse pp(pp_prefix);
46 
47  // Which LES closure?
48  std::string les_type_string = "None";
49  query_one_or_per_level(pp, "les_type", les_type, lev, max_level);
50 
51  // Handle 2-D Smag
52  if (les_type == LESType::Smagorinsky2D) {
53  les_type = LESType::Smagorinsky;
54  smag2d = true;
55  }
56 
57  // Which RANS closure?
58  std::string rans_type_string = "None";
59  query_one_or_per_level(pp, "rans_type", rans_type, lev, max_level);
60 
61  if ((rans_type != RANSType::None) && (les_type != LESType::None)) {
62  amrex::Error("Hybrid RANS-LES not implemented");
63  }
64 
65  // Which PBL Closure
66  static std::string pbl_type_string = "None";
67  query_one_or_per_level(pp, "pbl_type", pbl_type, lev, max_level);
68 
69  // Do some more stuff for PBL Modeling
70  if (pbl_type != PBLType::None) {
71  // Check for compatibility between PBL, LES, Molec Transport
72  if (les_type != LESType::None) {
73  amrex::Print() << "Selected a PBL model and an LES model: "
74  << "Using PBL for vertical transport, LES for horizontal"
75  << std::endl;
76  }
77  if (les_type == LESType::Smagorinsky) {
78  if (!smag2d)
79  amrex::Error("If using Smagorinsky with a PBL model, the 2-D "
80  "formulation should be used");
81  } else if (les_type == LESType::Deardorff) {
82  amrex::Error(
83  "It is not recommended to use Deardorff LES and a PBL model");
84  }
85 
86  if (pbl_type == PBLType::MYNN25 || pbl_type == PBLType::MYNNEDMF) {
87  query_one_or_per_level(pp, "pbl_mynn_A1", pbl_mynn.A1, lev, max_level);
88  query_one_or_per_level(pp, "pbl_mynn_A2", pbl_mynn.A2, lev, max_level);
89  query_one_or_per_level(pp, "pbl_mynn_B1", pbl_mynn.B1, lev, max_level);
90  query_one_or_per_level(pp, "pbl_mynn_B2", pbl_mynn.B2, lev, max_level);
91  query_one_or_per_level(pp, "pbl_mynn_C1", pbl_mynn.C1, lev, max_level);
92  query_one_or_per_level(pp, "pbl_mynn_C2", pbl_mynn.C2, lev, max_level);
93  query_one_or_per_level(pp, "pbl_mynn_C3", pbl_mynn.C3, lev, max_level);
94  query_one_or_per_level(pp, "pbl_mynn_C4", pbl_mynn.C4, lev, max_level);
95  query_one_or_per_level(pp, "pbl_mynn_C5", pbl_mynn.C5, lev, max_level);
100  pp, "pbl_mynn_diffuse_moistvars", pbl_mynn.diffuse_moistvars, lev,
101  max_level);
103  pp, "pbl_mynn_SMmin", pbl_mynn.SMmin, lev, max_level);
105  pp, "pbl_mynn_SMmax", pbl_mynn.SMmax, lev, max_level);
107  pp, "pbl_mynn_SHmin", pbl_mynn.SHmin, lev, max_level);
109  pp, "pbl_mynn_SHmax", pbl_mynn.SHmax, lev, max_level);
111  pp, "pbl_mynn_SQfactor", pbl_mynn.SQfac, lev, max_level);
112  } else if (pbl_type == PBLType::YSU) {
114  pp, "pbl_ysu_coriolis_freq", pbl_ysu_coriolis_freq, lev, max_level);
116  pp, "pbl_ysu_use_consistent_coriolis",
117  pbl_ysu_use_consistent_coriolis, lev, max_level);
119  pp, "pbl_ysu_force_over_water", pbl_ysu_force_over_water, lev,
120  max_level);
122  pp, "pbl_ysu_land_Ribcr", pbl_ysu_land_Ribcr, lev, max_level);
124  pp, "pbl_ysu_unst_Ribcr", pbl_ysu_unst_Ribcr, lev, max_level);
125  } else if (pbl_type == PBLType::MRF) {
127  pp, "pbl_mrf_coriolis_freq", pbl_mrf_coriolis_freq, lev, max_level);
129  pp, "pbl_mrf_Ribcr", pbl_mrf_Ribcr, lev, max_level);
131  pp, "pbl_mrf_const_b", pbl_mrf_const_b, lev, max_level);
132  query_one_or_per_level(pp, "pbl_mrf_sf", pbl_mrf_sf, lev, max_level);
134  pp, "mrf_moistvars", mrf_moistvars, lev, max_level);
135  }
136  }
137 
138  // Flags for QKE/TKE euation
139  if (pbl_type == PBLType::MYJ || pbl_type == PBLType::MYNN25 || pbl_type == PBLType::MYNNEDMF) {
140  // Add sources/sinks to QKE/TKE? (MYJ does this inline)
141  if (pbl_type == PBLType::MYNN25 || pbl_type == PBLType::MYNNEDMF) {
142  use_pbl_tke = true;
143  }
144  // Advect QKE/TKE?
145  query_one_or_per_level(pp, "advect_tke" , advect_tke , lev, max_level);
146  // Apply numerical diffusion to QKE/TKE?
147  query_one_or_per_level(pp, "diffuse_tke_3D", diffuse_tke_3D, lev, max_level);
148  }
149 
150  // LES constants...
151  query_one_or_per_level(pp, "Cs", Cs, lev, max_level);
152 
153  query_one_or_per_level(pp, "Pr_t", Pr_t, lev, max_level);
154  query_one_or_per_level(pp, "Sc_t", Sc_t, lev, max_level);
155 
156  // Compute relevant forms of diffusion parameters
157  Pr_t_inv = amrex::Real(1.0) / Pr_t;
158  Sc_t_inv = amrex::Real(1.0) / Sc_t;
159 
160  if (les_type == LESType::Deardorff) {
161  query_one_or_per_level(pp, "Ck", Ck, lev, max_level);
162  query_one_or_per_level(pp, "Ce", Ce, lev, max_level);
163  query_one_or_per_level(pp, "Ce_wall", Ce_wall, lev, max_level);
164  }
165 
166  // To quantify atmospheric stability for subgrid modeling
167  query_one_or_per_level(pp, "thermal_stratification", strat_type, lev, max_level);
168  if (strat_type == StratType::theta) {
169  amrex::Print() << "Thermal stratification based on gradient of potential temperature" << std::endl;
170  } else if (strat_type == StratType::thetav) {
171  amrex::Print() << "Thermal stratification based on gradient of virtual potential temperature" << std::endl;
172  } else if (strat_type == StratType::thetal) {
173  amrex::Print() << "Thermal stratification based on gradient of linearized liquid-water potential temperature" << std::endl;
174  }
175 
176  // k-eqn constants
177  query_one_or_per_level(pp, "Cmu0", Cmu0, lev, max_level);
178  query_one_or_per_level(pp, "Cb", Cb, lev, max_level);
179  query_one_or_per_level(pp, "Rt_crit", Rt_crit, lev, max_level);
180  query_one_or_per_level(pp, "Rt_min", Rt_min, lev, max_level);
181  query_one_or_per_level(pp, "max_geom_lscale", l_g_max, lev, max_level);
182 
183  // Common inputs (LES or RANS)
184  query_one_or_per_level(pp, "sigma_k", sigma_k, lev, max_level);
185  query_one_or_per_level(pp, "theta_ref", theta_ref, lev, max_level);
186 
187  query_one_or_per_level(pp, "mix_isotropic", mix_isotropic, lev, max_level);
188  query_one_or_per_level(pp, "use_smag_stratification", use_smag_stratification, lev, max_level);
189 
190  // Set common flags
191  use_kturb =
192  ((les_type != LESType::None) || (rans_type != RANSType::None) ||
193  (pbl_type != PBLType::None));
194  use_keqn =
195  ((les_type == LESType::Deardorff) || (rans_type == RANSType::kEqn));
196  use_tke =
197  ((les_type == LESType::Deardorff) || (rans_type == RANSType::kEqn) ||
198  (pbl_type == PBLType::MYJ) || (pbl_type == PBLType::MYNN25) ||
199  (pbl_type == PBLType::MYNNEDMF));
200 
201  // Validate inputs
202  if (les_type == LESType::Smagorinsky) {
203  if (Cs == 0) {
204  amrex::Error("Need to specify Cs for Smagorsinky LES");
205  }
206  if (smag2d && mix_isotropic) {
207  amrex::Print() << "Turning off mix_isotropic for 2-D Smagorinsky" << std::endl;
208  mix_isotropic = false;
209  }
210  }
211  }
212 
213  void display (int lev)
214  {
215  amrex::Print() << "Turbulence Settings at level " << lev << std::endl;
216 
217  if (
218  les_type == LESType::None && rans_type == RANSType::None &&
219  pbl_type == PBLType::None) {
220  amrex::Print() << " Using DNS model at level " << lev << std::endl;
221  } else if (les_type == LESType::Smagorinsky) {
222  if (smag2d) {
223  amrex::Print() << " Using 2D Smagorinsky LES model at level " << lev << std::endl;
224  } else {
225  amrex::Print() << " Using Smagorinsky LES model at level " << lev << std::endl;
226  }
228  amrex::Print() << " Smagorinsky accounts for moist stratification" << std::endl;
229  }
230  } else if (les_type == LESType::Deardorff) {
231  amrex::Print() << " Using Deardorff LES model at level " << lev << std::endl;
232  } else if (rans_type == RANSType::kEqn) {
233  amrex::Print()
234  << " Using Axell & Liungman one-equation RANS k model at level " << lev << std::endl;
235  } else if (pbl_type == PBLType::MYJ) {
236  amrex::Print() << " Using MYJ PBL model at level " << lev << std::endl;
237  } else if (pbl_type == PBLType::MYNN25) {
238  amrex::Print() << " Using MYNN2.5 PBL model at level " << lev << std::endl;
239  } else if (pbl_type == PBLType::MYNNEDMF) {
240  amrex::Print() << " Using MYNNEDMF PBL model at level " << lev << std::endl;
241  } else if (pbl_type == PBLType::YSU) {
242  amrex::Print() << " Using YSU PBL model at level " << lev << std::endl;
243  } else if (pbl_type == PBLType::MRF) {
244  amrex::Print() << " Using MRF PBL model at level " << lev << std::endl;
245  } else {
246  amrex::Error("Unknown turbulence model");
247  }
248 
249  if (les_type != LESType::None) {
250  if (les_type == LESType::Smagorinsky) {
251  amrex::Print() << " Cs : " << Cs << std::endl;
252  }
253  if (les_type == LESType::Deardorff) {
254  amrex::Print() << " Ce : " << Ce << std::endl;
255  amrex::Print() << " Ce at wall : " << Ce_wall << std::endl;
256  amrex::Print() << " Ck : " << Ck << std::endl;
257  amrex::Print() << " sigma_k : " << sigma_k << std::endl;
258 
259  // Sullivan et al 1994, Eqn 14
260  amrex::Real Cs_equiv = std::sqrt(Ck * std::sqrt(Ck / Ce));
261  amrex::Print() << " equivalent Cs : " << Cs_equiv
262  << std::endl;
263  }
264  amrex::Print() << " isotropic mixing : " << mix_isotropic
265  << std::endl;
266  }
267 
268  if (rans_type != RANSType::None) {
269  if (rans_type == RANSType::kEqn) {
270  amrex::Print() << "Cmu0 : " << Cmu0 << std::endl;
271  amrex::Print() << "sigma_k : " << sigma_k << std::endl;
272  amrex::Print() << "Cb : " << Cb << std::endl;
273  amrex::Print() << "Rt_crit : " << Rt_crit << std::endl;
274  amrex::Print() << "Rt_min : " << Rt_min << std::endl;
275  amrex::Print() << "max_geom_lscale : " << l_g_max << std::endl;
276  }
277  }
278 
279  if ((les_type == LESType::Deardorff) ||
280  (les_type == LESType::Smagorinsky && use_smag_stratification) ||
281  (rans_type == RANSType::kEqn)) {
282  if (theta_ref > 0) {
283  amrex::Print() << " reference theta : " << theta_ref << std::endl;
284  } else {
285  amrex::Print() << " reference theta : n/a" << std::endl;
286  }
287  }
288 
289  if ((les_type != LESType::None) || (rans_type != RANSType::None)) {
290  amrex::Print() << " Pr_t : " << Pr_t << std::endl;
291  amrex::Print() << " Sc_t : " << Sc_t << std::endl;
292  }
293 
294  if (pbl_type == PBLType::MYNN25 || pbl_type == PBLType::MYNNEDMF) {
295  amrex::Print() << " pbl_mynn_A1 : " << pbl_mynn.A1 << std::endl;
296  amrex::Print() << " pbl_mynn_A2 : " << pbl_mynn.A2 << std::endl;
297  amrex::Print() << " pbl_mynn_B1 : " << pbl_mynn.B1 << std::endl;
298  amrex::Print() << " pbl_mynn_B2 : " << pbl_mynn.B2 << std::endl;
299  amrex::Print() << " pbl_mynn_C1 : " << pbl_mynn.C1 << std::endl;
300  amrex::Print() << " pbl_mynn_C2 : " << pbl_mynn.C2 << std::endl;
301  amrex::Print() << " pbl_mynn_C3 : " << pbl_mynn.C3 << std::endl;
302  amrex::Print() << " pbl_mynn_C4 : " << pbl_mynn.C4 << std::endl;
303  amrex::Print() << " pbl_mynn_C5 : " << pbl_mynn.C5 << std::endl;
304  } else if (pbl_type == PBLType::YSU) {
305  amrex::Print() << " pbl_ysu_coriolis_freq : "
306  << pbl_ysu_coriolis_freq << std::endl;
307  amrex::Print() << " pbl_ysu_use_consistent_coriolis : "
308  << pbl_ysu_use_consistent_coriolis << std::endl;
309  amrex::Print() << " pbl_ysu_force_over_water : "
310  << pbl_ysu_force_over_water << std::endl;
311  amrex::Print() << " pbl_ysu_land_Ribcr : "
312  << pbl_ysu_land_Ribcr << std::endl;
313  amrex::Print() << " pbl_ysu_unst_Ribcr : "
314  << pbl_ysu_unst_Ribcr << std::endl;
315  } else if (pbl_type == PBLType::MRF) {
316  amrex::Print() << " pbl_mrf_coriolis_freq : " << pbl_mrf_coriolis_freq
317  << std::endl;
318  amrex::Print() << " pbl_mrf_Ribcr : " << pbl_mrf_Ribcr
319  << std::endl;
320  amrex::Print() << " pbl_mrf_const_b : " << pbl_mrf_const_b
321  << std::endl;
322  amrex::Print() << " pbl_mrf_sf : " << pbl_mrf_sf
323  << std::endl;
324  amrex::Print() << " mrf_moistvars : " << mrf_moistvars
325  << std::endl;
326  }
327  }
328 
329  // LES model
330  LESType les_type;
331 
332  // Turbulent Prandtl number
335 
336  // Turbulent Schmidt number
339 
340  // Smagorinsky
342  bool smag2d = false;
343 
344  // Deardorff
345  amrex::Real Ce = 0.93;
346  amrex::Real Ce_wall = 0.0; // if > 0, then set Ce to this at k=0
348 
349  // k-eqn RANS coefficients (Axell & Liungman 2001)
350  amrex::Real Cmu0 = 0.5562;
351  amrex::Real Cb = 0.35;
355 
356  // Deardorff or k-eqn RANS
357  // - diffusivity of tke is 1/sigma_k times the eddy viscosity
359  // - reference potential temperature used to quantify the stratification in
360  // a stable region
362  // - how to quantify stratification effects
364 
365  // Anistropic length scales
366  bool mix_isotropic = true;
367 
368  // Use the stratification adjustment for the Smagorinsky Scheme
370 
371  // RANS type
372  RANSType rans_type;
373 
374  // PBL model
375  PBLType pbl_type;
376 
378  MYNNLevel2 pbl_mynn_level2; // for filtering
379 
380  // Common Flags
381  bool use_kturb = false; // Any turbulence modeling?
382  bool use_keqn =
383  false; // Any microscale turbulence modeling (LES, RANS) with TKE closure?
384  // Then need to populate SmnSmn_lev for production term.
385  bool use_pbl_tke =
386  false; // Any mesoscale turbulence modeling (PBL) with TKE closure?
387  bool use_tke = false; // Any TKE closure (meso or microscale)?
388 
389  // Model coefficients - YSU
390  // TODO: Add parmparse for all of these above
392  1.0e-4; // 1e-4 is hardcoded in WRF, we let the user specify or tske the
393  // value from ERF coriolis forcing
395  false; // ignore input pbl_ysu_coriolis_freq, take value from ERF coriolis
396  // forcing instead
398  false; // Force YSU to act as if it is over water regardless of other inputs
399  // (for testing)
401  0.25; // Critical Bulk Richardson number of Land for stable conditions
403  0.0; // Critical Bulk Richardson number for unstable conditions
404  // Model coefficients - MRF
407  0.5; // Critical Bulk Richardson number for MRF PBL model
409  7.8; // Constant b in MRF PBL model, used to compute the PBL height
411  0.1; // MRF surface flux, used to compute the PBL height
412  bool mrf_moistvars = false; // if true, adds turbulence to moisture
413  // MYNN2.5 PBL model
414  // TKE/QKE stuff
415  bool advect_tke = true; // if MYNN2.5 PBL is used default is turb transport in
416  // Z-direction only
417  bool diffuse_tke_3D = true; // if numerical diffusion is turned on
418 };
419 #endif
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real pp(amrex::Real y)
Definition: ERF_MicrophysicsUtils.H:230
amrex::Real Real
Definition: ERF_ShocInterface.H:19
AMREX_ENUM(LESType, None, Smagorinsky, Smagorinsky2D, Deardorff)
void query_one_or_per_level(const amrex::ParmParse &pp, const char *query_string, T &query_var, const int lev, const int maxlev)
Definition: ERF_TurbStruct.H:16
@ theta
Definition: ERF_MM5.H:20
@ T
Definition: ERF_IndexDefines.H:110
Definition: ERF_MYNNStruct.H:9
amrex::Real SMmax
Definition: ERF_MYNNStruct.H:56
amrex::Real SHmax
Definition: ERF_MYNNStruct.H:58
amrex::Real SQfac
Definition: ERF_MYNNStruct.H:52
amrex::Real C4
Definition: ERF_MYNNStruct.H:48
amrex::Real C1
Definition: ERF_MYNNStruct.H:45
amrex::Real C3
Definition: ERF_MYNNStruct.H:47
amrex::Real C2
Definition: ERF_MYNNStruct.H:46
amrex::Real A2
Definition: ERF_MYNNStruct.H:42
amrex::Real SHmin
Definition: ERF_MYNNStruct.H:57
amrex::Real B1
Definition: ERF_MYNNStruct.H:43
amrex::Real B2
Definition: ERF_MYNNStruct.H:44
amrex::Real C5
Definition: ERF_MYNNStruct.H:49
amrex::Real SMmin
Definition: ERF_MYNNStruct.H:55
amrex::Real A1
Definition: ERF_MYNNStruct.H:41
bool diffuse_moistvars
Definition: ERF_MYNNStruct.H:63
Definition: ERF_MYNNStruct.H:66
void init_coeffs(amrex::Real A1_lvl25, amrex::Real A2_lvl25, amrex::Real B1, amrex::Real B2, amrex::Real C1, amrex::Real C2, amrex::Real C3, amrex::Real, amrex::Real C5)
Definition: ERF_MYNNStruct.H:67
Definition: ERF_TurbStruct.H:41
bool smag2d
Definition: ERF_TurbStruct.H:342
amrex::Real sigma_k
Definition: ERF_TurbStruct.H:358
MYNNLevel25 pbl_mynn
Definition: ERF_TurbStruct.H:377
PBLType pbl_type
Definition: ERF_TurbStruct.H:375
amrex::Real pbl_mrf_const_b
Definition: ERF_TurbStruct.H:408
bool use_pbl_tke
Definition: ERF_TurbStruct.H:385
amrex::Real Sc_t_inv
Definition: ERF_TurbStruct.H:338
amrex::Real Rt_min
Definition: ERF_TurbStruct.H:353
bool use_keqn
Definition: ERF_TurbStruct.H:382
amrex::Real pbl_mrf_coriolis_freq
Definition: ERF_TurbStruct.H:405
MYNNLevel2 pbl_mynn_level2
Definition: ERF_TurbStruct.H:378
RANSType rans_type
Definition: ERF_TurbStruct.H:372
StratType strat_type
Definition: ERF_TurbStruct.H:363
amrex::Real pbl_mrf_Ribcr
Definition: ERF_TurbStruct.H:406
bool advect_tke
Definition: ERF_TurbStruct.H:415
amrex::Real Ck
Definition: ERF_TurbStruct.H:347
amrex::Real Cmu0
Definition: ERF_TurbStruct.H:350
bool pbl_ysu_use_consistent_coriolis
Definition: ERF_TurbStruct.H:394
amrex::Real Cb
Definition: ERF_TurbStruct.H:351
amrex::Real pbl_mrf_sf
Definition: ERF_TurbStruct.H:410
void init_params(int lev, int max_level, std::string pp_prefix)
Definition: ERF_TurbStruct.H:43
amrex::Real pbl_ysu_land_Ribcr
Definition: ERF_TurbStruct.H:400
bool mrf_moistvars
Definition: ERF_TurbStruct.H:412
bool use_tke
Definition: ERF_TurbStruct.H:387
bool diffuse_tke_3D
Definition: ERF_TurbStruct.H:417
amrex::Real Cs
Definition: ERF_TurbStruct.H:341
amrex::Real Pr_t_inv
Definition: ERF_TurbStruct.H:334
amrex::Real Ce_wall
Definition: ERF_TurbStruct.H:346
amrex::Real Ce
Definition: ERF_TurbStruct.H:345
LESType les_type
Definition: ERF_TurbStruct.H:330
amrex::Real pbl_ysu_coriolis_freq
Definition: ERF_TurbStruct.H:391
bool pbl_ysu_force_over_water
Definition: ERF_TurbStruct.H:397
void display(int lev)
Definition: ERF_TurbStruct.H:213
bool use_smag_stratification
Definition: ERF_TurbStruct.H:369
amrex::Real Pr_t
Definition: ERF_TurbStruct.H:333
bool use_kturb
Definition: ERF_TurbStruct.H:381
amrex::Real theta_ref
Definition: ERF_TurbStruct.H:361
amrex::Real Sc_t
Definition: ERF_TurbStruct.H:337
amrex::Real Rt_crit
Definition: ERF_TurbStruct.H:352
amrex::Real pbl_ysu_unst_Ribcr
Definition: ERF_TurbStruct.H:402
bool mix_isotropic
Definition: ERF_TurbStruct.H:366
amrex::Real l_g_max
Definition: ERF_TurbStruct.H:354