Merge pull request #4029 from CliMA/sa/fix_vanishing_ar_anomalies

Fix updraft values for vanishing ar

Author: sajjadazimi

config/model_configs/prognostic_edmfx_rico_column.yml

@@ -27,7 +27,7 @@ z_elem: 100
 z_stretch: false
 perturb_initstate: false
 dt: "5secs"
-t_end: "12hours"
+t_end: "24hours"
 dt_save_state_to_disk: "60mins"
 toml: [toml/prognostic_edmfx_1M.toml]
 netcdf_interpolation_num_points: [8, 8, 100]

reproducibility_tests/ref_counter.jl

@@ -1,4 +1,4 @@
-263
+264
 
 # **README**
 #
@@ -20,6 +20,9 @@
 
 
 #=
+264
+- Allow entrainment and detrainment when updraft area fraction is negligible
+
 263
 - Make vertical diffusion of updrafts implicit
 

src/prognostic_equations/edmfx_entr_detr.jl

@@ -139,9 +139,11 @@ function entrainment(
     ::PiGroupsEntrainment,
 )
     FT = eltype(thermo_params)
-    # Entrainment is not well-defined or should be zero if updraft area is negligible.
+    entr_inv_tau = CAP.entr_inv_tau(turbconv_params)
+    # Entrainment is not well-defined if updraft area is negligible.
+    # Fix at entr_inv_tau to ensure some mixing with the environment.
     if ᶜaʲ <= eps(FT)
-        return 0
+        return entr_inv_tau
     end
 
     elev_above_sfc = ᶜz - z_sfc
@@ -168,7 +170,6 @@ function entrainment(
     )
 
     entr_param_vec = CAP.entr_param_vec(turbconv_params)
-    entr_inv_tau = CAP.entr_inv_tau(turbconv_params)
     pi_sum =
         entr_param_vec[1] * abs(Π₁) +
         entr_param_vec[2] * abs(Π₂) +
@@ -212,10 +213,10 @@ function entrainment(
     min_area_limiter_power = CAP.min_area_limiter_power(turbconv_params)
     a_min = CAP.min_area(turbconv_params)
 
-    # Entrainment is not well-defined or should be zero if updraft area is negligible,
+    # Entrainment is not well-defined if updraft area is negligible,
     # as some limiters depend on ᶜaʲ.
     if ᶜaʲ <= eps(FT) && min_area_limiter_scale == FT(0) # If no area and no base min_area_limiter
-        return 0
+        return entr_inv_tau
     end
 
     min_area_limiter =
@@ -386,12 +387,13 @@ function detrainment(
     ::PiGroupsDetrainment,
 )
     FT = eltype(thermo_params)
+    detr_inv_tau = CAP.detr_inv_tau(turbconv_params)
 
-    # If ᶜρaʲ (updraft effective density) is zero or negligible,
-    # detrainment is considered zero. This also protects division by ᶜρaʲ later.
+    # If ᶜρaʲ (updraft effective density) is zero or negligible, detrainment is 
+    # considered to be fixed at detr_inv_tau. This also protects division by ᶜρaʲ later.
     # This condition implies the updraft area (ᶜaʲ) is also likely negligible.
     if ᶜρaʲ <= eps(FT)
-        return 0
+        return detr_inv_tau
     end
 
     elev_above_sfc = ᶜz - z_sfc
@@ -466,9 +468,10 @@ function detrainment(
     max_area_limiter_power = CAP.max_area_limiter_power(turbconv_params)
     a_max = CAP.max_area(turbconv_params)
 
-    # If ᶜρaʲ (updraft effective density) is zero or negligible, detrainment is zero.
-    if ᶜρaʲ <= eps(FT) # Consistent check
-        return 0
+    # If ᶜρaʲ (updraft effective density) is zero or negligible, detrainment is not well defined.
+    # Fix at detr_inv_tau to ensure some mixing with the environment.
+    if ᶜρaʲ <= eps(FT)
+        return detr_inv_tau
     end
 
     max_area_limiter =
@@ -508,8 +511,12 @@ function detrainment(
     ::SmoothAreaDetrainment,
 )
     FT = eltype(thermo_params)
-    # If ᶜρaʲ is negligible or vertical velocity divergence term is non-negative, detrainment is zero.
-    if (ᶜρaʲ <= eps(FT)) || (ᶜw_vert_div >= 0) # Consistent check for ᶜρaʲ
+    detr_inv_tau = CAP.detr_inv_tau(turbconv_params)
+    # If ᶜρaʲ is negligible detrainment is fixed at detr_inv_tau.
+    if (ᶜρaʲ <= eps(FT)) # Consistent check for ᶜρaʲ
+        detr = detr_inv_tau
+        # If vertical velocity divergence term is non-negative detrainment is zero.
+    elseif (ᶜw_vert_div >= 0)
         detr = FT(0)
     else
         detr = ᶜentr - ᶜw_vert_div