Merge pull request #2940 from CliMA/zs/mixing_length

add different mixing lengths to cache and diagnostics

Author: szy21

config/model_configs/prognostic_edmfx_bomex_column.yml

@@ -33,5 +33,5 @@ diagnostics:
     period: 10mins
   - short_name: [arup, waup, taup, thetaaup, haup, husup, hurup, clwup, cliup, waen, taen, thetaaen, haen, husen, huren, clwen, clien, tke]
     period: 10mins
-  - short_name: [entr, detr, lmix, bgrad, strain, edt, evu]
+  - short_name: [entr, detr, lmix, lmixw, lmixtke, lmixb, bgrad, strain, edt, evu]
     period: 10mins

src/cache/diagnostic_edmf_precomputed_quantities.jl

@@ -766,7 +766,12 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_env_closures!
     (; ustar, obukhov_length) = p.precomputed.sfc_conditions
     (; ᶜρaʲs, ᶠu³ʲs, ᶜdetrʲs) = p.precomputed
     (; ᶜtke⁰, ᶠu³⁰, ᶜu⁰) = p.precomputed
-    (; ᶜlinear_buoygrad, ᶜstrain_rate_norm, ᶜmixing_length) = p.precomputed
+    (;
+        ᶜlinear_buoygrad,
+        ᶜstrain_rate_norm,
+        ᶜmixing_length_tuple,
+        ᶜmixing_length,
+    ) = p.precomputed
     (; ᶜK_h, ᶜK_u, ρatke_flux) = p.precomputed
     thermo_params = CAP.thermodynamics_params(params)
     ᶜlg = Fields.local_geometry_field(Y.c)
@@ -816,7 +821,7 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_env_closures!
 
     sfc_tke = Fields.level(ᶜtke⁰, 1)
     z_sfc = Fields.level(Fields.coordinate_field(Y.f).z, half)
-    @. ᶜmixing_length = mixing_length(
+    @. ᶜmixing_length_tuple = mixing_length(
         params,
         ustar,
         ᶜz,
@@ -830,6 +835,7 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_env_closures!
         ᶜprandtl_nvec,
         ᶜtke_exch,
     )
+    @. ᶜmixing_length = ᶜmixing_length_tuple.master
 
     turbconv_params = CAP.turbconv_params(params)
     c_m = CAP.tke_ed_coeff(turbconv_params)

src/cache/precomputed_quantities.jl

@@ -85,6 +85,7 @@ function precomputed_quantities(Y, atmos)
             ᶜq_tot⁰ = similar(Y.c, FT),
             ᶜts⁰ = similar(Y.c, TST),
             ᶜρ⁰ = similar(Y.c, FT),
+            ᶜmixing_length_tuple = similar(Y.c, MixingLength{FT}),
             ᶜK_u = similar(Y.c, FT),
             ᶜK_h = similar(Y.c, FT),
             ρatke_flux = similar(Fields.level(Y.f, half), C3{FT}),
@@ -126,6 +127,7 @@ function precomputed_quantities(Y, atmos)
             ᶜu⁰ = similar(Y.c, C123{FT}),
             ᶜK⁰ = similar(Y.c, FT),
             ᶜtke⁰ = similar(Y.c, FT),
+            ᶜmixing_length_tuple = similar(Y.c, MixingLength{FT}),
             ᶜK_u = similar(Y.c, FT),
             ᶜK_h = similar(Y.c, FT),
             ρatke_flux = similar(Fields.level(Y.f, half), C3{FT}),

src/cache/prognostic_edmf_precomputed_quantities.jl

@@ -196,6 +196,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
 
     (; ᶜtke⁰, ᶜu, ᶜp, ᶜρa⁰, ᶠu³⁰, ᶜts⁰, ᶜq_tot⁰) = p.precomputed
     (;
+        ᶜmixing_length_tuple,
         ᶜmixing_length,
         ᶜlinear_buoygrad,
         ᶜstrain_rate_norm,
@@ -307,7 +308,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
     end
 
     sfc_tke = Fields.level(ᶜtke⁰, 1)
-    @. ᶜmixing_length = mixing_length(
+    @. ᶜmixing_length_tuple = mixing_length(
         p.params,
         ustar,
         ᶜz,
@@ -322,6 +323,8 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_closures!(
         ᶜtke_exch,
     )
 
+    @. ᶜmixing_length = ᶜmixing_length_tuple.master
+
     turbconv_params = CAP.turbconv_params(params)
     c_m = CAP.tke_ed_coeff(turbconv_params)
     @. ᶜK_u = c_m * ᶜmixing_length * sqrt(max(ᶜtke⁰, 0))

src/diagnostics/edmfx_diagnostics.jl

@@ -882,6 +882,93 @@ add_diagnostic_variable!(
     compute! = compute_tke!,
 )
 
+###
+# Environment Wall Constrained Mixing Length (3d)
+###
+compute_lmixw!(out, state, cache, time) =
+    compute_lmixw!(out, state, cache, time, cache.atmos.turbconv_model)
+compute_lmixw!(_, _, _, _, turbconv_model::T) where {T} =
+    error_diagnostic_variable("lmixw", turbconv_model)
+
+function compute_lmixw!(
+    out,
+    state,
+    cache,
+    time,
+    turbconv_model::Union{PrognosticEDMFX, DiagnosticEDMFX},
+)
+    if isnothing(out)
+        return copy(cache.precomputed.ᶜmixing_length_tuple.wall)
+    else
+        out .= cache.precomputed.ᶜmixing_length_tuple.wall
+    end
+end
+
+add_diagnostic_variable!(
+    short_name = "lmixw",
+    long_name = "Environment Wall Constrained Mixing Length",
+    units = "m",
+    compute! = compute_lmixw!,
+)
+
+###
+# Environment TKE Balanced Mixing Length (3d)
+###
+compute_lmixtke!(out, state, cache, time) =
+    compute_lmixtke!(out, state, cache, time, cache.atmos.turbconv_model)
+compute_lmixtke!(_, _, _, _, turbconv_model::T) where {T} =
+    error_diagnostic_variable("lmixtke", turbconv_model)
+
+function compute_lmixtke!(
+    out,
+    state,
+    cache,
+    time,
+    turbconv_model::Union{PrognosticEDMFX, DiagnosticEDMFX},
+)
+    if isnothing(out)
+        return copy(cache.precomputed.ᶜmixing_length_tuple.tke)
+    else
+        out .= cache.precomputed.ᶜmixing_length_tuple.tke
+    end
+end
+
+add_diagnostic_variable!(
+    short_name = "lmixtke",
+    long_name = "Environment TKE Balanced Mixing Length",
+    units = "m",
+    compute! = compute_lmixtke!,
+)
+
+###
+# Environment Stability Mixing Length (3d)
+###
+compute_lmixb!(out, state, cache, time) =
+    compute_lmixb!(out, state, cache, time, cache.atmos.turbconv_model)
+compute_lmixb!(_, _, _, _, turbconv_model::T) where {T} =
+    error_diagnostic_variable("lmixb", turbconv_model)
+
+function compute_lmixb!(
+    out,
+    state,
+    cache,
+    time,
+    turbconv_model::Union{PrognosticEDMFX, DiagnosticEDMFX},
+)
+    if isnothing(out)
+        return copy(cache.precomputed.ᶜmixing_length_tuple.buoy)
+    else
+        out .= cache.precomputed.ᶜmixing_length_tuple.buoy
+    end
+end
+
+add_diagnostic_variable!(
+    short_name = "lmixb",
+    long_name = "Environment Static Stability Mixing Length",
+    units = "m",
+    compute! = compute_lmixb!,
+)
+
 ###
 # Diffusivity of heat (3d)
 ###

src/prognostic_equations/edmfx_closures.jl

@@ -251,7 +251,7 @@ function mixing_length(
     l_smin = lamb_smooth_minimum(l, smin_ub, smin_rm)
     l_limited = max(l_smag, min(l_smin, l_z))
 
-    return l_limited
+    return MixingLength{FT}(l_limited, l_W, l_TKE, l_N)
 end
 
 """

src/solver/types.jl

@@ -153,6 +153,12 @@ end
 Base.eltype(::EnvBuoyGradVars{FT}) where {FT} = FT
 Base.broadcastable(x::EnvBuoyGradVars) = tuple(x)
 
+struct MixingLength{FT}
+    master::FT
+    wall::FT
+    tke::FT
+    buoy::FT
+end
 
 abstract type AbstractEDMF end