From 305f977bf47dba57a16f9a34aaaff1a96d513436 Mon Sep 17 00:00:00 2001
From: Zhaoyi Shen <11598433+szy21@users.noreply.github.com>
Date: Wed, 20 Nov 2024 11:14:48 -0800
Subject: [PATCH] modify constant area detrainment

---
 .../prognostic_edmfx_bomex_column.yml           |  6 +++---
 .../prognostic_edmfx_simpleplume_column.yml     | 11 +++++------
 .../diagnostic_edmf_precomputed_quantities.jl   |  1 +
 .../prognostic_edmf_precomputed_quantities.jl   |  3 +++
 src/prognostic_equations/edmfx_entr_detr.jl     | 17 +++++++++++++++--
 src/surface_conditions/surface_setups.jl        |  2 +-
 toml/prognostic_edmfx_bomex.toml                |  4 ++--
 7 files changed, 30 insertions(+), 14 deletions(-)

diff --git a/config/model_configs/prognostic_edmfx_bomex_column.yml b/config/model_configs/prognostic_edmfx_bomex_column.yml
index a495d58f85..ea41aea576 100644
--- a/config/model_configs/prognostic_edmfx_bomex_column.yml
+++ b/config/model_configs/prognostic_edmfx_bomex_column.yml
@@ -4,17 +4,17 @@ edmf_coriolis: "Bomex"
 ls_adv: "Bomex"
 surface_setup: "Bomex"
 turbconv: "prognostic_edmfx"
-implicit_diffusion: true
+implicit_diffusion: false
 implicit_sgs_advection: false
 approximate_linear_solve_iters: 2
 edmfx_upwinding: first_order
 edmfx_entr_model: "Generalized"
-edmfx_detr_model: "Generalized"
+edmfx_detr_model: "ConstantArea"
 edmfx_sgs_mass_flux: true
 edmfx_sgs_diffusive_flux: true
 edmfx_nh_pressure: true
 edmfx_filter: true
-prognostic_tke: true
+prognostic_tke: false
 moist: "equil"
 cloud_model: "quadrature_sgs"
 call_cloud_diagnostics_per_stage: true
diff --git a/config/model_configs/prognostic_edmfx_simpleplume_column.yml b/config/model_configs/prognostic_edmfx_simpleplume_column.yml
index 2c21f5f097..f0cc6bffec 100644
--- a/config/model_configs/prognostic_edmfx_simpleplume_column.yml
+++ b/config/model_configs/prognostic_edmfx_simpleplume_column.yml
@@ -1,14 +1,13 @@
 initial_condition: "SimplePlume"
 surface_setup: "SimplePlume"
 turbconv: "prognostic_edmfx"
-implicit_diffusion: true
-implicit_sgs_advection: true
+implicit_diffusion: false
+implicit_sgs_advection: false
 approximate_linear_solve_iters: 2
 max_newton_iters_ode: 3
-zero_tendency: "grid_scale"
 edmfx_upwinding: first_order
 edmfx_entr_model: "Generalized"
-edmfx_detr_model: "Generalized"
+edmfx_detr_model: "ConstantArea"
 edmfx_sgs_mass_flux: false
 edmfx_sgs_diffusive_flux: false
 edmfx_nh_pressure: true
@@ -19,10 +18,10 @@ cloud_model: "quadrature_sgs"
 call_cloud_diagnostics_per_stage: true
 config: "column"
 z_max: 4e3
-z_elem: 80
+z_elem: 40
 z_stretch: false
 perturb_initstate: false
-dt: 20secs"
+dt: 30secs"
 t_end: "12hours"
 dt_save_state_to_disk: "10mins"
 toml: [toml/prognostic_edmfx_simpleplume.toml]
diff --git a/src/cache/diagnostic_edmf_precomputed_quantities.jl b/src/cache/diagnostic_edmf_precomputed_quantities.jl
index fbf43ecf2a..582b525f9e 100644
--- a/src/cache/diagnostic_edmf_precomputed_quantities.jl
+++ b/src/cache/diagnostic_edmf_precomputed_quantities.jl
@@ -617,6 +617,7 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_do_integral!(
                 entrʲ_prev_level,
                 vert_div_level,
                 FT(0), # mass flux divergence is not implemented for diagnostic edmf
+                FT(0), # w divergence is not implemented for diagnostic edmf
                 tke_prev_level,
                 p.atmos.edmfx_model.detr_model,
             )
diff --git a/src/cache/prognostic_edmf_precomputed_quantities.jl b/src/cache/prognostic_edmf_precomputed_quantities.jl
index 6e41f0f2e5..efb00326a4 100644
--- a/src/cache/prognostic_edmf_precomputed_quantities.jl
+++ b/src/cache/prognostic_edmf_precomputed_quantities.jl
@@ -215,6 +215,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
 
     ᶜvert_div = p.scratch.ᶜtemp_scalar
     ᶜmassflux_vert_div = p.scratch.ᶜtemp_scalar_2
+    ᶜw_vert_div = p.scratch.ᶜtemp_scalar_3
     for j in 1:n
         # entrainment/detrainment
         @. ᶜentrʲs.:($$j) = entrainment(
@@ -252,6 +253,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
         @. ᶜvert_div = ᶜdivᵥ(ᶠinterp(ᶜρʲs.:($$j)) * ᶠu³ʲs.:($$j)) / ᶜρʲs.:($$j)
         @. ᶜmassflux_vert_div =
             ᶜdivᵥ(ᶠinterp(Y.c.sgsʲs.:($$j).ρa) * ᶠu³ʲs.:($$j))
+        @. ᶜw_vert_div = ᶜdivᵥ(ᶠu³ʲs.:($$j))
         @. ᶜdetrʲs.:($$j) = detrainment(
             thermo_params,
             turbconv_params,
@@ -270,6 +272,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
             ᶜentrʲs.:($$j),
             ᶜvert_div,
             ᶜmassflux_vert_div,
+            ᶜw_vert_div,
             ᶜtke⁰,
             p.atmos.edmfx_model.detr_model,
         )
diff --git a/src/prognostic_equations/edmfx_entr_detr.jl b/src/prognostic_equations/edmfx_entr_detr.jl
index 092761c648..f48e20ef56 100644
--- a/src/prognostic_equations/edmfx_entr_detr.jl
+++ b/src/prognostic_equations/edmfx_entr_detr.jl
@@ -142,6 +142,7 @@ function detrainment_from_thermo_state(
     entrʲ_prev_level,
     vert_div_level,
     ᶜmassflux_vert_div, # mass flux divergence is not implemented for diagnostic edmf
+    ᶜw_vert_div, # mass flux divergence is not implemented for diagnostic edmf
     tke_prev_level,
     edmfx_detr_model,
 )
@@ -164,6 +165,7 @@ function detrainment_from_thermo_state(
         entrʲ_prev_level,
         vert_div_level,
         FT(0), # mass flux divergence is not implemented for diagnostic edmf
+        FT(0), # w divergence is not implemented for diagnostic edmf
         tke_prev_level,
         edmfx_detr_model,
     )
@@ -202,6 +204,7 @@ function detrainment(
     ᶜentr,
     ᶜvert_div,
     ᶜmassflux_vert_div,
+    ᶜw_vert_div,
     ᶜtke⁰,
     ::NoDetrainment,
 )
@@ -226,6 +229,7 @@ function detrainment(
     ᶜentr,
     ᶜvert_div,
     ᶜmassflux_vert_div,
+    ᶜw_vert_div,
     ᶜtke⁰,
     ::PiGroupsDetrainment,
 )
@@ -278,6 +282,7 @@ function detrainment(
     ᶜentr,
     ᶜvert_div,
     ᶜmassflux_vert_div,
+    ᶜw_vert_div,
     ᶜtke⁰,
     ::GeneralizedDetrainment,
 )
@@ -329,9 +334,15 @@ function detrainment(
     ᶜentr,
     ᶜvert_div,
     ᶜmassflux_vert_div,
+    ᶜw_vert_div,
+    ᶜtke⁰,
     ::ConstantAreaDetrainment,
 )
-    detr = ᶜentr - ᶜvert_div
+    if ᶜρaʲ <= 0
+        detr = 0
+    else
+        detr = (ᶜw_vert_div < 0) * (ᶜentr - ᶜw_vert_div)
+    end
     return max(detr, 0)
 end
 
@@ -374,8 +385,10 @@ limit_entrainment(entr::FT, a, dt) where {FT} = max(
     min(entr, FT(0.9) * (1 - a) / max(a, eps(FT)) / dt, FT(0.9) * 1 / dt),
     0,
 )
+# limit_detrainment(detr::FT, a, dt) where {FT} =
+#     max(min(detr, FT(0.9) * 1 / dt), 0)
 limit_detrainment(detr::FT, a, dt) where {FT} =
-    max(min(detr, FT(0.9) * 1 / dt), 0)
+    max(detr, 0)
 
 function limit_turb_entrainment(dyn_entr::FT, turb_entr, dt) where {FT}
     return max(min((FT(0.9) * 1 / dt) - dyn_entr, turb_entr), 0)
diff --git a/src/surface_conditions/surface_setups.jl b/src/surface_conditions/surface_setups.jl
index 3ad6d88162..0bb451130b 100644
--- a/src/surface_conditions/surface_setups.jl
+++ b/src/surface_conditions/surface_setups.jl
@@ -247,7 +247,7 @@ function (::SimplePlume)(params)
     T = FT(310)
     p = FT(101500)
     q_vap = FT(0.02245)
-    θ_flux = FT(8)
+    θ_flux = FT(8e-2)
     q_flux = FT(0)
     z0 = FT(1e-4)
     ustar = FT(0.28)
diff --git a/toml/prognostic_edmfx_bomex.toml b/toml/prognostic_edmfx_bomex.toml
index 3d6a8c492c..677a221502 100644
--- a/toml/prognostic_edmfx_bomex.toml
+++ b/toml/prognostic_edmfx_bomex.toml
@@ -8,10 +8,10 @@ value = 1.0e-5
 value = 0.7
 
 [entr_inv_tau]
-value = 0.0001
+value = 0.0002
 
 [entr_coeff]
-value = 0.3
+value = 0
 
 [min_area_limiter_scale]
 value = 0