Skip to content

Commit c48de1a

Browse files
jmsexton03jmsexton03
authored
MYYN-EDMF initial code (#2131)
* Copying tridiag compiles and runs * Duplicating the behavior of MYNN25 with MYNNEDMF flag * float transform * First pass compiles * Update from mynn_changes from ewquon PR6 * Update from mynn_changes from ewquon PR8 * Style fixes * More style fixes * Fix workarray size * Add if guard for unused functions * codespell fixes * More if guards
1 parent ac3e595 commit c48de1a

23 files changed

+8777
-16
lines changed

‎.codespell-ignore-words

+3-1
Original file line numberDiff line numberDiff line change
@@ -10,4 +10,6 @@ ABL
1010
abl
1111
indx
1212
aer
13-
13+
delt
14+
dum
15+
adn

‎CMake/BuildERFExe.cmake

+1
Original file line numberDiff line numberDiff line change
@@ -166,6 +166,7 @@ function(build_erf_lib erf_lib_name)
166166
${SRC_DIR}/IO/ERF_WriteJobInfo.cpp
167167
${SRC_DIR}/IO/ERF_ConsoleIO.cpp
168168
${SRC_DIR}/PBL/ERF_ComputeDiffusivityMYNN25.cpp
169+
${SRC_DIR}/PBL/ERF_ComputeDiffusivityMYNNEDMF.cpp
169170
${SRC_DIR}/PBL/ERF_ComputeDiffusivityYSU.cpp
170171
${SRC_DIR}/SourceTerms/ERF_ApplySpongeZoneBCs.cpp
171172
${SRC_DIR}/SourceTerms/ERF_ApplySpongeZoneBCs_ReadFromFile.cpp

‎Exec/ABL/inputs_most_mynnedmf

+63
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,63 @@
1+
# ------------------ INPUTS TO MAIN PROGRAM -------------------
2+
max_step = 3000
3+
4+
amrex.fpe_trap_invalid = 1
5+
6+
fabarray.mfiter_tile_size = 1024 1024 1024
7+
8+
# PROBLEM SIZE & GEOMETRY
9+
geometry.prob_extent = 1024 1024 1024
10+
amr.n_cell = 4 4 64
11+
12+
#erf.grid_stretching_ratio = 1.025
13+
#erf.initial_dz = 16.0
14+
15+
geometry.is_periodic = 1 1 0
16+
17+
# MOST BOUNDARY (DEFAULT IS ADIABATIC FOR THETA)
18+
zlo.type = "Most"
19+
erf.most.z0 = 0.1
20+
erf.most.zref = 8.0
21+
22+
zhi.type = "SlipWall"
23+
24+
# TIME STEP CONTROL
25+
erf.fixed_dt = 0.1 # fixed time step depending on grid resolution
26+
27+
# DIAGNOSTICS & VERBOSITY
28+
erf.sum_interval = 1 # timesteps between computing mass
29+
erf.v = 1 # verbosity in ERF.cpp
30+
amr.v = 1 # verbosity in Amr.cpp
31+
32+
# REFINEMENT / REGRIDDING
33+
amr.max_level = 0 # maximum level number allowed
34+
35+
# CHECKPOINT FILES
36+
erf.check_file = chk # root name of checkpoint file
37+
erf.check_int = -1 # number of timesteps between checkpoints
38+
39+
# PLOTFILES
40+
erf.plot_file_1 = plt # prefix of plotfile name
41+
erf.plot_int_1 = 100 # number of timesteps between plotfiles
42+
erf.plot_vars_1 = density rhoadv_0 x_velocity y_velocity z_velocity pressure temp theta
43+
44+
# SOLVER CHOICE
45+
erf.use_gravity = false
46+
47+
erf.molec_diff_type = "None"
48+
erf.les_type = "Smagorinsky"
49+
erf.Cs = 0.1
50+
51+
erf.pbl_type = "MYNNEDMF"
52+
53+
erf.init_type = "uniform"
54+
55+
# PROBLEM PARAMETERS
56+
prob.rho_0 = 1.0
57+
prob.A_0 = 1.0
58+
59+
prob.U_0 = 10.0
60+
prob.V_0 = 0.0
61+
prob.W_0 = 0.0
62+
prob.T_0 = 300.0
63+

‎Exec/ABL/module_bl_mynn_1d_cc.cpp

+4,180
Large diffs are not rendered by default.

‎Exec/ABL/transform_text_sed_real.sh

+21
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,21 @@
1+
#!/bin/bash
2+
3+
sed -i s/a2fac/afac2/g module_bl_mynn_1d_cc.cpp
4+
max=10
5+
for ((i = 0 ; i < max ; i++ )); do
6+
sed -i 's/\([0-9]f\)/\1_rt/' module_bl_mynn_1d_cc.cpp
7+
sed -i s/f_rt/_rt/g module_bl_mynn_1d_cc.cpp
8+
sed -i 's/\(\.f\)/\1_rt/' module_bl_mynn_1d_cc.cpp
9+
sed -i s/f_rt/_rt/g module_bl_mynn_1d_cc.cpp
10+
sed -i s/'float\*'/'Real\*'/g module_bl_mynn_1d_cc.cpp
11+
sed -i s/'float('/'Real('/g module_bl_mynn_1d_cc.cpp
12+
sed -i s/'float '/'Real '/g module_bl_mynn_1d_cc.cpp
13+
sed -i s/'float\&'/'Real\&'/g module_bl_mynn_1d_cc.cpp
14+
sed -i s/'<float>'/'<Real>'/g module_bl_mynn_1d_cc.cpp
15+
sed -i s/' abs'/' std::abs'/g module_bl_mynn_1d_cc.cpp
16+
sed -i s/'=abs'/'=std::abs'/g module_bl_mynn_1d_cc.cpp
17+
sed -i s/'(abs'/'(std::abs'/g module_bl_mynn_1d_cc.cpp
18+
done
19+
sed -i s/afac2/a2fac/g module_bl_mynn_1d_cc.cpp
20+
#grep [0-9]f ~/codes/WRF_mynn_changes/phys/module_bl_mynn_1d_cc.cpp | sed 's/\([0-9]f\)/\1_rt/' | sed s/f_rt/_rt/g | sed 's/\([0-9]f\)/\1_rt/' | sed s/f_rt/_rt/g | sed 's/\([0-9]f\)/\1_rt/' | sed s/f_rt/_rt/g | sed 's/\([0-9]f\)/\1_rt/' | sed s/f_rt/_rt/g | sed 's/\([0-9]f\)/\1_rt/' | sed s/f_rt/_rt/g | sed 's/\([0-9]f\)/\1_rt/' | sed s/f_rt/_rt/g | sed 's/\([0-9]f\)/\1_rt/' | sed s/f_rt/_rt/g | sed 's/\([0-9]f\)/\1_rt/' | sed s/f_rt/_rt/g
21+

‎Source/BoundaryConditions/ERF_ABLMost.H

+2
Original file line numberDiff line numberDiff line change
@@ -79,6 +79,8 @@ public:
7979
pblh_type = PBLHeightCalcType::None;
8080
} else if (pblh_string == "MYNN25") {
8181
pblh_type = PBLHeightCalcType::MYNN25;
82+
} else if (pblh_string == "MYNNEDMF") {
83+
pblh_type = PBLHeightCalcType::MYNN25;
8284
} else if (pblh_string == "YSU") {
8385
pblh_type = PBLHeightCalcType::YSU;
8486
} else {

‎Source/DataStructs/ERF_TurbStruct.H

+6-4
Original file line numberDiff line numberDiff line change
@@ -7,7 +7,7 @@ AMREX_ENUM(LESType, None, Smagorinsky, Deardorff);
77

88
AMREX_ENUM(RANSType, None, kEqn);
99

10-
AMREX_ENUM(PBLType, None, MYNN25, YSU);
10+
AMREX_ENUM(PBLType, None, MYNN25, MYNNEDMF, YSU);
1111

1212
template <typename T>
1313
void query_one_or_per_level(const amrex::ParmParse& pp, const char* query_string, T& query_var, const int lev, const int maxlev)
@@ -60,7 +60,7 @@ struct TurbChoice {
6060
amrex::Error("It is not recommended to use Deardorff LES and a PBL model");
6161
}
6262

63-
if (pbl_type == PBLType::MYNN25) {
63+
if (pbl_type == PBLType::MYNN25 || pbl_type == PBLType::MYNNEDMF) {
6464
query_one_or_per_level(pp, "pbl_mynn_A1", pbl_mynn.A1, lev, max_level);
6565
query_one_or_per_level(pp, "pbl_mynn_A2", pbl_mynn.A2, lev, max_level);
6666
query_one_or_per_level(pp, "pbl_mynn_B1", pbl_mynn.B1, lev, max_level);
@@ -87,7 +87,7 @@ struct TurbChoice {
8787
}
8888

8989
// Right now, solving the QKE equation is only supported when MYNN PBL is turned on
90-
if (pbl_type == PBLType::MYNN25) {
90+
if (pbl_type == PBLType::MYNN25 || pbl_type == PBLType::MYNNEDMF) {
9191
use_KE = true;
9292
query_one_or_per_level(pp, "advect_KE" , advect_KE, lev, max_level);
9393
query_one_or_per_level(pp, "diffuse_KE_3D", diffuse_KE_3D, lev, max_level);
@@ -142,6 +142,8 @@ struct TurbChoice {
142142
amrex::Print() << " Using Axell & Liungman one-equation RANS k model at level " << lev << std::endl;
143143
} else if (pbl_type == PBLType::MYNN25) {
144144
amrex::Print() << " Using MYNN2.5 PBL model at level " << lev << std::endl;
145+
} else if (pbl_type == PBLType::MYNNEDMF) {
146+
amrex::Print() << " Using MYNNEDMF PBL model at level " << lev << std::endl;
145147
} else if (pbl_type == PBLType::YSU) {
146148
amrex::Print() << " Using YSU PBL model at level " << lev << std::endl;
147149
}
@@ -179,7 +181,7 @@ struct TurbChoice {
179181
amrex::Print() << "Sc_t : " << Sc_t << std::endl;
180182
}
181183

182-
if (pbl_type == PBLType::MYNN25) {
184+
if (pbl_type == PBLType::MYNN25 || pbl_type == PBLType::MYNNEDMF) {
183185
amrex::Print() << " pbl_mynn_A1 : " << pbl_mynn.A1 << std::endl;
184186
amrex::Print() << " pbl_mynn_A2 : " << pbl_mynn.A2 << std::endl;
185187
amrex::Print() << " pbl_mynn_B1 : " << pbl_mynn.B1 << std::endl;

‎Source/Diffusion/ERF_ComputeTurbulentViscosity.cpp

+7
Original file line numberDiff line numberDiff line change
@@ -579,6 +579,7 @@ void ComputeTurbulentViscosity (const MultiFab& xvel , const MultiFab& yvel ,
579579
turbChoice.les_type == LESType::Deardorff ||
580580
turbChoice.rans_type == RANSType::kEqn ||
581581
turbChoice.pbl_type == PBLType::MYNN25 ||
582+
turbChoice.pbl_type == PBLType::MYNNEDMF ||
582583
turbChoice.pbl_type == PBLType::YSU );
583584
AMREX_ALWAYS_ASSERT_WITH_MESSAGE(l_use_turb,
584585
"A turbulence model must be utilized with MOST boundaries to compute the turbulent viscosity");
@@ -617,6 +618,12 @@ void ComputeTurbulentViscosity (const MultiFab& xvel , const MultiFab& yvel ,
617618
level, bc_ptr, vert_only, z_phys_nd,
618619
solverChoice.RhoQv_comp, solverChoice.RhoQc_comp,
619620
solverChoice.RhoQr_comp);
621+
} else if (turbChoice.pbl_type == PBLType::MYNNEDMF) {
622+
ComputeDiffusivityMYNNEDMF(xvel, yvel, cons_in, eddyViscosity,
623+
geom, turbChoice, most, use_moisture,
624+
level, bc_ptr, vert_only, z_phys_nd,
625+
solverChoice.RhoQv_comp, solverChoice.RhoQc_comp,
626+
solverChoice.RhoQr_comp);
620627
} else if (turbChoice.pbl_type == PBLType::YSU) {
621628
ComputeDiffusivityYSU(xvel, yvel, cons_in, eddyViscosity,
622629
geom, turbChoice, most, use_moisture,

‎Source/Diffusion/ERF_DiffusionSrcForState_N.cpp

+2-1
Original file line numberDiff line numberDiff line change
@@ -84,13 +84,14 @@ DiffusionSrcForState_N (const Box& bx, const Box& domain,
8484

8585
bool l_use_keqn = ( (turbChoice.les_type == LESType::Deardorff) ||
8686
(turbChoice.rans_type == RANSType::kEqn) );
87-
bool l_use_mynn = (turbChoice.pbl_type == PBLType::MYNN25);
87+
bool l_use_mynn = ( (turbChoice.pbl_type == PBLType::MYNN25) || (turbChoice.pbl_type == PBLType::MYNNEDMF) ) ;
8888

8989
bool l_consA = (diffChoice.molec_diff_type == MolecDiffType::ConstantAlpha);
9090
bool l_turb = ( (turbChoice.les_type == LESType::Smagorinsky) ||
9191
(turbChoice.les_type == LESType::Deardorff ) ||
9292
(turbChoice.rans_type == RANSType::kEqn ) ||
9393
(turbChoice.pbl_type == PBLType::MYNN25 ) ||
94+
(turbChoice.pbl_type == PBLType::MYNNEDMF ) ||
9495
(turbChoice.pbl_type == PBLType::YSU ) );
9596

9697
const Box xbx = surroundingNodes(bx,0);

‎Source/Diffusion/ERF_DiffusionSrcForState_T.cpp

+2-1
Original file line numberDiff line numberDiff line change
@@ -96,13 +96,14 @@ DiffusionSrcForState_T (const Box& bx, const Box& domain,
9696

9797
bool l_use_keqn = ( (turbChoice.les_type == LESType::Deardorff) ||
9898
(turbChoice.rans_type == RANSType::kEqn) );
99-
bool l_use_mynn = (turbChoice.pbl_type == PBLType::MYNN25);
99+
bool l_use_mynn = ( (turbChoice.pbl_type == PBLType::MYNN25) || (turbChoice.pbl_type == PBLType::MYNNEDMF) ) ;
100100

101101
bool l_consA = (diffChoice.molec_diff_type == MolecDiffType::ConstantAlpha);
102102
bool l_turb = ( (turbChoice.les_type == LESType::Smagorinsky) ||
103103
(turbChoice.les_type == LESType::Deardorff ) ||
104104
(turbChoice.rans_type == RANSType::kEqn ) ||
105105
(turbChoice.pbl_type == PBLType::MYNN25 ) ||
106+
(turbChoice.pbl_type == PBLType::MYNNEDMF ) ||
106107
(turbChoice.pbl_type == PBLType::YSU ) );
107108

108109
const Box xbx = surroundingNodes(bx,0);

‎Source/ERF.cpp

+2-1
Original file line numberDiff line numberDiff line change
@@ -653,9 +653,10 @@ ERF::InitData_pre ()
653653
// Verify BCs are compatible with solver choice
654654
for (int lev(0); lev <= max_level; ++lev) {
655655
if ( ( (solverChoice.turbChoice[lev].pbl_type == PBLType::MYNN25) ||
656+
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNNEDMF) ||
656657
(solverChoice.turbChoice[lev].pbl_type == PBLType::YSU) ) &&
657658
phys_bc_type[Orientation(Direction::z,Orientation::low)] != ERF_BC::MOST ) {
658-
Abort("MYNN2.5/YSU PBL Model requires MOST at lower boundary");
659+
Abort("MYNN2.5/MYNNEDMF/YSU PBL Model requires MOST at lower boundary");
659660
}
660661

661662
if ( (solverChoice.turbChoice[lev].les_type == LESType::Deardorff) &&

‎Source/IO/ERF_Checkpoint.cpp

+2-1
Original file line numberDiff line numberDiff line change
@@ -408,7 +408,8 @@ ERF::ReadCheckpointFile ()
408408
MultiFab::Copy(vars_new[lev][Vars::cons],cons,0,0,(RhoKE_comp+1),0);
409409

410410
// Only if we have a PBL model do we need to copy QKE is src to KE in dst
411-
if (solverChoice.turbChoice[lev].pbl_type == PBLType::MYNN25) {
411+
if ( (solverChoice.turbChoice[lev].pbl_type == PBLType::MYNN25) ||
412+
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNNEDMF) ) {
412413
MultiFab::Copy(vars_new[lev][Vars::cons],cons,(RhoKE_comp+1),RhoKE_comp,1,0);
413414
vars_new[lev][Vars::cons].mult(0.5,RhoKE_comp,1,0);
414415
}

‎Source/IO/ERF_Write1DProfiles.cpp

+2-2
Original file line numberDiff line numberDiff line change
@@ -213,8 +213,8 @@ ERF::derive_diag_profiles(Real /*time*/,
213213
(solverChoice.turbChoice[lev].pbl_type != PBLType::None));
214214
bool l_use_KE = ( (solverChoice.turbChoice[lev].les_type == LESType::Deardorff) ||
215215
(solverChoice.turbChoice[lev].rans_type == RANSType::kEqn) ||
216-
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNN25) );
217-
216+
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNN25) ||
217+
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNNEDMF) );
218218
// This will hold rho, theta, ksgs, Kmh, Kmv, uu, uv, uw, vv, vw, ww, uth, vth, wth,
219219
// 0 1 2 3 4 5 6 7 8 9 10 11 12 13
220220
// thth, uiuiu, uiuiv, uiuiw, p, pu, pv, pw, qv, qc, qr, wqv, wqc, wqr,

‎Source/IO/ERF_Write1DProfiles_stag.cpp

+2-2
Original file line numberDiff line numberDiff line change
@@ -321,8 +321,8 @@ ERF::derive_diag_profiles_stag (Real /*time*/,
321321
(solverChoice.turbChoice[lev].pbl_type != PBLType::None));
322322
bool l_use_KE = ( (solverChoice.turbChoice[lev].les_type == LESType::Deardorff) ||
323323
(solverChoice.turbChoice[lev].rans_type == RANSType::kEqn) ||
324-
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNN25) );
325-
324+
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNN25) ||
325+
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNNEDMF) );
326326
// Note: "uiui" == u_i*u_i = u*u + v*v + w*w
327327
// This will hold rho, theta, ksgs, Kmh, Kmv, uu, uv, vv, uth, vth,
328328
// indices: 0 1 2 3 4 5 6 7 8 9

‎Source/Initialization/ERF_InitCustom.cpp

+2-1
Original file line numberDiff line numberDiff line change
@@ -85,7 +85,8 @@ ERF::init_custom (int lev)
8585
// RhoKE is relevant if using Deardorff with LES, k-equation for RANS, or MYNN with PBL
8686
if ((solverChoice.turbChoice[lev].les_type == LESType::Deardorff) ||
8787
(solverChoice.turbChoice[lev].rans_type == RANSType::kEqn) ||
88-
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNN25)) {
88+
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNN25) ||
89+
(solverChoice.turbChoice[lev].pbl_type == PBLType::MYNNEDMF) ) {
8990
MultiFab::Add(lev_new[Vars::cons], cons_pert, RhoKE_comp, RhoKE_comp, 1, cons_pert.nGrow());
9091
}
9192

‎Source/PBL/ERF_ComputeDiffusivityMYNN25.cpp

+1
Original file line numberDiff line numberDiff line change
@@ -24,6 +24,7 @@ ComputeDiffusivityMYNN25 (const MultiFab& xvel,
2424
const int RhoQc_comp,
2525
const int RhoQr_comp)
2626
{
27+
2728
const bool use_terrain = (z_phys_nd != nullptr);
2829
const bool use_most = (most != nullptr);
2930

0 commit comments

Comments
 (0)