reduce allocations in canopy II

src/shared_utilities/drivers.jl

@@ -23,6 +23,7 @@ export AbstractAtmosphericDrivers,
     snow_precipitation,
     vapor_pressure_deficit,
     displacement_height,
+    relative_humidity,
     make_update_drivers
 
 """
@@ -648,6 +649,33 @@ function vapor_pressure_deficit(T_air, P_air, q_air, thermo_params)
     return es - ea
 end
 
+"""
+    relative_humidity(T_air, P_air, q_air, thermo_params)
+
+Computes the vapor pressure deficit for air with temperature T_air,
+pressure P_air, and specific humidity q_air, using thermo_params,
+a Thermodynamics.jl param set.
+"""
+function relative_humidity(
+    T_air::FT,
+    P_air::FT,
+    q_air::FT,
+    thermo_params,
+) where {FT}
+    es = Thermodynamics.saturation_vapor_pressure(
+        thermo_params,
+        T_air,
+        Thermodynamics.Liquid(),
+    )
+    ea = Thermodynamics.partial_pressure_vapor(
+        thermo_params,
+        P_air,
+        Thermodynamics.PhasePartition(q_air),
+    )
+    return es / ea
+end
+
+
 """
     initialize_drivers(a::PrescribedAtmosphere{FT}, coords) where {FT}
 

src/standalone/Vegetation/Canopy.jl

@@ -410,6 +410,7 @@ function ClimaLand.make_update_aux(
         fa = p.canopy.hydraulics.fa
         inc_par = p.canopy.radiative_transfer.inc_par
         inc_nir = p.canopy.radiative_transfer.inc_nir
+        frac_diff = p.canopy.radiative_transfer.frac_diff
 
         # Current atmospheric conditions
         θs = p.drivers.θs
@@ -445,26 +446,21 @@ function ClimaLand.make_update_aux(
             canopy.radiative_transfer.parameters.ϵ_canopy *
             (1 - exp(-(LAI + SAI))) #from CLM 5.0, Tech note 4.20
         RT = canopy.radiative_transfer
-        inc_par .= compute_PAR(RT, canopy.radiation, p, t)
-        inc_nir .= compute_NIR(RT, canopy.radiation, p, t)
+        compute_PAR!(inc_par, RT, canopy.radiation, p, t)
+        compute_NIR!(inc_nir, RT, canopy.radiation, p, t)
         K = extinction_coeff.(G_Function, θs)
-        e_sat =
-            Thermodynamics.saturation_vapor_pressure.(
-                thermo_params,
-                T_air,
-                Ref(Thermodynamics.Liquid()),
-            )
-        e =
-            Thermodynamics.partial_pressure_vapor.(
-                thermo_params,
-                P_air,
-                PhasePartition.(q_air),
-            )
-        rel_hum = e ./ e_sat
         DOY = Dates.dayofyear(
             canopy.atmos.ref_time + Dates.Second(floor(Int64, t)),
         )
-        frac_diff = @. diffuse_fraction(DOY, T_air, p.drivers.SW_d, rel_hum, θs)
+        @. frac_diff = diffuse_fraction(
+            DOY,
+            T_air,
+            P_air,
+            q_air,
+            p.drivers.SW_d,
+            θs,
+            thermo_params,
+        )
 
         compute_absorbances!(
             p,

src/standalone/Vegetation/canopy_parameterizations.jl

@@ -406,11 +406,11 @@ function plant_absorbed_pfd(
 end
 
 """
-    diffuse_fraction(td::FT, T::FT, SW_IN::FT, RH::FT, θs::FT) where {FT}
+    diffuse_fraction(td::FT, T::FT, P, q, SW_IN::FT, θs::FT, thermo_params) where {FT}
 
 Computes the fraction of diffuse radiation (`diff_frac`) as a function
 of the solar zenith angle (`θs`), the total surface incident shortwave radiation (`SW_IN`),
-the air temperature (`T`), the relative humidity (`RH`), and the day of the year
+the air temperature (`T`), air pressure (`P`), specific humidity (`q`), and the day of the year
 (`td`). 
 
 See Appendix A of Braghiere, "Evaluation of turbulent fluxes of CO2, sensible heat,
@@ -427,7 +427,16 @@ can become large negative numbers.
 
 Because of that, we cap the returned value to lie within [0,1].
 """
-function diffuse_fraction(td, T::FT, SW_IN::FT, RH::FT, θs::FT) where {FT}
+function diffuse_fraction(
+    td,
+    T::FT,
+    P::FT,
+    q::FT,
+    SW_IN::FT,
+    θs::FT,
+    thermo_params,
+) where {FT}
+    RH = ClimaLand.relative_humidity(T, P, q, thermo_params)
     k₀ = FT(1370 * (1 + 0.033 * cos(2π * td / 365))) * cos(θs)
     kₜ = SW_IN / k₀
     if kₜ ≤ 0.3

src/standalone/Vegetation/radiation.jl

@@ -128,55 +128,57 @@ function TwoStreamModel{FT}(
 end
 
 """
-    compute_PAR(
+    compute_PAR!(par,
         model::AbstractRadiationModel,
         solar_radiation::ClimaLand.PrescribedRadiativeFluxes,
         p,
         t,
     )
 
-Returns the estimated PAR (W/,m^2) given the input solar radiation
+Updates `par` with the estimated PAR (W/,m^2) given the input solar radiation
 for a radiative transfer model.
 
 The estimated PAR is half of the incident shortwave radiation.
 """
-function compute_PAR(
+function compute_PAR!(
+    par,
     model::AbstractRadiationModel,
     solar_radiation::ClimaLand.PrescribedRadiativeFluxes,
     p,
     t,
 )
-    return @. p.drivers.SW_d / 2
+    @. par = p.drivers.SW_d / 2
 end
 
 """
-    compute_NIR(
+    compute_NIR!(nir,
         model::AbstractRadiationModel,
         solar_radiation::ClimaLand.PrescribedRadiativeFluxes,
         p,
         t,
     )
 
-Returns the estimated NIR (W/m^2) given the input solar radiation
+Update `nir` with the estimated NIR (W/m^2) given the input solar radiation
 for a radiative transfer model.
 
 The estimated PNIR is half of the incident shortwave radiation.
 """
-function compute_NIR(
+function compute_NIR!(
+    nir,
     model::AbstractRadiationModel,
     solar_radiation::ClimaLand.PrescribedRadiativeFluxes,
     p,
     t,
 )
-    return @. p.drivers.SW_d / 2
+    @. nir = p.drivers.SW_d / 2
 end
 
 # Make radiation models broadcastable
 Base.broadcastable(RT::AbstractRadiationModel) = tuple(RT)
 
 ClimaLand.name(model::AbstractRadiationModel) = :radiative_transfer
 ClimaLand.auxiliary_vars(model::Union{BeerLambertModel, TwoStreamModel}) =
-    (:inc_nir, :inc_par, :nir, :par, :LW_n, :SW_n, :ϵ)
+    (:inc_nir, :inc_par, :nir, :par, :LW_n, :SW_n, :ϵ, :frac_diff)
 ClimaLand.auxiliary_types(
     model::Union{BeerLambertModel{FT}, TwoStreamModel{FT}},
 ) where {FT} = (
@@ -187,9 +189,18 @@ ClimaLand.auxiliary_types(
     FT,
     FT,
     FT,
+    FT,
+)
+ClimaLand.auxiliary_domain_names(::Union{BeerLambertModel, TwoStreamModel}) = (
+    :surface,
+    :surface,
+    :surface,
+    :surface,
+    :surface,
+    :surface,
+    :surface,
+    :surface,
 )
-ClimaLand.auxiliary_domain_names(::Union{BeerLambertModel, TwoStreamModel}) =
-    (:surface, :surface, :surface, :surface, :surface, :surface, :surface)
 
 """
     canopy_radiant_energy_fluxes!(p::NamedTuple,