ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
TurbChoice Struct Reference

#include <ERF_TurbStruct.H>

Collaboration diagram for TurbChoice:

Public Member Functions

void init_params (int lev, int max_level, std::string pp_prefix)
 
void display (int lev)
 

Public Attributes

LESType les_type
 
amrex::Real Pr_t = amrex::Real(1.0) / amrex::Real(3.0)
 
amrex::Real Pr_t_inv = amrex::Real(3.0)
 
amrex::Real Sc_t = 1.0
 
amrex::Real Sc_t_inv = 1.0
 
amrex::Real Cs = 0.0
 
bool smag2d = false
 
amrex::Real Ce = 0.93
 
amrex::Real Ce_wall = 0.0
 
amrex::Real Ck = 0.1
 
amrex::Real Cmu0 = 0.5562
 
amrex::Real Cb = 0.35
 
amrex::Real Rt_crit = -1.0
 
amrex::Real Rt_min = -3.0
 
amrex::Real l_g_max = 1e99
 
amrex::Real sigma_k = 0.5
 
amrex::Real theta_ref = 300.0
 
bool mix_isotropic = true
 
RANSType rans_type
 
PBLType pbl_type
 
MYNNLevel25 pbl_mynn
 
MYNNLevel2 pbl_mynn_level2
 
bool use_kturb = false
 
bool use_keqn
 
bool use_pbl_tke
 
bool use_tke = false
 
amrex::Real pbl_ysu_coriolis_freq
 
bool pbl_ysu_use_consistent_coriolis
 
bool pbl_ysu_force_over_water
 
amrex::Real pbl_ysu_land_Ribcr
 
amrex::Real pbl_ysu_unst_Ribcr
 
amrex::Real pbl_mrf_coriolis_freq = 1.0e-4
 
amrex::Real pbl_mrf_Ribcr
 
amrex::Real pbl_mrf_const_b
 
amrex::Real pbl_mrf_sf
 
bool mrf_moistvars = false
 
bool advect_tke = true
 
bool diffuse_tke_3D = true
 

Detailed Description

Container holding quantities related to turbulence parametrizations

Member Function Documentation

◆ display()

void TurbChoice::display ( int  lev)
inline
198  {
199  amrex::Print() << "Turbulence Settings at level " << lev << std::endl;
200 
201  if (
202  les_type == LESType::None && rans_type == RANSType::None &&
203  pbl_type == PBLType::None) {
204  amrex::Print() << " Using DNS model at level " << lev << std::endl;
205  } else if (les_type == LESType::Smagorinsky) {
206  if (smag2d) {
207  amrex::Print() << " Using 2D Smagorinsky LES model at level " << lev
208  << std::endl;
209  } else {
210  amrex::Print() << " Using Smagorinsky LES model at level " << lev
211  << std::endl;
212  }
213  } else if (les_type == LESType::Deardorff) {
214  amrex::Print() << " Using Deardorff LES model at level " << lev
215  << std::endl;
216  } else if (rans_type == RANSType::kEqn) {
217  amrex::Print()
218  << " Using Axell & Liungman one-equation RANS k model at level "
219  << lev << std::endl;
220  } else if (pbl_type == PBLType::MYNN25) {
221  amrex::Print() << " Using MYNN2.5 PBL model at level " << lev
222  << std::endl;
223  } else if (pbl_type == PBLType::MYNNEDMF) {
224  amrex::Print() << " Using MYNNEDMF PBL model at level " << lev
225  << std::endl;
226  } else if (pbl_type == PBLType::YSU) {
227  amrex::Print() << " Using YSU PBL model at level " << lev << std::endl;
228  } else if (pbl_type == PBLType::MRF) {
229  amrex::Print() << " Using MRF PBL model at level " << lev << std::endl;
230  } else {
231  amrex::Error("Unknown turbulence model");
232  }
233 
234  if (les_type != LESType::None) {
235  if (les_type == LESType::Smagorinsky) {
236  amrex::Print() << "Cs : " << Cs << std::endl;
237  }
238  if (les_type == LESType::Deardorff) {
239  amrex::Print() << "Ce : " << Ce << std::endl;
240  amrex::Print() << "Ce at wall : " << Ce_wall
241  << std::endl;
242  amrex::Print() << "Ck : " << Ck << std::endl;
243  amrex::Print() << "sigma_k : " << sigma_k
244  << std::endl;
245  amrex::Print() << "reference theta : " << theta_ref
246  << std::endl;
247  // Sullivan et al 1994, Eqn 14
248  amrex::Real Cs_equiv = std::sqrt(Ck * std::sqrt(Ck / Ce));
249  amrex::Print() << "equivalent Cs : " << Cs_equiv
250  << std::endl;
251  }
252  amrex::Print() << "isotropic mixing : " << mix_isotropic
253  << std::endl;
254  }
255 
256  if (rans_type != RANSType::None) {
257  if (rans_type == RANSType::kEqn) {
258  amrex::Print() << "Cmu0 : " << Cmu0 << std::endl;
259  amrex::Print() << "sigma_k : " << sigma_k
260  << std::endl;
261  amrex::Print() << "Cb : " << Cb << std::endl;
262  amrex::Print() << "Rt_crit : " << Rt_crit
263  << std::endl;
264  amrex::Print() << "Rt_min : " << Rt_min
265  << std::endl;
266  amrex::Print() << "max_geom_lscale : " << l_g_max
267  << std::endl;
268  amrex::Print() << "reference theta : " << theta_ref
269  << std::endl;
270  }
271  }
272 
273  if ((les_type != LESType::None) || (rans_type != RANSType::None)) {
274  amrex::Print() << "Pr_t : " << Pr_t << std::endl;
275  amrex::Print() << "Sc_t : " << Sc_t << std::endl;
276  }
277 
278  if (pbl_type == PBLType::MYNN25 || pbl_type == PBLType::MYNNEDMF) {
279  amrex::Print() << " pbl_mynn_A1 : " << pbl_mynn.A1 << std::endl;
280  amrex::Print() << " pbl_mynn_A2 : " << pbl_mynn.A2 << std::endl;
281  amrex::Print() << " pbl_mynn_B1 : " << pbl_mynn.B1 << std::endl;
282  amrex::Print() << " pbl_mynn_B2 : " << pbl_mynn.B2 << std::endl;
283  amrex::Print() << " pbl_mynn_C1 : " << pbl_mynn.C1 << std::endl;
284  amrex::Print() << " pbl_mynn_C2 : " << pbl_mynn.C2 << std::endl;
285  amrex::Print() << " pbl_mynn_C3 : " << pbl_mynn.C3 << std::endl;
286  amrex::Print() << " pbl_mynn_C4 : " << pbl_mynn.C4 << std::endl;
287  amrex::Print() << " pbl_mynn_C5 : " << pbl_mynn.C5 << std::endl;
288  } else if (pbl_type == PBLType::YSU) {
289  amrex::Print() << " pbl_ysu_coriolis_freq : "
290  << pbl_ysu_coriolis_freq << std::endl;
291  amrex::Print() << " pbl_ysu_use_consistent_coriolis : "
292  << pbl_ysu_use_consistent_coriolis << std::endl;
293  amrex::Print() << " pbl_ysu_force_over_water : "
294  << pbl_ysu_force_over_water << std::endl;
295  amrex::Print() << " pbl_ysu_land_Ribcr : "
296  << pbl_ysu_land_Ribcr << std::endl;
297  amrex::Print() << " pbl_ysu_unst_Ribcr : "
298  << pbl_ysu_unst_Ribcr << std::endl;
299  } else if (pbl_type == PBLType::MRF) {
300  amrex::Print() << " pbl_mrf_coriolis_freq : " << pbl_mrf_coriolis_freq
301  << std::endl;
302  amrex::Print() << " pbl_mrf_Ribcr : " << pbl_mrf_Ribcr
303  << std::endl;
304  amrex::Print() << " pbl_mrf_const_b : " << pbl_mrf_const_b
305  << std::endl;
306  amrex::Print() << " pbl_mrf_sf : " << pbl_mrf_sf
307  << std::endl;
308  amrex::Print() << " mrf_moistvars : " << mrf_moistvars
309  << std::endl;
310  }
311  }
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 B1
Definition: ERF_MYNNStruct.H:43
amrex::Real B2
Definition: ERF_MYNNStruct.H:44
amrex::Real C5
Definition: ERF_MYNNStruct.H:49
amrex::Real A1
Definition: ERF_MYNNStruct.H:41
bool smag2d
Definition: ERF_TurbStruct.H:326
amrex::Real sigma_k
Definition: ERF_TurbStruct.H:342
MYNNLevel25 pbl_mynn
Definition: ERF_TurbStruct.H:356
PBLType pbl_type
Definition: ERF_TurbStruct.H:354
amrex::Real pbl_mrf_const_b
Definition: ERF_TurbStruct.H:387
amrex::Real Rt_min
Definition: ERF_TurbStruct.H:337
amrex::Real pbl_mrf_coriolis_freq
Definition: ERF_TurbStruct.H:384
RANSType rans_type
Definition: ERF_TurbStruct.H:351
amrex::Real pbl_mrf_Ribcr
Definition: ERF_TurbStruct.H:385
amrex::Real Ck
Definition: ERF_TurbStruct.H:331
amrex::Real Cmu0
Definition: ERF_TurbStruct.H:334
bool pbl_ysu_use_consistent_coriolis
Definition: ERF_TurbStruct.H:373
amrex::Real Cb
Definition: ERF_TurbStruct.H:335
amrex::Real pbl_mrf_sf
Definition: ERF_TurbStruct.H:389
amrex::Real pbl_ysu_land_Ribcr
Definition: ERF_TurbStruct.H:379
bool mrf_moistvars
Definition: ERF_TurbStruct.H:391
amrex::Real Cs
Definition: ERF_TurbStruct.H:325
amrex::Real Ce_wall
Definition: ERF_TurbStruct.H:330
amrex::Real Ce
Definition: ERF_TurbStruct.H:329
LESType les_type
Definition: ERF_TurbStruct.H:314
amrex::Real pbl_ysu_coriolis_freq
Definition: ERF_TurbStruct.H:370
bool pbl_ysu_force_over_water
Definition: ERF_TurbStruct.H:376
amrex::Real Pr_t
Definition: ERF_TurbStruct.H:317
amrex::Real theta_ref
Definition: ERF_TurbStruct.H:345
amrex::Real Sc_t
Definition: ERF_TurbStruct.H:321
amrex::Real Rt_crit
Definition: ERF_TurbStruct.H:336
amrex::Real pbl_ysu_unst_Ribcr
Definition: ERF_TurbStruct.H:381
bool mix_isotropic
Definition: ERF_TurbStruct.H:348
amrex::Real l_g_max
Definition: ERF_TurbStruct.H:338

◆ init_params()

void TurbChoice::init_params ( int  lev,
int  max_level,
std::string  pp_prefix 
)
inline
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);
108  } else if (pbl_type == PBLType::YSU) {
110  pp, "pbl_ysu_coriolis_freq", pbl_ysu_coriolis_freq, lev, max_level);
112  pp, "pbl_ysu_use_consistent_coriolis",
113  pbl_ysu_use_consistent_coriolis, lev, max_level);
115  pp, "pbl_ysu_force_over_water", pbl_ysu_force_over_water, lev,
116  max_level);
118  pp, "pbl_ysu_land_Ribcr", pbl_ysu_land_Ribcr, lev, max_level);
120  pp, "pbl_ysu_unst_Ribcr", pbl_ysu_unst_Ribcr, lev, max_level);
121  } else if (pbl_type == PBLType::MRF) {
123  pp, "pbl_mrf_coriolis_freq", pbl_mrf_coriolis_freq, lev, max_level);
125  pp, "pbl_mrf_Ribcr", pbl_mrf_Ribcr, lev, max_level);
127  pp, "pbl_mrf_const_b", pbl_mrf_const_b, lev, max_level);
128  query_one_or_per_level(pp, "pbl_mrf_sf", pbl_mrf_sf, lev, max_level);
130  pp, "mrf_moistvars", mrf_moistvars, lev, max_level);
131  }
132  }
133 
134  // Right now, solving the QKE equation is only supported when MYNN PBL is
135  // turned on
136  if (pbl_type == PBLType::MYNN25 || pbl_type == PBLType::MYNNEDMF) {
137  use_pbl_tke = true;
138 
139  query_one_or_per_level(pp, "advect_tke", advect_tke, lev, max_level);
141  pp, "diffuse_tke_3D", diffuse_tke_3D, lev, max_level);
142  }
143 
144  // LES constants...
145  query_one_or_per_level(pp, "Cs", Cs, lev, max_level);
146 
147  query_one_or_per_level(pp, "Pr_t", Pr_t, lev, max_level);
148  query_one_or_per_level(pp, "Sc_t", Sc_t, lev, max_level);
149 
150  // Compute relevant forms of diffusion parameters
151  Pr_t_inv = amrex::Real(1.0) / Pr_t;
152  Sc_t_inv = amrex::Real(1.0) / Sc_t;
153 
154  query_one_or_per_level(pp, "Ce", Ce, lev, max_level);
155  query_one_or_per_level(pp, "Ce_wall", Ce_wall, lev, max_level);
156 
157  if (les_type == LESType::Deardorff) {
158  query_one_or_per_level(pp, "Ck", Ck, lev, max_level);
159  }
160 
161  // k-eqn constants
162  query_one_or_per_level(pp, "Cmu0", Cmu0, lev, max_level);
163  query_one_or_per_level(pp, "Cb", Cb, lev, max_level);
164  query_one_or_per_level(pp, "Rt_crit", Rt_crit, lev, max_level);
165  query_one_or_per_level(pp, "Rt_min", Rt_min, lev, max_level);
166  query_one_or_per_level(pp, "max_geom_lscale", l_g_max, lev, max_level);
167 
168  // Common inputs (LES or RANS)
169  query_one_or_per_level(pp, "sigma_k", sigma_k, lev, max_level);
170  query_one_or_per_level(pp, "theta_ref", theta_ref, lev, max_level);
171 
172  query_one_or_per_level(pp, "mix_isotropic", mix_isotropic, lev, max_level);
173 
174  // Set common flags
175  use_kturb =
176  ((les_type != LESType::None) || (rans_type != RANSType::None) ||
177  (pbl_type != PBLType::None));
178  use_keqn =
179  ((les_type == LESType::Deardorff) || (rans_type == RANSType::kEqn));
180  use_tke =
181  ((les_type == LESType::Deardorff) || (rans_type == RANSType::kEqn) ||
182  (pbl_type == PBLType::MYNN25) || (pbl_type == PBLType::MYNNEDMF));
183 
184  // Validate inputs
185  if (les_type == LESType::Smagorinsky) {
186  if (Cs == 0) {
187  amrex::Error("Need to specify Cs for Smagorsinky LES");
188  }
189  if (smag2d && mix_isotropic) {
190  amrex::Print() << "Turning off mix_isotropic for 2-D Smagorinsky"
191  << std::endl;
192  mix_isotropic = false;
193  }
194  }
195  }
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real pp(amrex::Real y)
Definition: ERF_MicrophysicsUtils.H:219
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
amrex::Real SMmax
Definition: ERF_MYNNStruct.H:53
amrex::Real SHmax
Definition: ERF_MYNNStruct.H:55
amrex::Real SHmin
Definition: ERF_MYNNStruct.H:54
amrex::Real SMmin
Definition: ERF_MYNNStruct.H:52
bool diffuse_moistvars
Definition: ERF_MYNNStruct.H:60
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:64
bool use_pbl_tke
Definition: ERF_TurbStruct.H:364
amrex::Real Sc_t_inv
Definition: ERF_TurbStruct.H:322
bool use_keqn
Definition: ERF_TurbStruct.H:361
MYNNLevel2 pbl_mynn_level2
Definition: ERF_TurbStruct.H:357
bool advect_tke
Definition: ERF_TurbStruct.H:394
bool use_tke
Definition: ERF_TurbStruct.H:366
bool diffuse_tke_3D
Definition: ERF_TurbStruct.H:396
amrex::Real Pr_t_inv
Definition: ERF_TurbStruct.H:318
bool use_kturb
Definition: ERF_TurbStruct.H:360
Here is the call graph for this function:

Member Data Documentation

◆ advect_tke

bool TurbChoice::advect_tke = true

Referenced by erf_slow_rhs_post(), and init_params().

◆ Cb

amrex::Real TurbChoice::Cb = 0.35

◆ Ce

amrex::Real TurbChoice::Ce = 0.93

◆ Ce_wall

amrex::Real TurbChoice::Ce_wall = 0.0

◆ Ck

amrex::Real TurbChoice::Ck = 0.1

◆ Cmu0

amrex::Real TurbChoice::Cmu0 = 0.5562

◆ Cs

amrex::Real TurbChoice::Cs = 0.0

◆ diffuse_tke_3D

bool TurbChoice::diffuse_tke_3D = true

Referenced by init_params(), and make_sources().

◆ l_g_max

amrex::Real TurbChoice::l_g_max = 1e99

◆ les_type

◆ mix_isotropic

bool TurbChoice::mix_isotropic = true

◆ mrf_moistvars

bool TurbChoice::mrf_moistvars = false

Referenced by display(), and init_params().

◆ pbl_mrf_const_b

amrex::Real TurbChoice::pbl_mrf_const_b
Initial value:
=
7.8

Referenced by ComputeDiffusivityMRF(), display(), and init_params().

◆ pbl_mrf_coriolis_freq

amrex::Real TurbChoice::pbl_mrf_coriolis_freq = 1.0e-4

Referenced by display(), and init_params().

◆ pbl_mrf_Ribcr

amrex::Real TurbChoice::pbl_mrf_Ribcr
Initial value:
=
0.5

Referenced by ComputeDiffusivityMRF(), display(), and init_params().

◆ pbl_mrf_sf

amrex::Real TurbChoice::pbl_mrf_sf
Initial value:
=
0.1

Referenced by ComputeDiffusivityMRF(), display(), and init_params().

◆ pbl_mynn

◆ pbl_mynn_level2

MYNNLevel2 TurbChoice::pbl_mynn_level2

◆ pbl_type

◆ pbl_ysu_coriolis_freq

amrex::Real TurbChoice::pbl_ysu_coriolis_freq
Initial value:
=
1.0e-4

Referenced by ComputeDiffusivityYSU(), display(), and init_params().

◆ pbl_ysu_force_over_water

bool TurbChoice::pbl_ysu_force_over_water
Initial value:
=
false

Referenced by ComputeDiffusivityYSU(), display(), and init_params().

◆ pbl_ysu_land_Ribcr

amrex::Real TurbChoice::pbl_ysu_land_Ribcr
Initial value:
=
0.25

Referenced by ComputeDiffusivityYSU(), display(), and init_params().

◆ pbl_ysu_unst_Ribcr

amrex::Real TurbChoice::pbl_ysu_unst_Ribcr
Initial value:
=
0.0

Referenced by ComputeDiffusivityYSU(), display(), and init_params().

◆ pbl_ysu_use_consistent_coriolis

bool TurbChoice::pbl_ysu_use_consistent_coriolis
Initial value:
=
false

Referenced by display(), and init_params().

◆ Pr_t

amrex::Real TurbChoice::Pr_t = amrex::Real(1.0) / amrex::Real(3.0)

Referenced by display(), and init_params().

◆ Pr_t_inv

amrex::Real TurbChoice::Pr_t_inv = amrex::Real(3.0)

◆ rans_type

◆ Rt_crit

amrex::Real TurbChoice::Rt_crit = -1.0

◆ Rt_min

amrex::Real TurbChoice::Rt_min = -3.0

◆ Sc_t

amrex::Real TurbChoice::Sc_t = 1.0

Referenced by display(), and init_params().

◆ Sc_t_inv

amrex::Real TurbChoice::Sc_t_inv = 1.0

◆ sigma_k

amrex::Real TurbChoice::sigma_k = 0.5

◆ smag2d

bool TurbChoice::smag2d = false

◆ theta_ref

amrex::Real TurbChoice::theta_ref = 300.0

◆ use_keqn

bool TurbChoice::use_keqn
Initial value:
=
false

Referenced by erf_slow_rhs_pre(), and init_params().

◆ use_kturb

bool TurbChoice::use_kturb = false

◆ use_pbl_tke

bool TurbChoice::use_pbl_tke
Initial value:
=
false

Referenced by init_params().

◆ use_tke

bool TurbChoice::use_tke = false

The documentation for this struct was generated from the following file: