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