Optimal LAI

ZhouOptimalLAIModel{FT, OLPT <: OptimalLAIParameters{FT}, GD, RDTH, HTH} <: AbstractBiomassModel{FT}

An implementation of the optimal LAI model from Zhou et al. (2025) as a biomass model.

This model computes LAI dynamically based on optimality principles, balancing energy and water constraints. LAI is stored in p.canopy.biomass.area_index.leaf, consistent with PrescribedBiomassModel.

Fields

• parameters: Required parameters for the optimal LAI model
• optimal_lai_inputs: NamedTuple with spatially varying GSL (growing season length in days), A0annual (annual potential GPP in mol CO2 m^-2 yr^-1), precipannual (mol H2O m^-2 yr^-1), vpdgs (Pa), laiinit (initial LAI from MODIS), and f0 (fraction of precip for transpiration), typically created using optimal_lai_initial_conditions.
• SAI: Prescribed stem area index (m^2 m^-2)
• RAI: Prescribed root area index (m^2 m^-2)
• rooting_depth: Rooting depth parameter (m) - a characteristic depth below which 1/e of the root mass lies
• height: Canopy height (m) - can be scalar (uniform) or spatially-varying Field

References

Zhou et al. (2025) "A General Model for the Seasonal to Decadal Dynamics of Leaf Area" Global Change Biology. https://onlinelibrary.wiley.com/doi/pdf/10.1111/gcb.70125

OptimalLAIParameters{FT<:AbstractFloat}

The required parameters for the optimal LAI model based on Zhou et al. (2025).

Water limitation is handled through the f0*P/A0 term following Zhou et al. (2025) Equation 11, where P is annual precipitation and A0 is annual potential GPP.

References

Zhou et al. (2025) "A General Model for the Seasonal to Decadal Dynamics of Leaf Area" Global Change Biology. https://onlinelibrary.wiley.com/doi/pdf/10.1111/gcb.70125

• k: Light extinction coefficient (dimensionless), typically 0.5

• z: Unit cost of constructing and maintaining leaves (mol m^-2 yr^-1), globally fitted as 12.227 mol m^-2 yr^-1

• sigma: Dimensionless parameter representing departure from square-wave LAI dynamics, globally fitted as 0.771

• alpha: Smoothing factor for exponential moving average (dimensionless, 0-1). Set to 0.067 for ~15 days of memory

• f0: Fraction of annual precipitation available for transpiration (dimensionless, 0-1). Following Zhou et al. (2025), f0 = 0.65 at the energy-water limitation transition. In arid regions, f0 can be lower: f0 = 0.65 * exp(-0.604 * ln^2(AI/1.9)) where AI is aridity index. Default value 0.65 assumes optimal water use efficiency.

OptimalLAIParameters{FT}(toml_dict::CP.ParamDict) where {FT}

Creates an OptimalLAIParameters object from a TOML parameter dictionary.

update_optimal_LAI(local_noon_mask, A0_daily, L, k, A0_annual, z, GSL, sigma, alpha, precip_annual, f0, ca_pa, chi, vpd_gs)

Update LAI using the optimal LAI model with precomputed daily and annual potential GPP.

Arguments

• local_noon_mask::FT: Mask (0 or 1) indicating if it's local noon
• A0_daily::FT: Daily potential GPP (mol CO2 m^-2 day^-1), with moisture stress factor beta
• L::FT: Current LAI (m^2 m^-2)
• k::FT: Light extinction coefficient
• A0_annual::FT: Annual potential GPP (mol CO2 m^-2 yr^-1)
• z::FT: Unit cost of constructing and maintaining leaves (mol m^-2 yr^-1)
• GSL::FT: Growing season length (days)
• sigma::FT: Dimensionless parameter for LAI dynamics
• alpha::FT: Smoothing factor for exponential moving average (~15-day memory)
• precip_annual::FT: Mean annual precipitation (mol H2O m^-2 yr^-1)
• f0::FT: Fraction of precipitation available for transpiration (dimensionless)
• ca_pa::FT: Ambient CO2 partial pressure (Pa)
• chi::FT: Optimal ratio of intercellular to ambient CO2 (dimensionless)
• vpd_gs::FT: Mean vapor pressure deficit during growing season (Pa)

Returns

Updated LAI value.

Notes

Following Zhou et al. (2025):

• A0_daily uses moisture stress factor beta to drive daily LAI dynamics
• A0annual uses actual moisture stress factor beta for LAImax computation
• Water limitation enters LAI_max through the f0P/A0 * (ca(1-chi))/(1.6D) term (Equation 11)

compute_A0_daily(is_c3, parameters, constants, T_air, P_air, VPD, ca, PPFD, βm)

Compute daily potential GPP (A0) used in the optimal LAI model (Zhou et al. 2025).

This function computes the potential GPP assuming fAPAR = 1 (full light absorption), which represents the maximum carbon assimilation possible under given environmental conditions. The result is in units of kg C m^-2 s^-1.

Arguments

• is_c3::FT: Photosynthesis mechanism (1 for C3, 0 for C4)
• parameters: PModelParameters containing cstar, beta, etc.
• constants: PModelConstants containing physical constants
• T_air::FT: Air temperature (K)
• P_air::FT: Atmospheric pressure (Pa)
• VPD::FT: Vapor pressure deficit (Pa)
• ca::FT: Ambient CO2 concentration (mol/mol)
• PPFD::FT: Photosynthetic photon flux density (mol photons m^-2 s^-1)
• βm::FT: Soil moisture stress factor (dimensionless, 0-1)

Returns

• A0_daily::FT: Potential GPP with fAPAR=1 (kg C m^-2 s^-1)

Notes

Used by the optimal LAI model to drive LAI dynamics. The moisture stress factor βm is included to capture water limitation effects on photosynthetic capacity. Air temperature is used (rather than canopy temperature) because A0 represents potential GPP under reference conditions, independent of energy balance feedbacks.

compute_LAI(LAI_prev, L_steady, alpha, local_noon_mask)

Compute updated LAI using exponential weighted moving average to represent lag between carbon allocation and steady-state LAI.

This implements Equation 16 from Zhou et al. (2025). The exponential moving average represents the time lag (days to months) for photosynthate allocation to leaves and leaf development.

Arguments

• LAI_prev::FT: LAI from previous time step (m^2 m^-2, dimensionless)
• L_steady::FT: Current steady-state LAI (m^2 m^-2, dimensionless), from compute_steady_state_LAI()
• alpha::FT: Smoothing factor (dimensionless, 0-1). Set to 0.067 for ~15 days of memory, meaning LAI[t] = 0.067 * L_steady[t] + 0.933 * LAI[t-1]
• local_noon_mask::FT: A mask (0 or 1) indicating whether the current time is within the local noon window.

Returns

• LAI_new::FT: Updated leaf area index (m^2 m^-2, dimensionless). Always >= 0.

Notes

The parameter alpha controls the response time:

• alpha = 0.067 ~ 15 days of memory (paper default)
• alpha = 0.1 ~ 10 days of memory (faster response)
• alpha = 0.033 ~ 30 days of memory (slower response)

The time scale tau ~ 1/alpha days.

References

Zhou et al. (2025) Global Change Biology, Equation 16

make_OptimalLAI_callback(::Type{FT}, t0::ITime, dt, canopy; longitude) where {FT <: AbstractFloat}

This constructs an IntervalBasedCallback for the optimal LAI model that:

1. Computes and accumulates potential GPP (A0) at each timestep
2. Updates LAI using an exponential moving average at local noon
3. Tracks daily and annual A0 sums

We check for local noon using the provided longitude every dt. The time of local noon is expressed in seconds UTC and neglects the effects of obliquity and eccentricity, so it is constant throughout the year.

Arguments

• FT: The floating-point type used in the model (e.g., Float32, Float64).
• t0: ITime, with epoch in UTC.
• dt: timestep
• canopy: the canopy object containing the optimal LAI model parameters.
• longitude: optional longitude in degrees for local noon calculation (default is nothing, which means that it will be inferred from the canopy domain).

Notes

• Daily A0 is computed with fAPAR=1 and moisture stress factor beta - drives L_steady
• Annual A0 is computed with fAPAR=1 and actual moisture stress factor beta - used for LAI_max
• Water limitation enters LAI_max through the f0P/A0 * (ca(1-chi))/(1.6D) term (Equation 11)
• Daily A0 is accumulated over each day and finalized at local noon
• Annual A0 is accumulated and reset on January 1
• GSL (Growing Season Length) is read from p.canopy.biomass.GSL, which supports spatially varying values initialized via sethistoricalcache!.

call_update_optimal_LAI(p, Y, t, current_date; canopy, dt, local_noon)

Updates LAI and accumulates potential GPP (A0) at each timestep.

Every timestep: accumulates instantaneous potential GPP into the daily accumulator. At local noon: finalizes daily A0, adds it to annual accumulator, updates LAI. Every 365 days: finalizes annual A0 from the accumulator.

Uses air temperature (not canopy temperature) for A0 computation, since canopy temperature includes energy balance feedbacks that should not affect potential GPP.

GSL (Growing Season Length) is read from p.canopy.biomass.GSL, which supports spatially varying values initialized via sethistoricalcache!.

compute_L_max(Ao_annual, k, z, precip_annual, f0, ca_pa, chi, vpd_gs)

Compute seasonal maximum leaf area index (LAI_max) based on annual potential GPP and water availability, following Zhou et al. (2025) Equation 11.

LAI_max is determined by the minimum of energy-limited and water-limited fAPAR:

• Energy-limited: fAPAR_energy = 1 - z/(k*A0)
• Water-limited: fAPAR_water = f0P/A0 * (ca(1-chi))/(1.6D)

Arguments

• Ao_annual::FT: Annual total potential GPP (mol CO2 m^-2 yr^-1).
• k::FT: Light extinction coefficient (dimensionless), typically 0.5
• z::FT: Unit cost of constructing and maintaining leaves (mol m^-2 yr^-1), 12.227
• precip_annual::FT: Mean annual precipitation (mol H2O m^-2 yr^-1)
• f0::FT: Fraction of precipitation available for transpiration (dimensionless), 0.65
• ca_pa::FT: Ambient CO2 partial pressure (Pa), typically ~40 Pa at 400 ppm
• chi::FT: Optimal ratio of intercellular to ambient CO2 (dimensionless), typically 0.7-0.8
• vpd_gs::FT: Mean vapor pressure deficit during growing season (Pa)

Returns

• LAI_max::FT: Seasonal maximum leaf area index (m^2 m^-2)

Notes

Following Zhou et al. (2025) Equation 11:

fAPAR_max = min{1 - z/(k*A0), f0*P/A0 * (ca(1-chi))/(1.6*D)}

The first term is energy-limited (carbon gain vs leaf cost trade-off). The second term is water-limited (precipitation constrains transpiration, scaled by intrinsic water use efficiency iWUE = ca(1-chi)/(1.6*D)).

The iWUE factor converts water flux to carbon flux:

• ca(1-chi): CO2 drawdown from ambient to intercellular (Pa)
• 1.6*D: VPD adjusted for CO2/H2O diffusivity ratio (Pa)

References

Zhou et al. (2025) Global Change Biology, Equation 11

compute_m(GSL, LAI_max, Ao_annual, sigma, k)

Compute the parameter m, which represents the ratio of steady-state LAI to steady-state GPP.

This implements Equation 20 from Zhou et al. (2025). The parameter m quantifies the relationship between LAI and GPP dynamics, representing the extent to which seasonal LAI dynamics depart from a "square wave" (where maximum LAI would be maintained throughout the growing season).

Arguments

• GSL::FT: Growing season length (days). Defined as the length of continuous period above 0C longer than 5 days.
• LAI_max::FT: Seasonal maximum leaf area index (m^2 m^-2, dimensionless)
• Ao_annual::FT: Annual total potential GPP (mol m^-2 yr^-1). This is the integral of daily A0 over the year.
• sigma::FT: Dimensionless parameter representing departure from square-wave LAI dynamics. Globally fitted as sigma = 0.771
• k::FT: Light extinction coefficient (dimensionless)

Returns

• m::FT: Parameter relating steady-state LAI to steady-state GPP (dimensionless, units work out as: days * m^2 m^-2 / (mol m^-2 yr^-1 * dimensionless) with implicit conversion)

References

Zhou et al. (2025) Global Change Biology, Equation 20

compute_steady_state_LAI(Ao_daily, m, k, LAI_max)

Compute steady-state LAI from daily potential GPP using the Lambert W function solution.

This implements Equations 13-15 from Zhou et al. (2025). The steady-state LAI (L_s) is the LAI that would be in equilibrium with GPP if weather conditions were held constant. Given daily meteorological conditions, this is computed on a daily basis.

Arguments

• Ao_daily::FT: Daily potential GPP (mol m^-2 day^-1). This is the GPP that would be achieved if fAPAR = 1, calculated from LUE * PPFD.
• m::FT: Parameter relating steady-state LAI to steady-state GPP (dimensionless), from compute_m()
• k::FT: Light extinction coefficient (dimensionless), typically 0.5
• LAI_max::FT: Seasonal maximum LAI constraint (m^2 m^-2, dimensionless)

Returns

• L_steady::FT: Steady-state leaf area index (m^2 m^-2, dimensionless). Always >= 0.

Notes

The solution uses the Lambert W0 function: Ls = min{mu + (1/k)W0[-kmuexp(-k*mu)], LAImax} where mu = m * A0. The result is constrained to be non-negative and below LAI_max.

References

Zhou et al. (2025) Global Change Biology, Equations 13-15

lambertw0(x::T; maxiter::Int = 8) where {T<:AbstractFloat}

Compute the principal branch (W0) of the Lambert W function for x in [-1/e, Inf).

This is a GPU-device-friendly implementation using a fixed number of Halley iterations. The Lambert W function satisfies W(x)*exp(W(x)) = x.

Arguments

• x::T: Input value, must be >= -1/e ~ -0.36788
• maxiter::Int: Maximum number of Halley iterations (default: 8; Halley's method has cubic convergence, so 8 is generous)

Returns

• W::T: Lambert W0(x), the principal branch value, or NaN for invalid inputs

Algorithm

Uses Halley's method with a fixed number of iterations for GPU compatibility:

• No dynamic memory allocation
• No conditional breaks (runs all iterations)
• Broadcastable for use with CuArrays: lambertw0.(cuarray)

Device Compatibility

This implementation is designed to work on both CPU and GPU:

• All operations are scalar and supported on CUDA.jl
• No array allocations or dynamic loops
• Type-generic over AbstractFloat (Float32, Float64)

References

Corless et al. (1996) "On the Lambert W function"