Skip to content

Commit

Permalink
review comments
Browse files Browse the repository at this point in the history
  • Loading branch information
kmdeck authored and juliasloan25 committed Sep 6, 2024
1 parent ee2042f commit 7eb5cb7
Show file tree
Hide file tree
Showing 3 changed files with 21 additions and 19 deletions.
26 changes: 13 additions & 13 deletions experiments/standalone/Vegetation/timestep_test.jl
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,7 @@ using ClimaLand.Canopy
using ClimaLand.Canopy.PlantHydraulics
import ClimaLand
import ClimaLand.Parameters as LP
const FT = Float32;
const FT = Float64;
earth_param_set = LP.LandParameters(FT);
f_root_to_shoot = FT(3.5)
plant_ν = FT(2.46e-4) # kg/m^2
Expand Down Expand Up @@ -94,13 +94,12 @@ AR_model = AutotrophicRespirationModel{FT}(AR_params);

f_root_to_shoot = FT(3.5)
SAI = FT(1.0)
RAI = FT(3 * f_root_to_shoot)
RAI = FT(3f_root_to_shoot)
ai_parameterization = PrescribedSiteAreaIndex{FT}(LAIfunction, SAI, RAI)

function root_distribution(z::T) where {T}
rooting_depth = T(1.0)
return T(1.0 / rooting_depth) * exp(z / T(rooting_depth))
end;
rooting_depth = FT(1.0)
function root_distribution(z::T; rooting_depth = rooting_depth) where {T}
return T(1.0 / rooting_depth) * exp(z / T(rooting_depth)) # 1/m
end

K_sat_plant = FT(1.8e-6)
ψ63 = FT(-4 / 0.0098)
Expand Down Expand Up @@ -132,9 +131,9 @@ plant_hydraulics = PlantHydraulics.PlantHydraulicsModel{FT}(;
compartment_surfaces = compartment_surfaces,
compartment_midpoints = compartment_midpoints,
);

ac_canopy = FT(1e3)
energy_model = ClimaLand.Canopy.BigLeafEnergyModel{FT}(
BigLeafEnergyParameters{FT}(FT(1e3)),
BigLeafEnergyParameters{FT}(ac_canopy),
)

canopy = ClimaLand.Canopy.CanopyModel{FT}(;
Expand Down Expand Up @@ -175,10 +174,10 @@ Y.canopy.hydraulics.ϑ_l.:1 .= augmented_liquid_fraction.(ν, S_l_ini[1])
Y.canopy.hydraulics.ϑ_l.:2 .= augmented_liquid_fraction.(ν, S_l_ini[2])



t0 = 150 * 3600 * 24.0
seconds_per_day = 3600 * 24.0
t0 = 150seconds_per_day
N_days = 100
tf = t0 + 3600 * 24 * N_days
tf = t0 + N_days * seconds_per_day
evaluate!(Y.canopy.energy.T, atmos.T, t0)
set_initial_cache! = make_set_initial_cache(canopy)
set_initial_cache!(p, Y, t0);
Expand Down Expand Up @@ -216,8 +215,9 @@ ref_T = [parent(ref_sol.u[k].canopy.energy.T)[1] for k in 1:length(ref_sol.t)]
mean_err = []
p95_err = []
p99_err = []
dts = [225.0, 450.0, 900.0, 1800.0, 3600.0]
dts = [12.0, 24.0, 48.0, 100.0, 225.0, 450.0, 900.0, 1800.0, 3600.0]
for dt in dts
@info dt
saveat = Array(t0:(3 * 3600):tf)
evaluate!(Y.canopy.energy.T, atmos.T, t0)
updateat = Array(t0:(3600 * 3):tf)
Expand Down
4 changes: 2 additions & 2 deletions src/standalone/Vegetation/canopy_energy.jl
Original file line number Diff line number Diff line change
Expand Up @@ -188,8 +188,8 @@ function ClimaLand.make_compute_jacobian(
area_index = p.canopy.hydraulics.area_index
ac_canopy = model.parameters.ac_canopy
earth_param_set = canopy.parameters.earth_param_set
_T_freeze = FT(LP.T_freeze(earth_param_set))
= FT(LP.Stefan(earth_param_set))
_T_freeze = LP.T_freeze(earth_param_set)
= LP.Stefan(earth_param_set)
@. ∂LW_n∂Tc = -2 * 4 ** ϵ_c * Y.canopy.energy.T^3 # ≈ ϵ_ground = 1
@. ∂qc∂Tc = partial_q_sat_partial_T_liq(
p.drivers.P,
Expand Down
10 changes: 6 additions & 4 deletions src/standalone/Vegetation/component_models.jl
Original file line number Diff line number Diff line change
Expand Up @@ -129,8 +129,9 @@ function ClimaLand.make_compute_exp_tendency(
)
function compute_exp_tendency!(dY, Y, p, t)
vars = prognostic_vars(component)
if vars != ()
getproperty(getproperty(dY, name(canopy)), name(component)) .= 0
dY_canopy = getproperty(dY, name(canopy))
if !isempty(vars)
getproperty(dY_canopy, name(component)) .= 0
end
end
return compute_exp_tendency!
Expand All @@ -152,8 +153,9 @@ function ClimaLand.make_compute_imp_tendency(
)
function compute_imp_tendency!(dY, Y, p, t)
vars = prognostic_vars(component)
if vars != ()
getproperty(getproperty(dY, name(canopy)), name(component)) .= 0
dY_canopy = getproperty(dY, name(canopy))
if !isempty(vars)
getproperty(dY_canopy, name(component)) .= 0
end

end
Expand Down

0 comments on commit 7eb5cb7

Please sign in to comment.