limit area fraction in entrainment and detrainment limiter

CliMA · Oct 25, 2023 · ae3b983 · ae3b983
1 parent 768aeb0
commit ae3b983
Show file tree

Hide file tree

Showing 2 changed files with 9 additions and 6 deletions.
diff --git a/config/model_configs/diagnostic_edmfx_dycoms_rf01_box.yml b/config/model_configs/diagnostic_edmfx_dycoms_rf01_box.yml
@@ -23,7 +23,7 @@ y_elem: 2
 z_elem: 30 
 z_max: 1500 
 z_stretch: false
-dt: 6secs 
+dt: 10secs 
 t_end: 4hours 
 dt_save_to_disk: 2mins
 toml: [toml/diagnostic_edmfx_box.toml]
diff --git a/src/prognostic_equations/edmfx_entr_detr.jl b/src/prognostic_equations/edmfx_entr_detr.jl
@@ -112,7 +112,8 @@ function entrainment(
     turbconv_params = CAP.turbconv_params(params)
     a_min = TCP.min_area(turbconv_params)
     min_area_limiter =
-        min_area_limiter_scale * exp(-min_area_limiter_power * (ᶜaʲ - a_min))
+        min_area_limiter_scale *
+        exp(-min_area_limiter_power * (max(ᶜaʲ, 0) - a_min))
     entr = min(
         entr_inv_tau + entr_coeff * abs(ᶜwʲ) / (ᶜz - z_sfc) + min_area_limiter,
         1 / dt,
@@ -144,7 +145,8 @@ function entrainment(
     turbconv_params = CAP.turbconv_params(params)
     a_min = TCP.min_area(turbconv_params)
     min_area_limiter =
-        min_area_limiter_scale * exp(-min_area_limiter_power * (ᶜaʲ - a_min))
+        min_area_limiter_scale *
+        exp(-min_area_limiter_power * (max(ᶜaʲ, 0) - a_min))
     entr = min(
         entr_inv_tau + entr_coeff * abs(ᶜwʲ) / (ᶜz - z_sfc) + min_area_limiter,
         1 / dt,
@@ -262,12 +264,13 @@ function detrainment(
     turbconv_params = CAP.turbconv_params(params)
     a_max = TCP.max_area(turbconv_params)
     max_area_limiter =
-        max_area_limiter_scale * exp(-max_area_limiter_power * (a_max - ᶜaʲ))
+        max_area_limiter_scale *
+        exp(-max_area_limiter_power * (a_max - min(ᶜaʲ, 1)))
     detr = min(
         detr_inv_tau +
         detr_coeff * abs(ᶜwʲ) +
         detr_buoy_coeff * abs(min(ᶜbuoyʲ - ᶜbuoy⁰, 0)) /
-        max(eps(FT), abs((ᶜwʲ - ᶜw⁰))) +
+        max(eps(FT), abs(ᶜwʲ - ᶜw⁰)) +
         max_area_limiter,
         1 / dt,
     )
@@ -304,7 +307,7 @@ function detrainment(
         detr_inv_tau +
         detr_coeff * abs(ᶜwʲ) +
         detr_buoy_coeff * abs(min(ᶜbuoyʲ - ᶜbuoy⁰, 0)) /
-        max(eps(FT), abs((ᶜwʲ - ᶜw⁰))) +
+        max(eps(FT), abs(ᶜwʲ - ᶜw⁰)) +
         max_area_limiter,
         1 / dt,
     )