Change LES default imin index, expose surface area parameter in diag…

…nostic EDMF, set top cell forcing tendency to zero

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, "theta_mean")[:, 1],
             gcm_driven_profile(ds, "u_mean")[:, 1],
             gcm_driven_profile(ds, "v_mean")[:, 1],
             gcm_driven_profile(ds, "qt_mean")[:, 1],

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
 
@@ -182,5 +181,11 @@ function external_forcing_tendency!(Yₜ, Y, p, t, ::GCMForcing)
     )
     # <-- subsidence
 
+    ρ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

src/solver/model_getters.jl

@@ -426,7 +426,10 @@ 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,
+            turbconv_params.min_area,
+        )
     else
         nothing
     end

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]),

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