Skip to content

Commit 0107433

Browse files
costachriscostachris
committed
debug
1 parent 82462e9 commit 0107433

15 files changed

+90
-61
lines changed

‎src/cache/precomputed_quantities.jl

Lines changed: 9 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -3,6 +3,7 @@
33
#####
44
import Thermodynamics as TD
55
import ClimaCore: Spaces, Fields
6+
using Base.Broadcast: materialize
67

78
"""
89
implicit_precomputed_quantities(Y, atmos)
@@ -372,21 +373,20 @@ function thermo_state(
372373
end
373374

374375
function thermo_vars(moisture_model, precip_model, Y_c, K, Φ)
375-
# Compute specific quantities on-the-fly
376-
e_tot = @. lazy(specific(Y_c.ρe_tot, Y_c.ρ))
376+
e_tot = materialize(@. lazy(specific(Y_c.ρe_tot, Y_c.ρ)))
377377
energy_var = (; e_int = e_tot - K - Φ)
378-
378+
379379
moisture_var = if moisture_model isa DryModel
380380
(;)
381381
elseif moisture_model isa EquilMoistModel
382-
q_tot = @. lazy(specific(Y_c.ρq_tot, Y_c.ρ))
382+
q_tot = materialize(@. lazy(specific(Y_c.ρq_tot, Y_c.ρ)))
383383
(; q_tot)
384384
elseif moisture_model isa NonEquilMoistModel
385-
q_tot = @. lazy(specific(Y_c.ρq_tot, Y_c.ρ))
386-
q_liq = @. lazy(specific(Y_c.ρq_liq, Y_c.ρ))
387-
q_ice = @. lazy(specific(Y_c.ρq_ice, Y_c.ρ))
388-
q_rai = @. lazy(specific(Y_c.ρq_rai, Y_c.ρ))
389-
q_sno = @. lazy(specific(Y_c.ρq_sno, Y_c.ρ))
385+
q_tot = materialize(@. lazy(specific(Y_c.ρq_tot, Y_c.ρ)))
386+
q_liq = materialize(@. lazy(specific(Y_c.ρq_liq, Y_c.ρ)))
387+
q_ice = materialize(@. lazy(specific(Y_c.ρq_ice, Y_c.ρ)))
388+
q_rai = materialize(@. lazy(specific(Y_c.ρq_rai, Y_c.ρ)))
389+
q_sno = materialize(@. lazy(specific(Y_c.ρq_sno, Y_c.ρ)))
390390
q_pt_args = (q_tot, q_liq + q_rai, q_ice + q_sno)
391391
(; q_pt = TD.PhasePartition(q_pt_args...))
392392
end
@@ -439,7 +439,6 @@ NVTX.@annotate function set_implicit_precomputed_quantities!(Y, p, t)
439439
thermo_params = CAP.thermodynamics_params(p.params)
440440
thermo_args = (thermo_params, moisture_model, precip_model)
441441

442-
ᶜspecific .= ᶜspecific_gs_tracers(Y)
443442
@. ᶠuₕ³ = $compute_ᶠuₕ³(Y.c.uₕ, Y.c.ρ)
444443

445444
# TODO: We might want to move this to dss! (and rename dss! to something

‎src/cache/prognostic_edmf_precomputed_quantities.jl

Lines changed: 18 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -25,11 +25,15 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_environment!(
2525
(; ᶠu₃⁰, ᶜu⁰, ᶠu³⁰, ᶜK⁰, ᶜts⁰) = p.precomputed
2626

2727
ᶜρa⁰ = @.lazy(ρa⁰(Y.c))
28-
@. ᶜtke⁰ = specific_tke(Y.c.sgs⁰, Y.c, turbconv_model)
28+
ᶜtke⁰ = @. lazy(specific_tke(Y.c.sgs⁰, Y.c, turbconv_model))
2929
set_sgs_ᶠu₃!(u₃⁰, ᶠu₃⁰, Y, turbconv_model)
3030
set_velocity_quantities!(ᶜu⁰, ᶠu³⁰, ᶜK⁰, ᶠu₃⁰, Y.c.uₕ, ᶠuₕ³)
3131
# @. ᶜK⁰ += ᶜtke⁰
3232
ᶜq_tot⁰ = @.lazy(specific_env_value(:q_tot, Y.c, turbconv_model))
33+
34+
ᶜmse⁰ = p.scratch.ᶜtemp_scalar_2
35+
ᶜmse⁰ .= specific_env_mse(Y.c, p)
36+
3337
if p.atmos.moisture_model isa NonEquilMoistModel &&
3438
p.atmos.precip_model isa Microphysics1Moment
3539
ᶜq_liq⁰ = @.lazy(specific_env_value(:q_liq, Y.c, turbconv_model))
@@ -39,14 +43,15 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_environment!(
3943
@. ᶜts⁰ = TD.PhaseNonEquil_phq(
4044
thermo_params,
4145
ᶜp,
42-
specific_env_mse(Y.c, p) - ᶜΦ,
46+
ᶜmse⁰ - ᶜΦ,
4347
TD.PhasePartition(ᶜq_tot⁰, ᶜq_liq⁰ + ᶜq_rai⁰, ᶜq_ice⁰ + ᶜq_sno⁰),
4448
)
4549
else
50+
4651
@. ᶜts⁰ = TD.PhaseEquil_phq(
4752
thermo_params,
4853
ᶜp,
49-
specific_env_mse(Y.c, p) - ᶜΦ,
54+
ᶜmse⁰ - ᶜΦ,
5055
ᶜq_tot⁰,
5156
)
5257
end
@@ -132,7 +137,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_bottom_bc!(
132137
turbconv_params = CAP.turbconv_params(p.params)
133138

134139
(; ᶜΦ,) = p.core
135-
(; ᶜp, ᶜK, ᶜtsʲs, ᶜρʲs) = p.precomputed
140+
(; ᶜp, ᶜK, ᶜtsʲs, ᶜρʲs, ᶜts) = p.precomputed
136141
(; ustar, obukhov_length, buoyancy_flux) = p.precomputed.sfc_conditions
137142

138143
for j in 1:n
@@ -176,8 +181,8 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_bottom_bc!(
176181
# TODO: replace this with the actual surface area fraction when
177182
# using prognostic surface area
178183
@. ᶜaʲ_int_val = FT(turbconv_params.surface_area)
179-
ᶜh_tot = @. lazy(TD.total_specific_enthalpy(thermo_params, ᶜtsʲ, specific(Y.c.ρe_tot, Y.c.ρ)))
180-
ᶜh_tot_int_val = Fields.field_values(Fields.level(ᶜh_tot, 1))
184+
ᶜh_tot = @. lazy(TD.total_specific_enthalpy(thermo_params, ᶜts, specific(Y.c.ρe_tot, Y.c.ρ)))
185+
ᶜh_tot_int_val = Fields.field_values(Fields.level(Base.materialize(ᶜh_tot), 1))
181186
ᶜK_int_val = Fields.field_values(Fields.level(ᶜK, 1))
182187
ᶜmseʲ_int_val = Fields.field_values(Fields.level(ᶜmseʲ, 1))
183188
@. ᶜmseʲ_int_val = sgs_scalar_first_interior_bc(
@@ -193,8 +198,10 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_bottom_bc!(
193198
)
194199

195200
# ... and the first interior point for EDMFX ᶜq_totʲ.
201+
202+
ᶜq_tot = @. lazy(specific(Y.c.ρq_tot, Y.c.ρ))
196203
ᶜq_tot_int_val =
197-
Fields.field_values(Fields.level(specific(Y.c.ρq_tot, Y.c.ρ), 1))
204+
Fields.field_values(Fields.level(Base.materialize(ᶜq_tot), 1))
198205
ᶜq_totʲ_int_val = Fields.field_values(Fields.level(ᶜq_totʲ, 1))
199206
@. ᶜq_totʲ_int_val = sgs_scalar_first_interior_bc(
200207
ᶜz_int_val - z_sfc_val,
@@ -359,6 +366,8 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_explicit_clos
359366
ᶜdz = Fields.Δz_field(axes(Y.c))
360367
ᶜlg = Fields.local_geometry_field(Y.c)
361368
ᶠlg = Fields.local_geometry_field(Y.f)
369+
ᶜtke⁰ = @. lazy(specific_tke(Y.c.sgs⁰, Y.c, turbconv_model))
370+
ᶜρa⁰ = @. lazy(ρa⁰(Y.c))
362371

363372
ᶜvert_div = p.scratch.ᶜtemp_scalar
364373
ᶜmassflux_vert_div = p.scratch.ᶜtemp_scalar_2
@@ -480,7 +489,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_explicit_clos
480489
(1 / 2 * norm_sqr(ᶜinterp(ᶠu³⁰) - ᶜinterp(ᶠu³ʲs.:($$j))) - ᶜtke⁰)
481490
end
482491

483-
sfc_tke = Fields.level(ᶜtke⁰, 1)
492+
sfc_tke = Fields.level(Base.materialize(ᶜtke⁰), 1)
484493
@. ᶜmixing_length_tuple = mixing_length(
485494
p.params,
486495
ustar,
@@ -503,7 +512,7 @@ NVTX.@annotate function set_prognostic_edmf_precomputed_quantities_explicit_clos
503512
@. ᶜK_h = eddy_diffusivity(ᶜK_u, ᶜprandtl_nvec)
504513

505514
ρatke_flux_values = Fields.field_values(ρatke_flux)
506-
ρa_sfc_values = Fields.field_values(Fields.level(ᶜρa⁰, 1)) # TODO: replace by surface value
515+
ρa_sfc_values = Fields.field_values(Fields.level(Base.materialize(ᶜρa⁰), 1)) # TODO: replace by surface value
507516
ustar_values = Fields.field_values(ustar)
508517
sfc_local_geometry_values = Fields.field_values(
509518
Fields.level(Fields.local_geometry_field(Y.f), half),

‎src/diagnostics/Diagnostics.jl

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -58,6 +58,8 @@ import ..PrognosticSurfaceTemperature
5858
import ..draft_area
5959
import ..compute_gm_mixing_length!
6060
import ..horizontal_integral_at_boundary
61+
import ..ρa⁰
62+
import ..specific_tke
6163

6264
# We need the abbreviations for symbols like curl, grad, and so on
6365
include(joinpath("..", "utils", "abbreviations.jl"))

‎src/diagnostics/edmfx_diagnostics.jl

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -633,7 +633,7 @@ compute_aren!(_, _, _, _, turbconv_model::T) where {T} =
633633

634634
function compute_aren!(out, state, cache, time, turbconv_model::PrognosticEDMFX)
635635
thermo_params = CAP.thermodynamics_params(cache.params)
636-
ᶜρa⁰ = @.lazy(ρa⁰(state.c))
636+
ᶜρa⁰ = @. lazy(ρa⁰(state.c))
637637
if isnothing(out)
638638
return draft_area.(
639639
ᶜρa⁰,
@@ -1146,11 +1146,11 @@ function compute_tke!(
11461146
time,
11471147
turbconv_model::Union{EDOnlyEDMFX, PrognosticEDMFX, DiagnosticEDMFX},
11481148
)
1149-
1149+
ᶜtke = @. lazy(specific_tke(state.c.sgs⁰, state.c, turbconv_model))
11501150
if isnothing(out)
1151-
return specific_tke(state.c.sgs⁰, state.c, turbconv_model)
1151+
return Base.materialize(ᶜtke)
11521152
else
1153-
out .= specific_tke(state.c.sgs⁰, state.c, turbconv_model)
1153+
out .= ᶜtke
11541154
end
11551155
end
11561156

‎src/parameterized_tendencies/les_sgs_models/smagorinsky_lilly.jl

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -96,7 +96,7 @@ horizontal_smagorinsky_lilly_tendency!(Yₜ, Y, p, t, ::Nothing) = nothing
9696
vertical_smagorinsky_lilly_tendency!(Yₜ, Y, p, t, ::Nothing) = nothing
9797

9898
function horizontal_smagorinsky_lilly_tendency!(Yₜ, Y, p, t, ::SmagorinskyLilly)
99-
(; ᶜτ_smag, ᶠτ_smag, ᶜD_smag) = p.precomputed
99+
(; ᶜτ_smag, ᶠτ_smag, ᶜD_smag, ᶜts) = p.precomputed
100100

101101
## Momentum tendencies
102102
ᶠρ = @. p.scratch.ᶠtemp_scalar = ᶠinterp(Y.c.ρ)
@@ -128,6 +128,7 @@ function vertical_smagorinsky_lilly_tendency!(Yₜ, Y, p, t, ::SmagorinskyLilly)
128128
(; ᶜτ_smag, ᶠτ_smag, ᶠD_smag, sfc_conditions) =
129129
p.precomputed
130130
(; ρ_flux_uₕ, ρ_flux_h_tot) = sfc_conditions
131+
(; ᶜts) = p.precomputed
131132

132133
# Define operators
133134
ᶠgradᵥ = Operators.GradientC2F() # apply BCs to ᶜdivᵥ, which wraps ᶠgradᵥ

‎src/prognostic_equations/advection.jl

Lines changed: 4 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -36,7 +36,9 @@ Modifies `Yₜ.c.ρ`, `Yₜ.c.ρe_tot`, `Yₜ.c.uₕ`, and EDMFX-related fields
3636
NVTX.@annotate function horizontal_dynamics_tendency!(Yₜ, Y, p, t)
3737
n = n_mass_flux_subdomains(p.atmos.turbconv_model)
3838
(; ᶜΦ) = p.core
39-
(; ᶜu, ᶜK, ᶜp) = p.precomputed
39+
(; ᶜu, ᶜK, ᶜp, ᶜts) = p.precomputed
40+
(; params) = p
41+
thermo_params = CAP.thermodynamics_params(params)
4042

4143
if p.atmos.turbconv_model isa PrognosticEDMFX
4244
(; ᶜuʲs) = p.precomputed
@@ -174,7 +176,7 @@ NVTX.@annotate function explicit_vertical_advection_tendency!(Yₜ, Y, p, t)
174176
(; dt) = p
175177
ᶜJ = Fields.local_geometry_field(Y.c).J
176178
(; ᶜf³, ᶠf¹², ᶜΦ) = p.core
177-
(; ᶜu, ᶠu³, ᶜK) = p.precomputed
179+
(; ᶜu, ᶠu³, ᶜK, ᶜts) = p.precomputed
178180
(; edmfx_upwinding) = n > 0 || advect_tke ? p.atmos.numerics : all_nothing
179181
(; ᶜuʲs, ᶜKʲs, ᶠKᵥʲs) = n > 0 ? p.precomputed : all_nothing
180182
(; energy_upwinding, tracer_upwinding) = p.atmos.numerics

‎src/prognostic_equations/edmfx_entr_detr.jl

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -530,7 +530,8 @@ function edmfx_entr_detr_tendency!(Yₜ, Y, p, t, turbconv_model::PrognosticEDMF
530530
(; ᶜturb_entrʲs, ᶜentrʲs, ᶜdetrʲs) = p.precomputed
531531
(; ᶠu₃⁰) = p.precomputed
532532

533-
ᶜmse⁰ = @.lazy(specific_env_mse(Y.c, p))
533+
ᶜmse⁰ = p.scratch.ᶜtemp_scalar
534+
ᶜmse⁰ .= specific_env_mse(Y.c, p)
534535
if p.atmos.moisture_model isa NonEquilMoistModel &&
535536
p.atmos.precip_model isa Microphysics1Moment
536537
ᶜq_liq⁰ = @.lazy(specific_env_value(:q_liq, Y.c, turbconv_model))

‎src/prognostic_equations/edmfx_sgs_flux.jl

Lines changed: 14 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -42,7 +42,7 @@ function edmfx_sgs_mass_flux_tendency!(
4242
(; edmfx_sgsflux_upwinding) = p.atmos.numerics
4343
(; ᶠu³) = p.precomputed
4444
(; ᶠu³ʲs, ᶜKʲs, ᶜρʲs) = p.precomputed
45-
(; ᶠu³⁰, ᶜK⁰, ᶜts⁰) = p.precomputed
45+
(; ᶠu³⁰, ᶜK⁰, ᶜts⁰, ᶜts) = p.precomputed
4646
thermo_params = CAP.thermodynamics_params(p.params)
4747
ᶜρ⁰ = @. TD.air_density(thermo_params, ᶜts⁰)
4848
ᶜρa⁰ = @.lazy(ρa⁰(Y.c))
@@ -56,7 +56,7 @@ function edmfx_sgs_mass_flux_tendency!(
5656
# [best after removal of precomputed quantities]
5757
ᶠu³_diff = p.scratch.ᶠtemp_CT3
5858
ᶜa_scalar = p.scratch.ᶜtemp_scalar
59-
ᶜh_tot = @. lazy(TD.total_specific_enthalpy(thermo_params, ᶜtsʲ, specific(Y.c.ρe_tot, Y.c.ρ)))
59+
ᶜh_tot = @. lazy(TD.total_specific_enthalpy(thermo_params, ᶜts, specific(Y.c.ρe_tot, Y.c.ρ)))
6060
for j in 1:n
6161
@. ᶠu³_diff = ᶠu³ʲs.:($$j) - ᶠu³
6262
@. ᶜa_scalar =
@@ -73,7 +73,9 @@ function edmfx_sgs_mass_flux_tendency!(
7373
end
7474
# Add the environment fluxes
7575
@. ᶠu³_diff = ᶠu³⁰ - ᶠu³
76-
ᶜmse⁰ = @.lazy(specific_env_mse(Y.c, p))
76+
77+
ᶜmse⁰ = p.scratch.ᶜtemp_scalar_2
78+
ᶜmse⁰ .= specific_env_mse(Y.c, p)
7779
@. ᶜa_scalar = (ᶜmse⁰ + ᶜK⁰ - ᶜh_tot) * draft_area(ᶜρa⁰, ᶜρ⁰)
7880
vtt = vertical_transport(
7981
ᶜρ⁰,
@@ -92,7 +94,7 @@ function edmfx_sgs_mass_flux_tendency!(
9294
(Y.c.sgsʲs.:($$j).q_tot - specific(Y.c.ρq_tot, Y.c.ρ)) *
9395
draft_area(Y.c.sgsʲs.:($$j).ρa, ᶜρʲs.:($$j))
9496
vtt = vertical_transport(
95-
ᶜρʲs.:($$j),
97+
ᶜρʲs.:($j),
9698
ᶠu³_diff,
9799
ᶜa_scalar,
98100
dt,
@@ -254,14 +256,14 @@ function edmfx_sgs_mass_flux_tendency!(
254256
n = n_mass_flux_subdomains(turbconv_model)
255257
(; edmfx_sgsflux_upwinding) = p.atmos.numerics
256258
(; ᶠu³) = p.precomputed
257-
(; ᶜρaʲs, ᶜρʲs, ᶠu³ʲs, ᶜKʲs, ᶜmseʲs, ᶜq_totʲs) = p.precomputed
259+
(; ᶜρaʲs, ᶜρʲs, ᶠu³ʲs, ᶜKʲs, ᶜmseʲs, ᶜq_totʲs, ᶜts) = p.precomputed
258260
(; dt) = p
259261
ᶜJ = Fields.local_geometry_field(Y.c).J
260262
FT = eltype(Y)
261263

262264
if p.atmos.edmfx_model.sgs_mass_flux isa Val{true}
263265
# energy
264-
ᶜh_tot = @. lazy(TD.total_specific_enthalpy(thermo_params, ᶜtsʲ, specific(Y.c.ρe_tot, Y.c.ρ)))
266+
ᶜh_tot = @. lazy(TD.total_specific_enthalpy(thermo_params, ᶜts, specific(Y.c.ρe_tot, Y.c.ρ)))
265267
ᶠu³_diff = p.scratch.ᶠtemp_CT3
266268
ᶜa_scalar = p.scratch.ᶜtemp_scalar
267269
for j in 1:n
@@ -396,7 +398,7 @@ function edmfx_sgs_diffusive_flux_tendency!(
396398
turbconv_params = CAP.turbconv_params(params)
397399
c_d = CAP.tke_diss_coeff(turbconv_params)
398400
(; ᶜu⁰, ᶜK⁰, ᶜlinear_buoygrad, ᶜstrain_rate_norm,) = p.precomputed
399-
(; ρatke_flux) = p.precomputed
401+
(; ᶜmixing_length, ᶜK_u, ᶜK_h, ρatke_flux) = p.precomputed
400402
ᶠgradᵥ = Operators.GradientC2F()
401403
ᶜρa⁰ = @.lazy(ρa⁰(Y.c))
402404
ᶜtke⁰ = @.lazy(specific_tke(Y.c.sgs⁰, Y.c, turbconv_model))
@@ -412,7 +414,9 @@ function edmfx_sgs_diffusive_flux_tendency!(
412414
top = Operators.SetValue(C3(FT(0))),
413415
bottom = Operators.SetValue(C3(FT(0))),
414416
)
415-
ᶜmse⁰ = @.lazy(specific_env_mse(Y.c, p))
417+
# ᶜmse⁰ = @. lazy(specific_env_mse(Y.c, p))
418+
ᶜmse⁰ = p.scratch.ᶜtemp_scalar_2
419+
ᶜmse⁰ .= specific_env_mse(Y.c, p)
416420
@. Yₜ.c.ρe_tot -= ᶜdivᵥ_ρe_tot(-(ᶠρaK_h * ᶠgradᵥ(ᶜmse⁰ + ᶜK⁰)))
417421
if use_prognostic_tke(turbconv_model)
418422
# Turbulent TKE transport (diffusion)
@@ -503,7 +507,7 @@ function edmfx_sgs_diffusive_flux_tendency!(
503507
turbconv_params = CAP.turbconv_params(params)
504508
thermo_params = CAP.thermodynamics_params(params)
505509
c_d = CAP.tke_diss_coeff(turbconv_params)
506-
(; ᶜu, ᶜmixing_length) = p.precomputed
510+
(; ᶜu, ᶜmixing_length, ᶜts) = p.precomputed
507511
(; ᶜK_u, ᶜK_h, ρatke_flux) = p.precomputed
508512
ᶠgradᵥ = Operators.GradientC2F()
509513
ᶜtke⁰ = @.lazy(specific_tke(Y.c.sgs⁰, Y.c, turbconv_model))
@@ -519,7 +523,7 @@ function edmfx_sgs_diffusive_flux_tendency!(
519523
top = Operators.SetValue(C3(FT(0))),
520524
bottom = Operators.SetValue(C3(FT(0))),
521525
)
522-
ᶜh_tot = @. lazy(TD.total_specific_enthalpy(thermo_params, ᶜtsʲ, specific(Y.c.ρe_tot, Y.c.ρ)))
526+
ᶜh_tot = @. lazy(TD.total_specific_enthalpy(thermo_params, ᶜts, specific(Y.c.ρe_tot, Y.c.ρ)))
523527
@. Yₜ.c.ρe_tot -= ᶜdivᵥ_ρe_tot(-(ᶠρaK_h * ᶠgradᵥ(ᶜh_tot)))
524528

525529
if use_prognostic_tke(turbconv_model)

‎src/prognostic_equations/edmfx_tke.jl

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -48,7 +48,9 @@ function edmfx_tke_tendency!(
4848
ᶠu³⁰,
4949
ᶠu³,
5050
ᶜstrain_rate_norm,
51-
ᶜlinear_buoygrad
51+
ᶜlinear_buoygrad,
52+
ᶜK_u,
53+
ᶜK_h,
5254
) = p.precomputed
5355
turbconv_params = CAP.turbconv_params(p.params)
5456
FT = eltype(p.params)

‎src/prognostic_equations/forcing/subsidence.jl

Lines changed: 9 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -84,30 +84,30 @@ Arguments:
8484
- `Y`: The current state vector, used for density (`ρ`).
8585
- `p`: Cache containing parameters, and the subsidence model object.
8686
- `t`: Current simulation time.
87-
- `subsidence`: The subsidence model object, containing the prescribed vertical
88-
velocity profile `Dᵥ`.
87+
- `subsidence`: The subsidence model object.
8988
"""
90-
function subsidence_tendency!(Yₜ, Y, p, t, subsidence::Subsidence)
91-
(; Dᵥ) = subsidence
92-
ᶜρ = Y.c.ρ
89+
function subsidence_tendency!(Yₜ, Y, p, t, ::Subsidence)
9390
(; moisture_model) = p.atmos
9491
subsidence_profile = p.atmos.subsidence.prof
9592
thermo_params = CAP.thermodynamics_params(p.params)
93+
(; ᶜts) = p.precomputed
9694

9795
ᶠz = Fields.coordinate_field(axes(Y.f)).z
9896
ᶠlg = Fields.local_geometry_field(Y.f)
9997
ᶠsubsidence³ = p.scratch.ᶠtemp_CT3
10098
@. ᶠsubsidence³ =
10199
subsidence_profile(ᶠz) * CT3(unit_basis_vector_data(CT3, ᶠlg))
102-
100+
101+
ᶜρ = Y.c.ρ
103102
# LS Subsidence
104103
ᶜh_tot = @. lazy(TD.total_specific_enthalpy(thermo_params, ᶜts, specific(Y.c.ρe_tot, Y.c.ρ)))
104+
ᶜq_tot = @. lazy(specific(Y.c.ρq_tot, Y.c.ρ))
105105
subsidence!(Yₜ.c.ρe_tot, Y.c.ρ, ᶠsubsidence³, ᶜh_tot, Val{:first_order}())
106106
subsidence!(
107107
Yₜ.c.ρq_tot,
108108
Y.c.ρ,
109109
ᶠsubsidence³,
110-
specific(Y.c.ρq_tot, Y.c.ρ),
110+
ᶜq_tot,
111111
Val{:first_order}(),
112112
)
113113
if moisture_model isa NonEquilMoistModel
@@ -116,15 +116,15 @@ function subsidence_tendency!(Yₜ, Y, p, t, subsidence::Subsidence)
116116
Yₜ.c.ρq_liq,
117117
Y.c.ρ,
118118
ᶠsubsidence³,
119-
specific(Y.c.ρq_liq, Y.c.ρ),
119+
ᶜq_liq,
120120
Val{:first_order}(),
121121
)
122122
ᶜq_ice = @. lazy(specific(Y.c.ρq_ice, Y.c.ρ))
123123
subsidence!(
124124
Yₜ.c.ρq_ice,
125125
Y.c.ρ,
126126
ᶠsubsidence³,
127-
specific(Y.c.ρq_ice, Y.c.ρ),
127+
ᶜq_ice,
128128
Val{:first_order}(),
129129
)
130130
end

0 commit comments

Comments
 (0)