Merge pull request #3445 from CliMA/ck/detrainment_param_mem2

Pass subset of ClimaAtmosParameters to detrainment
CliMA · Nov 22, 2024 · 5f3da9f · 5f3da9f
2 parents fd0dbb6 + d9bbbc1
commit 5f3da9f
Show file tree

Hide file tree

Showing 4 changed files with 56 additions and 39 deletions.
diff --git a/src/cache/diagnostic_edmf_precomputed_quantities.jl b/src/cache/diagnostic_edmf_precomputed_quantities.jl
@@ -325,6 +325,7 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_do_integral!(
 
     thermo_params = CAP.thermodynamics_params(params)
     microphys_params = CAP.microphysics_precipitation_params(params)
+    turbconv_params = CAP.turbconv_params(params)
 
     ᶠΦ = p.scratch.ᶠtemp_scalar
     @. ᶠΦ = CAP.grav(params) * ᶠz
@@ -458,7 +459,8 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_do_integral!(
             tke_prev_level = Fields.field_values(Fields.level(ᶜtke⁰, i - 1))
 
             @. entrʲ_prev_level = entrainment(
-                params,
+                thermo_params,
+                turbconv_params,
                 z_prev_level,
                 z_sfc_halflevel,
                 p_prev_level,
@@ -469,7 +471,7 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_do_integral!(
                     local_geometry_prev_halflevel,
                 ),
                 TD.relative_humidity(thermo_params, tsʲ_prev_level),
-                ᶜphysical_buoyancy(params, ρ_prev_level, ρʲ_prev_level),
+                ᶜphysical_buoyancy(thermo_params, ρ_prev_level, ρʲ_prev_level),
                 get_physical_w(
                     u³_prev_halflevel,
                     local_geometry_prev_halflevel,
@@ -481,7 +483,7 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_do_integral!(
             )
 
             @. turb_entrʲ_prev_level = turbulent_entrainment(
-                params,
+                turbconv_params,
                 draft_area(ρaʲ_prev_level, ρʲ_prev_level),
             )
 
@@ -598,7 +600,8 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_do_integral!(
                 ) / local_geometry_level.J / ρ_level
 
             @. detrʲ_prev_level = detrainment_from_thermo_state(
-                params,
+                thermo_params,
+                turbconv_params,
                 z_prev_level,
                 z_sfc_halflevel,
                 p_prev_level,

diff --git a/src/cache/prognostic_edmf_precomputed_quantities.jl b/src/cache/prognostic_edmf_precomputed_quantities.jl
@@ -190,6 +190,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
     (; params) = p
     (; dt) = p
     thermo_params = CAP.thermodynamics_params(params)
+    turbconv_params = CAP.turbconv_params(params)
 
     FT = eltype(params)
     n = n_mass_flux_subdomains(turbconv_model)
@@ -217,15 +218,16 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
     for j in 1:n
         # entrainment/detrainment
         @. ᶜentrʲs.:($$j) = entrainment(
-            params,
+            thermo_params,
+            turbconv_params,
             ᶜz,
             z_sfc,
             ᶜp,
             Y.c.ρ,
             draft_area(Y.c.sgsʲs.:($$j).ρa, ᶜρʲs.:($$j)),
             get_physical_w(ᶜuʲs.:($$j), ᶜlg),
             TD.relative_humidity(thermo_params, ᶜtsʲs.:($$j)),
-            ᶜphysical_buoyancy(params, Y.c.ρ, ᶜρʲs.:($$j)),
+            ᶜphysical_buoyancy(thermo_params, Y.c.ρ, ᶜρʲs.:($$j)),
             get_physical_w(ᶜu, ᶜlg),
             TD.relative_humidity(thermo_params, ᶜts⁰),
             FT(0),
@@ -240,7 +242,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
         )
 
         @. ᶜturb_entrʲs.:($$j) = turbulent_entrainment(
-            params,
+            turbconv_params,
             draft_area(Y.c.sgsʲs.:($$j).ρa, ᶜρʲs.:($$j)),
         )
 
@@ -251,7 +253,8 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
         @. ᶜmassflux_vert_div =
             ᶜdivᵥ(ᶠinterp(Y.c.sgsʲs.:($$j).ρa) * ᶠu³ʲs.:($$j))
         @. ᶜdetrʲs.:($$j) = detrainment(
-            params,
+            thermo_params,
+            turbconv_params,
             ᶜz,
             z_sfc,
             ᶜp,
@@ -260,7 +263,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
             draft_area(Y.c.sgsʲs.:($$j).ρa, ᶜρʲs.:($$j)),
             get_physical_w(ᶜuʲs.:($$j), ᶜlg),
             TD.relative_humidity(thermo_params, ᶜtsʲs.:($$j)),
-            ᶜphysical_buoyancy(params, Y.c.ρ, ᶜρʲs.:($$j)),
+            ᶜphysical_buoyancy(thermo_params, Y.c.ρ, ᶜρʲs.:($$j)),
             get_physical_w(ᶜu, ᶜlg),
             TD.relative_humidity(thermo_params, ᶜts⁰),
             FT(0),

diff --git a/src/prognostic_equations/edmfx_closures.jl b/src/prognostic_equations/edmfx_closures.jl
@@ -3,6 +3,7 @@
 #####
 
 import StaticArrays as SA
+import Thermodynamics.Parameters as TDP
 import ClimaCore.Geometry as Geometry
 import ClimaCore.Fields as Fields
 
@@ -16,9 +17,9 @@ end
 """
     Return buoyancy on cell centers.
 """
-function ᶜphysical_buoyancy(params, ᶜρ_ref, ᶜρ)
+function ᶜphysical_buoyancy(thermo_params, ᶜρ_ref, ᶜρ)
     # TODO - replace by ᶜgradᵥᶠΦ when we move to deep atmosphere
-    g = CAP.grav(params)
+    g = TDP.grav(thermo_params)
     return (ᶜρ_ref - ᶜρ) / ᶜρ * g
 end
 """

diff --git a/src/prognostic_equations/edmfx_entr_detr.jl b/src/prognostic_equations/edmfx_entr_detr.jl
@@ -2,11 +2,14 @@
 ##### EDMF entrainment detrainment
 #####
 
+import Thermodynamics.Parameters as TDP
+
 """
    Return entrainment rate [1/s].
 
    Inputs (everything defined on cell centers):
-   - params set with model parameters
+   - `thermo_params` thermodynamics parameters
+   - `turbconv_params` turbulence convection parameters
    - ᶜz, z_sfc, ᶜp, ᶜρ, - grid-scale height, surface height, grid-scale pressure and density
    - buoy_flux_surface - buoyancy flux at the surface
    - ᶜaʲ, ᶜwʲ, ᶜRHʲ, ᶜbuoyʲ - updraft area, physical vertical velocity,
@@ -17,7 +20,8 @@
 """
 
 function entrainment(
-    params,
+    thermo_params,
+    turbconv_params,
     ᶜz,
     z_sfc,
     ᶜp,
@@ -32,11 +36,12 @@ function entrainment(
     ᶜtke⁰,
     ::NoEntrainment,
 )
-    return zero(eltype(params))
+    return zero(eltype(thermo_params))
 end
 
 function entrainment(
-    params,
+    thermo_params,
+    turbconv_params,
     ᶜz,
     z_sfc,
     ᶜp,
@@ -51,14 +56,14 @@ function entrainment(
     ᶜtke⁰,
     ::PiGroupsEntrainment,
 )
-    FT = eltype(params)
+    FT = eltype(thermo_params)
     if ᶜaʲ <= FT(0)
         return FT(0)
     else
-        g = CAP.grav(params)
+        g = TDP.grav(thermo_params)
         ref_H = ᶜp / (ᶜρ * g)
 
-        entr_param_vec = CAP.entr_param_vec(params)
+        entr_param_vec = CAP.entr_param_vec(turbconv_params)
 
         # non-dimensional pi-groups
         Π₁ = (ᶜz - z_sfc) * (ᶜbuoyʲ - ᶜbuoy⁰) / ((ᶜwʲ - ᶜw⁰)^2 + eps(FT)) / 100
@@ -85,7 +90,8 @@ function entrainment(
 end
 
 function entrainment(
-    params,
+    thermo_params,
+    turbconv_params,
     ᶜz,
     z_sfc,
     ᶜp,
@@ -100,8 +106,7 @@ function entrainment(
     ᶜtke⁰,
     ::GeneralizedEntrainment,
 )
-    FT = eltype(params)
-    turbconv_params = CAP.turbconv_params(params)
+    FT = eltype(thermo_params)
     entr_inv_tau = CAP.entr_tau(turbconv_params)
     entr_coeff = CAP.entr_coeff(turbconv_params)
     min_area_limiter_scale = CAP.min_area_limiter_scale(turbconv_params)
@@ -120,7 +125,8 @@ function entrainment(
 end
 
 function detrainment_from_thermo_state(
-    params,
+    thermo_params,
+    turbconv_params,
     z_prev_level,
     z_sfc_halflevel,
     p_prev_level,
@@ -139,10 +145,10 @@ function detrainment_from_thermo_state(
     tke_prev_level,
     edmfx_detr_model,
 )
-    FT = eltype(params)
-    thermo_params = CAP.thermodynamics_params(params)
+    FT = eltype(thermo_params)
     detrainment(
-        params,
+        thermo_params,
+        turbconv_params,
         z_prev_level,
         z_sfc_halflevel,
         p_prev_level,
@@ -151,7 +157,7 @@ function detrainment_from_thermo_state(
         draft_area(ρaʲ_prev_level, ρʲ_prev_level),
         get_physical_w(u³ʲ_prev_halflevel, local_geometry_prev_halflevel),
         TD.relative_humidity(thermo_params, tsʲ_prev_level),
-        ᶜphysical_buoyancy(params, ρ_prev_level, ρʲ_prev_level),
+        ᶜphysical_buoyancy(thermo_params, ρ_prev_level, ρʲ_prev_level),
         get_physical_w(u³_prev_halflevel, local_geometry_prev_halflevel),
         TD.relative_humidity(thermo_params, ts_prev_level),
         FT(0),
@@ -167,7 +173,8 @@ end
    Return detrainment rate [1/s].
 
    Inputs (everything defined on cell centers):
-   - params set with model parameters
+   - `thermo_params` thermodynamics parameters
+   - `turbconv_params` turbulence convection parameters
    - ᶜz, z_sfc, ᶜp, ᶜρ, - grid-scale height, surface height, grid-scale pressure and density
    - buoy_flux_surface - buoyancy flux at the surface
    - ᶜρaʲ, ᶜaʲ, ᶜwʲ, ᶜRHʲ, ᶜbuoyʲ - updraft effective density, updraft area, physical vertical velocity,
@@ -178,7 +185,8 @@ end
    - ᶜvert_div,ᶜmassflux_vert_div - vertical divergence, vertical mass flux divergence
 """
 function detrainment(
-    params,
+    thermo_params,
+    turbconv_params,
     ᶜz,
     z_sfc,
     ᶜp,
@@ -197,11 +205,12 @@ function detrainment(
     ᶜtke⁰,
     ::NoDetrainment,
 )
-    return zero(eltype(params))
+    return zero(eltype(thermo_params))
 end
 
 function detrainment(
-    params,
+    thermo_params,
+    turbconv_params,
     ᶜz,
     z_sfc,
     ᶜp,
@@ -220,14 +229,14 @@ function detrainment(
     ᶜtke⁰,
     ::PiGroupsDetrainment,
 )
-    FT = eltype(params)
+    FT = eltype(thermo_params)
     if ᶜaʲ <= FT(0)
         return FT(0)
     else
-        g = CAP.grav(params)
+        g = TDP.grav(thermo_params)
         ref_H = ᶜp / (ᶜρ * g)
 
-        entr_param_vec = CAP.entr_param_vec(params)
+        entr_param_vec = CAP.entr_param_vec(turbconv_params)
 
         # non-dimensional pi-groups
         Π₁ = (ᶜz - z_sfc) * (ᶜbuoyʲ - ᶜbuoy⁰) / ((ᶜwʲ - ᶜw⁰)^2 + eps(FT)) / 100
@@ -252,7 +261,8 @@ function detrainment(
 end
 
 function detrainment(
-    params,
+    thermo_params,
+    turbconv_params,
     ᶜz,
     z_sfc,
     ᶜp,
@@ -271,8 +281,7 @@ function detrainment(
     ᶜtke⁰,
     ::GeneralizedDetrainment,
 )
-    FT = eltype(params)
-    turbconv_params = CAP.turbconv_params(params)
+    FT = eltype(thermo_params)
     detr_inv_tau = CAP.detr_tau(turbconv_params)
     detr_coeff = CAP.detr_coeff(turbconv_params)
     detr_buoy_coeff = CAP.detr_buoy_coeff(turbconv_params)
@@ -303,7 +312,8 @@ function detrainment(
 end
 
 function detrainment(
-    params,
+    thermo_params,
+    turbconv_params,
     ᶜz,
     z_sfc,
     ᶜp,
@@ -325,8 +335,8 @@ function detrainment(
     return max(detr, 0)
 end
 
-function turbulent_entrainment(params, ᶜaʲ)
-    turb_entr_param_vec = CAP.turb_entr_param_vec(params)
+function turbulent_entrainment(turbconv_params, ᶜaʲ)
+    turb_entr_param_vec = CAP.turb_entr_param_vec(turbconv_params)
     return max(turb_entr_param_vec[1] * exp(-turb_entr_param_vec[2] * ᶜaʲ), 0)
 end