Cleanup and return a few options for MYNN25 with TropicalCyclone.

erf-model · Oct 24, 2024 · 4a8b2ff · 4a8b2ff
1 parent 50f416d
commit 4a8b2ff
Show file tree

Hide file tree

Showing 5 changed files with 27 additions and 23 deletions.
diff --git a/Source/DataStructs/ERF_TurbStruct.H b/Source/DataStructs/ERF_TurbStruct.H
@@ -76,6 +76,13 @@ struct TurbChoice {
             }
         }
 
+        // Right now, solving the QKE equation is only supported when MYNN PBL is turned on
+        if (pbl_type == PBLType::MYNN25) {
+            use_KE = true;
+            query_one_or_per_level(pp, "advect_KE"    , advect_KE, lev, max_level);
+            query_one_or_per_level(pp, "diffuse_KE_3D", diffuse_KE_3D, lev, max_level);
+        }
+
         // LES constants...
         query_one_or_per_level(pp, "Cs"  ,Cs, lev, max_level);
         query_one_or_per_level(pp, "CI"  ,CI, lev, max_level);
@@ -194,5 +201,10 @@ struct TurbChoice {
     bool        pbl_ysu_force_over_water = false; // Force YSU to act as if it is over water regardless of other inputs (for testing)
     amrex::Real pbl_ysu_land_Ribcr = 0.25; // Critical Bulk Richardson number of Land for stable conditions
     amrex::Real pbl_ysu_unst_Ribcr = 0.0;  // Critical Bulk Richardson number for unstable conditions
+
+    // QKE stuff - default is to use it, if MYNN2.5 PBL is used default is turb transport in Z-direction only
+    bool use_KE = true;
+    bool diffuse_KE_3D = true;
+    bool advect_KE = true;
 };
 #endif
diff --git a/Source/TimeIntegration/ERF_TI_slow_headers.H b/Source/TimeIntegration/ERF_TI_slow_headers.H
@@ -122,8 +122,6 @@ void erf_slow_rhs_post (int level, int finest_level, int nrk,
                        amrex::Vector<amrex::MultiFab>& S_data,
                        const amrex::MultiFab& S_prim,
                        amrex::Vector<amrex::MultiFab >& S_scratch,
-                       const amrex::MultiFab& xvel,
-                       const amrex::MultiFab& yvel,
                        const amrex::MultiFab& zvel,
                        const amrex::MultiFab& source,
                        const amrex::MultiFab* SmnSmn,

diff --git a/Source/TimeIntegration/ERF_TI_slow_rhs_fun.H b/Source/TimeIntegration/ERF_TI_slow_rhs_fun.H
@@ -342,8 +342,7 @@
         if ( solverChoice.use_terrain && (solverChoice.terrain_type == TerrainType::Moving) ) {
             erf_slow_rhs_post(level, finest_level, nrk, slow_dt, n_qstate,
                               S_rhs, S_old, S_new, S_data, S_prim, S_scratch,
-                              xvel_new, yvel_new, zvel_new,
-                              cc_src, SmnSmn, eddyDiffs,
+                              zvel_new, cc_src, SmnSmn, eddyDiffs,
                               Hfx1, Hfx2, Hfx3, Q1fx1, Q1fx2, Q1fx3, Q2fx3, Diss,
                               fine_geom, solverChoice, m_most, domain_bcs_type_d, domain_bcs_type,
                               z_phys_nd[level], ax[level], ay[level], az[level], detJ_cc[level], detJ_cc_new[level],
@@ -360,8 +359,7 @@
         } else {
             erf_slow_rhs_post(level, finest_level, nrk, slow_dt, n_qstate,
                               S_rhs, S_old, S_new, S_data, S_prim, S_scratch,
-                              xvel_new, yvel_new, zvel_new,
-                              cc_src, SmnSmn, eddyDiffs,
+                              zvel_new, cc_src, SmnSmn, eddyDiffs,
                               Hfx1, Hfx2, Hfx3, Q1fx1, Q1fx2, Q1fx3, Q2fx3, Diss,
                               fine_geom, solverChoice, m_most, domain_bcs_type_d, domain_bcs_type,
                               z_phys_nd[level], ax[level], ay[level], az[level], detJ_cc[level], detJ_cc[level],

diff --git a/Source/TimeIntegration/ERF_slow_rhs_post.cpp b/Source/TimeIntegration/ERF_slow_rhs_post.cpp
@@ -56,8 +56,6 @@ void erf_slow_rhs_post (int level, int finest_level,
                         Vector<MultiFab>& S_data,
                         const MultiFab& S_prim,
                         Vector<MultiFab>& S_scratch,
-                        const MultiFab& xvel,
-                        const MultiFab& yvel,
                         const MultiFab& /*zvel*/,
                         const MultiFab& source,
                         const MultiFab* SmnSmn,
@@ -111,9 +109,6 @@ void erf_slow_rhs_post (int level, int finest_level,
     DiffChoice dc = solverChoice.diffChoice;
     TurbChoice tc = solverChoice.turbChoice[level];
 
-    const MultiFab* t_mean_mf = nullptr;
-    if (most) t_mean_mf = most->get_mac_avg(0,2);
-
     const bool l_use_terrain      = solverChoice.use_terrain;
     const bool l_reflux = (solverChoice.coupling_type != CouplingType::OneWay);
     const bool l_moving_terrain   = (solverChoice.terrain_type == TerrainType::Moving);
@@ -122,6 +117,7 @@ void erf_slow_rhs_post (int level, int finest_level,
     const bool l_use_mono_adv   = solverChoice.use_mono_adv;
     const bool l_use_KE         = ( (tc.les_type == LESType::Deardorff) ||
                                     (tc.pbl_type == PBLType::MYNN25) );
+    const bool l_advect_KE      = (tc.use_KE && tc.advect_KE);
     const bool l_use_diff       = ((dc.molec_diff_type != MolecDiffType::None) ||
                                    (tc.les_type        !=       LESType::None) ||
                                    (tc.pbl_type        !=       PBLType::None) );
@@ -266,9 +262,6 @@ void erf_slow_rhs_post (int level, int finest_level,
         Array4<Real> avg_ymom = S_scratch[IntVars::ymom].array(mfi);
         Array4<Real> avg_zmom = S_scratch[IntVars::zmom].array(mfi);
 
-        const Array4<const Real> & u = xvel.array(mfi);
-        const Array4<const Real> & v = yvel.array(mfi);
-
         const Array4<Real const>& mu_turb = l_use_turb ? eddyDiffs->const_array(mfi) : Array4<const Real>{};
 
         const Array4<const Real>& z_nd         = l_use_terrain    ? z_phys_nd->const_array(mfi) : Array4<const Real>{};
@@ -393,11 +386,19 @@ void erf_slow_rhs_post (int level, int finest_level,
                     num_comp = 1;
                 }
 
-                if (l_use_diff) {
-
-                    const Array4<const Real> tm_arr = t_mean_mf ? t_mean_mf->const_array(mfi) :
-                                                                  Array4<const Real>{};
+                if (( ivar != RhoKE_comp                 ) ||
+                    ((ivar == RhoKE_comp) && l_advect_KE))
+                {
+                    AdvectionSrcForScalars(dt, tbx, start_comp, num_comp, avg_xmom, avg_ymom, avg_zmom,
+                                           cur_cons, cur_prim, cell_rhs,
+                                           l_use_mono_adv, max_s_ptr, min_s_ptr,
+                                           detJ_arr, dxInv, mf_m,
+                                           horiz_adv_type, vert_adv_type,
+                                           horiz_upw_frac, vert_upw_frac,
+                                           flx_arr, flx_tmp_arr, domain, bc_ptr_h);
+                }
 
+                if (l_use_diff) {
                     if (l_use_terrain) {
                         DiffusionSrcForState_T(tbx, domain, start_comp, num_comp, exp_most, rot_most,
                                                new_cons, cur_prim, cell_rhs,

diff --git a/Source/TimeIntegration/ERF_slow_rhs_pre.cpp b/Source/TimeIntegration/ERF_slow_rhs_pre.cpp
@@ -129,9 +129,6 @@ void erf_slow_rhs_pre (int level, int finest_level,
     DiffChoice dc = solverChoice.diffChoice;
     TurbChoice tc = solverChoice.turbChoice[level];
 
-    const MultiFab* t_mean_mf = nullptr;
-    if (most) t_mean_mf = most->get_mac_avg(0,2);
-
     int start_comp = 0;
     int   num_comp = 2;
     int   end_comp = start_comp + num_comp - 1;
@@ -492,8 +489,6 @@ void erf_slow_rhs_pre (int level, int finest_level,
             Array4<Real> q2fx_z = (Q2fx3) ? Q2fx3->array(mfi) : Array4<Real>{};
             Array4<Real> diss  = Diss->array(mfi);
 
-            const Array4<const Real> tm_arr = t_mean_mf ? t_mean_mf->const_array(mfi) : Array4<const Real>{};
-
             // NOTE: No diffusion for continuity, so n starts at 1.
             int n_start = amrex::max(start_comp,RhoTheta_comp);
             int n_comp  = end_comp - n_start + 1;