From c2daac01f79fed2d991cb66dcf1749dc9daa1d88 Mon Sep 17 00:00:00 2001 From: costachris Date: Tue, 25 Jun 2024 21:02:49 -0700 Subject: [PATCH] Change LES default imin index, expose surface area parameter in diagnostic EDMF, set top cell forcing tendency to zero --- src/cache/diagnostic_edmf_precomputed_quantities.jl | 7 ++++--- src/initial_conditions/initial_conditions.jl | 2 +- src/prognostic_equations/forcing/external_forcing.jl | 9 +++++++-- src/solver/model_getters.jl | 2 +- src/solver/types.jl | 5 ++--- src/surface_conditions/surface_setups.jl | 5 ++--- src/utils/read_gcm_driven_scm_data.jl | 4 ++-- 7 files changed, 19 insertions(+), 15 deletions(-) diff --git a/src/cache/diagnostic_edmf_precomputed_quantities.jl b/src/cache/diagnostic_edmf_precomputed_quantities.jl index 4496ebb380a..80222e7c59f 100644 --- a/src/cache/diagnostic_edmf_precomputed_quantities.jl +++ b/src/cache/diagnostic_edmf_precomputed_quantities.jl @@ -104,6 +104,7 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_bottom_bc!( (; ᶠu³⁰, ᶜK⁰) = p.precomputed thermo_params = CAP.thermodynamics_params(p.params) + turbconv_params = CAP.turbconv_params(params) ρ_int_level = Fields.field_values(Fields.level(Y.c.ρ, 1)) uₕ_int_level = Fields.field_values(Fields.level(Y.c.uₕ, 1)) @@ -155,7 +156,7 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_bottom_bc!( @. mseʲ_int_level = sgs_scalar_first_interior_bc( z_int_level - z_sfc_halflevel, ρ_int_level, - turbconv_model.a_int, + FT(turbconv_params.surface_area), h_tot_int_level - K_int_level, buoyancy_flux_sfc_halflevel, ρ_flux_h_tot_sfc_halflevel, @@ -166,7 +167,7 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_bottom_bc!( @. q_totʲ_int_level = sgs_scalar_first_interior_bc( z_int_level - z_sfc_halflevel, ρ_int_level, - turbconv_model.a_int, + FT(turbconv_params.surface_area), q_tot_int_level, buoyancy_flux_sfc_halflevel, ρ_flux_q_tot_sfc_halflevel, @@ -191,7 +192,7 @@ NVTX.@annotate function set_diagnostic_edmf_precomputed_quantities_bottom_bc!( p_int_level, Φ_int_level, ) - @. ρaʲ_int_level = ρʲ_int_level * turbconv_model.a_int + @. ρaʲ_int_level = ρʲ_int_level * FT(turbconv_params.surface_area) end ρaʲs_int_level = Fields.field_values(Fields.level(ᶜρaʲs, 1)) diff --git a/src/initial_conditions/initial_conditions.jl b/src/initial_conditions/initial_conditions.jl index 241e421eb60..b0034378e13 100644 --- a/src/initial_conditions/initial_conditions.jl +++ b/src/initial_conditions/initial_conditions.jl @@ -1177,7 +1177,7 @@ end function gcm_initial_conditions(external_forcing_file) NC.NCDataset(external_forcing_file) do ds ( # TODO: Cast to CuVector for GPU compatibility - gcm_driven_profile(ds, "thetali_mean")[:, 1], # 1 is initial time index + gcm_driven_profile(ds, "thetali_mean")[:, 1], # 1 is initial time index gcm_driven_profile(ds, "u_mean")[:, 1], gcm_driven_profile(ds, "v_mean")[:, 1], gcm_driven_profile(ds, "qt_mean")[:, 1], diff --git a/src/prognostic_equations/forcing/external_forcing.jl b/src/prognostic_equations/forcing/external_forcing.jl index 452d5295175..4353f122060 100644 --- a/src/prognostic_equations/forcing/external_forcing.jl +++ b/src/prognostic_equations/forcing/external_forcing.jl @@ -52,14 +52,13 @@ function external_forcing_cache(Y, external_forcing::GCMForcing) ᶜls_subsidence = similar(Y.c, FT) (; external_forcing_file) = external_forcing - imin = 100 # TODO: move into `GCMForcing` (and `parsed_args`) NC.Dataset(external_forcing_file, "r") do ds function setvar!(cc_field, varname, colidx, zc_gcm, zc_les) parent(cc_field[colidx]) .= interp_vertical_prof( zc_gcm, zc_les, - gcm_driven_profile_tmean(ds, varname; imin), # TODO: time-varying tendencies + gcm_driven_profile_tmean(ds, varname), # TODO: time-varying tendencies ) end @@ -181,6 +180,12 @@ function external_forcing_tendency!(Yₜ, Y, p, t, ::GCMForcing) Val{:first_order}(), ) # <-- subsidence + # needed to address top boundary condition for forcings. Otherwise upper portion of domain is anomalously warm + ρe_tot_top = Fields.level(Yₜ.c.ρe_tot, Spaces.nlevels(axes(Y.c))) + @. ρe_tot_top = 0.0 + + ρq_tot_top = Fields.level(Yₜ.c.ρq_tot, Spaces.nlevels(axes(Y.c))) + @. ρq_tot_top = 0.0 return nothing end diff --git a/src/solver/model_getters.jl b/src/solver/model_getters.jl index 7faf2042ffb..297fde28c3c 100644 --- a/src/solver/model_getters.jl +++ b/src/solver/model_getters.jl @@ -426,7 +426,7 @@ function get_turbconv_model(FT, parsed_args, turbconv_params) elseif turbconv == "diagnostic_edmfx" N = parsed_args["updraft_number"] TKE = parsed_args["prognostic_tke"] - DiagnosticEDMFX{N, TKE}(FT(0.1), turbconv_params.min_area) + DiagnosticEDMFX{N, TKE}(turbconv_params.surface_area) else nothing end diff --git a/src/solver/types.jl b/src/solver/types.jl index fe18df1fbd3..e467f96335b 100644 --- a/src/solver/types.jl +++ b/src/solver/types.jl @@ -176,11 +176,10 @@ PrognosticEDMFX{N, TKE}(a_half::FT) where {N, TKE, FT} = PrognosticEDMFX{N, TKE, FT}(a_half) struct DiagnosticEDMFX{N, TKE, FT} <: AbstractEDMF - a_int::FT # area fraction of the first interior cell above the surface a_half::FT # WARNING: this should never be used outside of divide_by_ρa end -DiagnosticEDMFX{N, TKE}(a_int::FT, a_half::FT) where {N, TKE, FT} = - DiagnosticEDMFX{N, TKE, FT}(a_int, a_half) +DiagnosticEDMFX{N, TKE}(a_half::FT) where {N, TKE, FT} = + DiagnosticEDMFX{N, TKE, FT}(a_half) n_mass_flux_subdomains(::PrognosticEDMFX{N}) where {N} = N n_mass_flux_subdomains(::DiagnosticEDMFX{N}) where {N} = N diff --git a/src/surface_conditions/surface_setups.jl b/src/surface_conditions/surface_setups.jl index a343a316e5d..8db78896dfe 100644 --- a/src/surface_conditions/surface_setups.jl +++ b/src/surface_conditions/surface_setups.jl @@ -252,14 +252,13 @@ end function (surface_setup::GCMDriven)(params) FT = eltype(params) (; external_forcing_file) = surface_setup - imin = 100 - T, lhf, shf = FT.(gcm_surface_conditions(external_forcing_file, imin)) + T, lhf, shf = FT.(gcm_surface_conditions(external_forcing_file)) z0 = FT(1e-4) # zrough parameterization = MoninObukhov(; z0, fluxes = HeatFluxes(; lhf, shf)) return SurfaceState(; parameterization, T) end -function gcm_surface_conditions(external_forcing_file, imin) +function gcm_surface_conditions(external_forcing_file; imin = 793) NC.NCDataset(external_forcing_file) do ds ( mean(gcm_driven_timeseries(ds, "surface_temperature")[imin:end]), diff --git a/src/utils/read_gcm_driven_scm_data.jl b/src/utils/read_gcm_driven_scm_data.jl index 07e9b579eb4..e65f5ee5887 100644 --- a/src/utils/read_gcm_driven_scm_data.jl +++ b/src/utils/read_gcm_driven_scm_data.jl @@ -3,7 +3,7 @@ import NCDatasets as NC import StatsBase: mean """ - gcm_driven_profile_tmean(ds, varname; imin = 100) + gcm_driven_profile_tmean(ds, varname; imin = 793) Extract time-averaged data (`imin:end`) for `varname` from the "profile" group in the GCM-driven dataset `ds` @@ -14,7 +14,7 @@ Returns a 1D ("z",) `Vector{FT}` of the time-averaged data. This method currently assumes the underlying data is `Float64`. If this is not the case, "garbage" data may be returned with no warning. """ -function gcm_driven_profile_tmean(ds, varname; imin = 100) +function gcm_driven_profile_tmean(ds, varname; imin = 793) vec(mean(gcm_driven_profile(ds, varname)[:, imin:end]; dims = 2)) end