Process-based upgrades to the WRF multi-layer green-roof scheme (WRF-MLGR v2.0) and evaluation against field observations

Green roofs can moderate urban heat by increasing latent heat flux and reducing sensible heat flux. However, capturing these effects in models depends on accurate representation of key green-roof processes, including substrate heat and moisture transport, soil-vegetation-atmosphere energy and moisture exchanges, interactions with the underlying roof, and drainage. Here we introduce targeted, process-based updates to the green-roof scheme (hereafter, MLGR) within the multi-layer urban canopy model (BEP-BEM) in the WRF mesoscale model to address key limitations in the original formulation. The updates include a non-linear dependence of soil thermal conductivity on moisture, vegetation-modulated surface thermal conductivity, explicit soil-surface evaporation, multi-layer root water uptake for transpiration, and canopy interception with evaporation and dew formation. We evaluate the original and modified MLGR schemes using hourly observations from an extensive sedum roof in London, Canada, for “summer” (1 July–31 August 2014) and “fall” (1 September–31 October 2014) periods. We also analyze 11–18 October 2025, when green roof modules were placed directly on the roof deck – which corresponds to the model's lower boundary assumption. Following implementation of the process-based improvements, model–measurement agreement for the conductive heat flux is markedly improved: RMSE is reduced from 105.9 to 24.0 W m⁻² in summer and from 94.2 to 24.0 W m⁻² in fall and the model produces more realistic overall green roof energy partitioning. The modified model better captures post-rain increases in latent heat flux (Q_E) and improves the timing and magnitude of daytime turbulent latent and sensible heat flux peaks (Q_E and Q_H). Drainage is reduced relative to the original scheme; however, it remains slightly underestimated during summer and slightly overestimated during fall, and biases persist in surface temperature (warm during the day and cool at night) and in the magnitude and variability of Q_E. Overall, the revised MLGR physics improves surface-flux realism, and future development should focus on developing a more realistic vegetation canopy submodule.

1 Introduction

Rapid urbanization has intensified heat exposure as cities alter the surface energy balance through suppressed evapotranspiration, enhanced heat storage in impervious materials, and complex radiative geometry (Oke, 1982; Grimmond and Oke, 1999; Oke et al., 2017), producing urban heat islands (UHIs) and compounding heat-wave risks (Arnfield, 2003; Li and Bou-Zeid, 2013). Green roofs, engineered vegetated systems atop buildings, have emerged as a practical intervention capable of reducing air temperatures because they repartition turbulent heat exchange, decreasing sensible heat flux (Q_H) and enhancing (Q_E) (Oberndorfer et al., 2007; Berardi et al., 2014; Santamouris, 2014; Tabares-Velasco and Srebric, 2012), act as a sink for air pollutants by trapping particulate matter on vegetated surfaces (Rowe, 2011), and provide hydrologic co-benefits when deployed at scale (Berndtsson, 2010; Stovin, 2010; Li and Babcock, 2014). The magnitude and persistence of these benefits depend on climate, roof composition, vegetation traits, and especially substrate moisture (Sailor et al., 2008; Tabares-Velasco and Srebric, 2012; Tabares-Velasco et al., 2012) and irrigation management, which together control the partitioning of net radiation into heat storage, sensible heat and evapotranspiration (Heusinger et al., 2018).

At the building scale, the EnergyPlus “EcoRoof” module introduced a physically based coupling between a plant canopy and an engineered soil substrate, while highlighting the need for more mechanistic hydrology and canopy microclimate formulations (Sailor, 2008; Sailor et al., 2008; Heusinger et al., 2018). Subsequent work showed that sensible-heat reduction depends strongly on water availability, with post-rain evaporation playing a key role in shifting energy partitioning toward latent heat flux and enhanced cooling (Heusinger et al., 2018). At the city scale, urban canopy parameterizations (UCPs) embedded in mesoscale models, single-layer schemes (Kusaka et al., 2001; Yang et al., 2015), and multi-layer approaches such as TEB (Masson, 2000; Hamdi and Masson, 2008) and BEP (Martilli et al., 2002; Salamanca et al., 2011) represent radiative geometry, turbulent exchange, storage, and building–atmosphere feedbacks, enabling scenario testing of rooftop interventions at neighbourhood and city scales (Zonato et al., 2021). These coupled model frameworks have been used to assess how green roofs influence urban temperature and heat mitigation, and how these changes affect building energy demand, in both idealized and real-city applications (Li et al., 2014; Sun et al., 2016; de Munck et al., 2018; Jiang et al., 2025).

Despite this progress, there are two persistent gaps that limit fidelity at urban scales. First, the hydrology and thermophysics of thin, engineered substrates are often simplified in mesoscale implementations relative to detailed building-scale models (Yang et al., 2015; Ribeiro et al., 2021). Second, ensuring strict conservation of energy and water across multi-layer substrate and canopy elements remains challenging in offline tests and coupled WRF runs (Tsiringakis et al., 2019; Lipson et al., 2024). Addressing these issues is essential if we are to translate building-scale realism into robust city-scale modelling of green roof scenarios.

This study advances the state of green-roof representation in the WRF multi-layer urban canopy model (BEP+BEM). A recent process-level evaluation of the WRF-BEP+BEM green roof scheme (Zonato et al., 2021) has identified structural deficiencies in its representation of vegetation and soil processes, including the absence of soil evaporation and lack of canopy radiative buffering, which lead to systematic biases in energy partitioning, pronounced daytime warm biases in surface temperature (T_s), and failure to reproduce the observed diurnal phase and magnitude of Q_E and Q_G (Martinez Mendoza et al., 2026). These findings motivate the targeted physical upgrades and evaluation presented here. We further develop and evaluate this model (Zonato et al., 2021) by (i) applying a moisture-dependent thermal conductivity that improves soil heat storage and release during high soil water content, (ii) modulating effective surface thermal conductivity as a function of vegetation cover to better capture diurnal heat storage, (iii) adding direct soil evaporation alongside plant transpiration to represent total evapotranspiration, (iv) allowing plant water uptake from all substrate layers rather than the top layer alone, and (v) including interception storage, evaporation, and dew formation on leaves. Together, these changes target the two gaps above: improving flux partitioning and the phase and magnitude of diurnal cycles, including the timing of Q_H, Q_E, and Q_G, while remaining computationally compatible with WRF's BEP+BEM framework and conserving energy and mass.

Our objectives are twofold. First, we document the implementation of new formulations in the green-roof scheme implemented in WRF's multi-layer urban canopy model (BEP+BEM). Second, we evaluate the modified model against observations from a monitored extensive green roof in London, Canada, assessing surface temperature, turbulent latent heat flux, storage heat flux, substrate moisture, and runoff.

2 Methods

2.1 Experimental design

2.1.1 Site and simulation periods

Field observations were collected from a green-roof test array on the rooftop of the Talbot College building at Western University (London, Ontario, Canada). The analysis focuses on the instrumented central portion of the array, which is elevated above the surrounding roof surface to accommodate lysimeters and drainage units. The monitored modules were planted with Sedum spurium and had a total substrate depth of 0.15 m.

This paper evaluates the WRF multi-layer green-roof scheme using observations from 2014 and a focused period in 2025. We analyze two continuous periods in 2014: “summer” (1 July–31 August 2014) and “fall” (1 September–31 October 2014). We also examine 11–18 October 2025, when the heat-flux module was moved from the upper array to a layer of styrofoam insulation above the roof deck, providing a better physical match to the modeled roof–substrate setup for clearer interpretation of measured ground heat flux.

Raw measurements were quality controlled, and time periods affected by instrument malfunction, maintenance, calibration events, or missing data were excluded. For this paper, both the observational time series and the corresponding model output were processed to hourly resolution (hourly means for state variables and fluxes, and hourly totals for precipitation). Further details on the site instrumentation, observational processing, and quality-control procedures are provided in Martinez Mendoza et al. (2026); readers are referred there for a complete description of the observational dataset.

2.1.2 Simulation development

The simulations were performed using the multi-layer green roof scheme from WRF v4.3.3 with BEP+BEM in offline mode. The full atmospheric WRF model was not run. Instead, the green roof model was driven directly by observed meteorological forcing.

Atmospheric forcing data were prescribed from observations collected on the roof at, or near, the green-roof array. The forcing time series includes incoming shortwave radiation (W m⁻²), incoming longwave radiation (W m⁻²), near-surface air temperature (°C), relative humidity (%), wind speed (m s⁻¹), precipitation (mm; hourly accumulation), and surface air pressure (Pa). Pressure was taken from the nearby London A meteorological station at London International Airport (WMO 71623).

Initial conditions for prognostic substrate temperature and water-content variables were prescribed using the available observations at the initialization time. Geometric and material parameters (e.g., layer depths) were prescribed from site measurements and prior characterization (Martinez Mendoza et al., 2026; Perelli, 2014).

2.1.3 Evaluation metrics

Model performance was assessed against measurements of latent heat flux Q_E (W m⁻²), surface ground heat flux Q_G (W m⁻²), surface temperature T_s (°C), and substrate volumetric water content θ (m³ m⁻³).

Performance statistics reported in Table 1 include mean absolute error (MAE), root-mean-square error (RMSE), systematic RMSE (RMSEs), unsystematic RMSE (RMSEu), and the refined index of agreement (d) following Willmott et al. (2012).

Table 1 Performance statistics for the original and modified WRF multi-layer green-roof models during summer and fall. Q_E and Q_G in W m⁻²; T_s in °C; θ in m³ m⁻³.

Note: MAE = mean absolute error; RMSEu = unbiased RMSE; RMSEs = systematic RMSE; d = refined index of agreement per Willmott et al. (2012).

The coefficient of determination (R²) shown in the scatter plots was calculated using the residual-based form $R^{\mathrm{2}}=\mathrm{1}-\sum (M_i-O_i{)}^{\mathrm{2}}/\sum (O_i-\stackrel{\mathrm{‾}}{O}{)}^{\mathrm{2}}$ (Onyutha, 2022), where M_i and O_i are the modeled and observed values, respectively, and $\stackrel{\mathrm{‾}}{O}$ is the mean of the observations. This metric evaluates model performance relative to the observed mean as a reference predictor; values below zero indicate that the model performs worse than simply using the observed mean (Nash and Sutcliffe, 1970; Onyutha, 2022).

2.2 Original WRF multi-layer urban canopy model green roof scheme

Following the formulations of de Munck et al. (2013), as summarized by Zonato et al. (2021), the offline WRF-MLGR configuration contains 10 thermal layers in total. Six layers represent the green-roof substrate and have a combined depth of 0.15 m, consistent with the observed substrate depth of the experimental modules, while the remaining four layers represent approximately 0.15 m of underlying roof construction. Thus, the complete modeled thermal column is approximately 0.30 m deep, but only the upper 0.15 m represents the green-roof substrate. Temperature is calculated in all thermal layers, while volumetric water content is calculated only for the six soil layers. Layer-resolved soil temperature and moisture diagnostics are provided in Figs. S2 and S3.

During the offline evaluation, the indoor temperature is not simulated. Instead, the model uses a fixed temperature of 23 °C at the bottom of the green-roof system. representing a typical indoor/building temperature and providing the lower thermal boundary condition for the offline simulations. This temperature affects the amount of heat transferred between the green roof and the underlying roof. This model computes the roof's heat and water balances by accounting for net radiation, water added by rainfall and irrigation, evapotranspiration from the vegetation, convective heat exchange with the atmosphere, and the vertical diffusion of heat and moisture within the soil layers. In the original scheme, the model makes an important simplification: it does not distinguish a separate vegetation canopy in terms of radiation exchange, or a distinct canopy air layer for heat and moisture exchange. Instead, the green roof has a single surface and associated temperature which is effectively a soil surface with the ability to transpire (and with substantive thermal admittance, unlike a vegetation canopy). Therefore, the reported surface temperature T_s in the current WRF-MLGR formulation should be interpreted as a bulk or aggregated green-roof skin temperature, rather than as a separately resolved soil-surface temperature, vegetation-canopy temperature, or canopy-air temperature. It is this surface that exchanges heat and moisture with the atmosphere and conducts heat downward into the roof/substrate column. Heat transfer within the layers is then calculated based on Fourier diffusion, and for the soil-based roof layers the model uses a moisture-dependent thermal diffusivity (Zonato et al., 2021):

$$\begin{array}{}\text{(1)}& \mathit{\lambda }=\left\{\begin{array}{ll}\frac{e^{-\left({\log }_{\mathrm{10}}\right|\mathrm{\Psi }|+\mathrm{2.7})}}{C_{\mathrm{s}}}\cdot \mathrm{4.186}\times {\mathrm{10}}^{\mathrm{7}},& \text{if }{\log }_{\mathrm{10}}\left|\mathrm{\Psi }\right|\le \mathrm{5.1},\\ \frac{\mathrm{4.1}\times {\mathrm{10}}^{-\mathrm{5}}}{C_{\mathrm{s}}}\cdot \mathrm{4.186}\times {\mathrm{10}}^{\mathrm{7}},& \text{if }{\log }_{\mathrm{10}}\left|\mathrm{\Psi }\right|>\mathrm{5.1},\end{array}\right.\end{array}$$

where Ψ is the moisture potential (matric suction; negative under unsaturated conditions), and $C_{\mathrm{s}}=(\mathrm{1}-\mathit{\theta })C_{\mathrm{d}}+\mathit{\theta }{C}_{\mathrm{w}}$ is the volumetric heat capacity of wet soil ($\mathrm{J}{\mathrm{m}}^{-\mathrm{3}}{\mathrm{K}}^{-\mathrm{1}}$). It is computed as a weighted average of the volumetric heat capacity of dry soil (C_d) and water (C_w), with θ the volumetric water content Zonato et al. (2021).

For hydrology, rainfall (and optional irrigation) is added at the green roof (i.e., soil) surface, and water moves only vertically through the substrate in a one-dimensional column; therefore, lateral subsurface flow is not represented. Any excess water that reaches the drainage layer is removed as drainage. In the baseline scheme, moisture loss to the atmosphere is represented through plant transpiration (controlled by canopy and atmospheric resistances), and the model does not include a separate explicit term for direct soil-surface evaporation. The scheme also does not explicitly track a separate canopy-interception water store.

As discussed by Martinez Mendoza et al. (2026), the original green roof scheme in the WRF multi-layer urban canopy model exhibits important limitations, particularly in its representation of soil heat flux (Q_g). The model does not correctly reproduce the timing and magnitude of Q_g, and these errors propagate through the surface energy balance, leading to inaccuracies in both the strength and timing of latent (Q_e) and sensible (Q_h) heat fluxes. This issue is traced to the fact that there is no radiative canopy, and the soil layer (which has high thermal admittance relative to vegetation) is the radiatively active surface in the model instead of the vegetation. In the next section we detail several updates to this original model version that correct this problem as well as other model issues, and also add additional functionality and realism. To clarify the structural and process-level differences between the original and modified WRF multi-layer green-roof schemes, Fig. 1 provides a schematic comparison of the two model formulations. The figure highlights the main layer structure, the key heat and water flux pathways, and the locations of the five process-based modifications introduced in Sect. 2.3.

Figure 1 Schematic comparison of the original and modified WRF multi-layer green-roof schemes.

2.3 Modifications to the WRF multi-layer urban canopy model green roof scheme

To improve the representation of heat and water exchange in the WRF multi-layer green-roof scheme, we implemented five physically based modifications. The two thermal modifications mainly change how heat moves through the green roof substrate and how much energy is stored as Q_G. The revised soil thermal conductivity provides a more realistic representation of heat transfer as soil moisture changes, while the vegetation-dependent surface conductivity reduces excessive heat storage at the soil surface to represent the effect of vegetation shading on heat storage.

The three hydrological modifications mainly affect evapotranspiration and soil water. Soil evaporation allows water to evaporate directly from the substrate surface. Multilayer root uptake allows plants to take up water from several soil layers rather than only from the top layer. Canopy interception allows rainfall to be temporarily stored on the vegetation, evaporated back to the atmosphere, and also allows dew formation. Together, these modifications change Q_G and Q_E and therefore also affect how the remaining available energy is partitioned into Q_H.

2.3.1 Revised soil thermal conductivity (“Ksoil”)

In the original scheme, the substrate soil thermal conductivity, k, was defined as a simple linear function of volumetric water content. This approach tends to overestimate heat conductivity when the soil is very wet, which can make the simulated ground heat flux too large and reduce the diurnal range of surface temperature under wet conditions. To address this, we replaced the linear relationship with a Johansen-type formulation (Peters-Lidard et al., 1998). In this method, the effective soil thermal conductivity, k, is calculated as a weighted combination of the dry and saturated conductivities using the Kersten number, K_e:

$$\begin{array}{}\text{(2)}& k=k_{\text{dry}}+K_e(k_{\text{sat}}-k_{\text{dry}}).\end{array}$$

Here k_dry is the dry soil thermal conductivity, k_sat is the saturated soil thermal conductivity, and K_e is a dimensionless function of the degree of saturation. For dry soils we follow Peters-Lidard et al. (1998), who recommend the following semi-empirical expression for k_dry as a function of dry bulk density ρ_b (kg m⁻³) and mineral density ρ_s (kg m⁻³):

$$\begin{array}{}\text{(3)}& k_{\text{dry}}={\displaystyle \frac{(\mathrm{0.135}{\mathit{\rho }}_{\mathrm{b}}+\mathrm{64.7})}{({\mathit{\rho }}_{\mathrm{s}}-\mathrm{0.947}{\mathit{\rho }}_{\mathrm{b}})}}\end{array}$$

with ρ_s=2700 kg m⁻³ for mineral soils. The saturated thermal conductivity k_sat is calculated from the volume fractions and thermal conductivities of solids, water, air, and (when present) ice. The Kersten number K_e is parameterized as a smooth function of liquid water saturation, with different coefficients for coarse and fine substrates, which reduces the sensitivity of k at very low moisture contents. Overall, the updated formulation gives a more realistic shift from dry to wet conditions and reduces the unrealistically high soil heat conduction produced by the original linear method.

2.3.2 Vegetation-dependent surface thermal conductivity (“Heat Storage”)

In the original green-roof scheme, the model does not include an explicit vegetation canopy, so it does not represent vegetation cover or the shading effect of plants. This means the soil surface receives the full incoming radiation, which can lead to unrealistically large surface heat storage, particularly in the morning, with associated impacts on morning sensible and latent heat fluxes. To make the simulated surface heat storage more realistic, we added a vegetation dependent effective surface conductivity that acts as a simplified representation of the vegetation effect on surface heat storage. Following Peters-Lidard et al. (1997) and Ek et al. (2003), we used a vegetation-dependent effective thermal conductivity at the top soil layer, k_veg, which scales the bare-soil conductivity k_s1 by an exponential function of vegetation cover σ_f:

$$\begin{array}{}\text{(4)}& k_{\text{veg}}=k_{\text{s1}}\exp (-{\mathit{\beta }}_{\text{veg}}{\mathit{\sigma }}_{\mathrm{f}}),\end{array}$$

where σ_f is the fractional vegetative coverage of the roof and β_veg is an empirical coefficient. Consistent with offline Noah LSM tests (Ek et al., 2003), β_veg was set to 2.0. Corresponding effective thermal conductivities were used within the upper 10 cm of soil for the calculation of heat conduction.

To quantify the effect of the revised thermal treatment, we compared the time-varying average soil thermal conductivity produced by the original and modified green-roof schemes during summer 2014 (Fig. S1 in the Supplement). The original formulation produced unrealistically large and highly variable conductivity values, with peaks exceeding 49 $\mathrm{W}{\mathrm{m}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$, indicating excessive conductive heat transfer under wet conditions. In contrast, the modified formulation produced conductivity values mainly within the range 0.32–0.98 $\mathrm{W}{\mathrm{m}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$, which is physically more reasonable for engineered green-roof growing media, for which thermal conductivity has been shown to vary strongly with substrate moisture content and composition (Barozzi et al., 2017; Shao et al., 2021).

2.3.3 Addition of soil surface evaporation (“Soil Evap”)

In the original WRF multi-layer green-roof model, latent heat flux was derived only from plant transpiration, and direct evaporation from the soil surface was not included. This approach is not applicable to scenarios where green roof vegetation fraction drops significantly below 1.0. We used the bulk-transfer approach of Viterbo and Beljaars (1995) to quantify evaporation from the soil/substrate surface. The evaporation rate from the top substrate layer, QE_E ($\mathrm{kg}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}$), is given by:

$$\begin{array}{}\text{(5)}& {\text{QE}}_{\mathrm{E}}={\displaystyle \frac{\mathit{\rho }}{r_{\mathrm{a}}}}[q_{\mathrm{L}}-\mathit{\alpha }({\mathit{\theta }}_{\mathrm{1}}\left)q_{\text{sat}}\right(T_{\text{sk}},p_{\mathrm{s}}\left)\right].\end{array}$$

where ρ is air density, r_a is the aerodynamic resistance between the surface and the lowest atmospheric model level, q_L is the specific humidity at that level, and q_sat(T_sk,p_s) is the saturation specific humidity at the surface skin temperature T_sk and surface pressure p_s. α(θ₁) is a soil-wetness factor that scales the soil-surface vapor pressure (and thus soil evaporation) from dry to saturated conditions. It is computed from the top-layer water content θ₁, with α≈1 for wet soil and smaller α for drier soil defined as:

$$\begin{array}{}\text{(6)}& \begin{array}{rlrl}& \mathit{\alpha }\left({\mathit{\theta }}_{\mathrm{1}}\right)=\mathrm{0.5}\left[\mathrm{1}-\cos \left({\displaystyle \frac{\mathit{\pi }{\mathit{\theta }}_{\mathrm{1}}}{\left({\mathit{\theta }}_{\text{cap}}\right)}}\right)\right],& & {\mathit{\theta }}_{\mathrm{1}}<{\mathit{\theta }}_{\text{cap}}\\ & \mathit{\alpha }\left({\mathit{\theta }}_{\mathrm{1}}\right)=\mathrm{1},& & {\mathit{\theta }}_{\mathrm{1}}\ge {\mathit{\theta }}_{\text{cap}},\end{array}\end{array}$$

where θ₁ is the volumetric water content of the upper soil layer and θ_cap the field-capacity value, Viterbo and Beljaars (1995).

2.3.4 Multi-layer transpiration and root water uptake (“Root Uptake”)

In the original scheme, plant transpiration draws water only from the thin uppermost soil layer, which can make the surface layer dry too quickly. As a result, when the top layer becomes dry, the model limits transpiration, leading to reduced transpiration and an underestimation of latent heat flux. To address this issue, we introduced a multi-layer root-uptake approach that distributes the transpiration demand across the soil layers based on their available water, assuming that plant roots have good access to all soil layers. Total transpiration E_t ($\mathrm{kg}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}$) is still computed using the original canopy-resistance formulation, but the fraction taken from layer i is:

$$\begin{array}{}\text{(7)}& w_i={\displaystyle \frac{d{z}_i}{{\sum }_{j=\mathrm{1}}^{N}d{z}_j}},\end{array}$$

where dz_i is the thickness of layer i, and the sum is over the N soil layers from which transpiration is allowed to draw water. The transpiration extraction from layer i is then:

$$\begin{array}{}\text{(8)}& E_{\mathrm{t},i}=w_i{E}_{\mathrm{t}}.\end{array}$$

2.3.5 Interception storage, evaporation, and dew formation on leaves (“Intercepted Water”)

Finally, we implemented a canopy interception sub-model following Heusinger et al. (2018), whose formulation is based on Cammalleri et al. (2010). This sub-model tracks the canopy water store on the leaves and accounts for both evaporation of intercepted water and dew formation (condensation) on the leaf surfaces, depending on atmospheric conditions. The maximum amount of water that can remain on the leaves (interception storage capacity) is assumed to be proportional to leaf area index (LAI):

$$\begin{array}{}\text{(9)}& S_{\text{max}}=F_{\text{Smax}}\times \text{LAI}.\end{array}$$

Here, F_Smax (m) is a parameter that depends on leaf inclination and morphological properties. A value of F_smax=0.2 mm was adopted based on measurements conducted for herbaceous species on a mountain meadow in Austria (Wohlfahrt et al., 2006). At each time step, a fraction of the incident precipitation P is intercepted in proportion to vegetation cover σ_f:

$$\begin{array}{}\text{(10)}& I={\mathit{\sigma }}_{\mathrm{f}}\times P.\end{array}$$

The intercepted water depth on the canopy, S, is updated by adding I. If S exceeds S_max, the excess water drips to the top soil layer as D:

$$\begin{array}{}\text{(11)}& \begin{array}{rlrl}& D={\displaystyle \frac{(S-S_{\text{max}})}{\mathrm{\Delta }t}},& & \text{if }S>S_{\text{max}}\\ & D=\mathrm{0},& & \text{otherwise},\end{array}\end{array}$$

where Δt is the model time step. Water stored on the leaves exchanges moisture with the canopy air through a bulk-transfer formulation of the latent heat flux from intercepted water, L_I (W m⁻²):

$$\begin{array}{}\text{(12)}& L_{\mathrm{I}}=\text{LAI}\cdot C_{\mathrm{f}}\cdot L_{\mathrm{v}}\cdot u_{\text{af}}\cdot {\mathit{\rho }}_{\text{af}}\cdot (q_{\text{af}}-q_{\text{f,sat}})\cdot r.\end{array}$$

In this expression, C_f is a bulk transfer coefficient, L_v is the latent heat of vaporization, u_af is wind speed at canopy height, ρ_af is air density at canopy height, q_af is the specific humidity in the canopy air space, and q_f,sat is the saturation specific humidity at foliage temperature. The factor r accounts for the fullness of the interception store:

$$\begin{array}{}\text{(13)}& r={\left({\displaystyle \frac{S}{S_{\text{max}}}}\right)}^{\frac{\mathrm{2}}{\mathrm{3}}}.\end{array}$$

When intercepted water is present on the canopy, a fraction of the leaf area is assumed to be water-covered (r), so transpiration is limited to the remaining dry fraction (1−r); as the interception store evaporates and r decreases, transpiration progressively returns toward its dry-canopy rate. So, transpiration (Zonato et al., 2021) is defined as:

$$\begin{array}{}\text{(14)}& L_{\mathrm{t}}=\text{LAI}\cdot {\mathit{\rho }}_{\text{af}}\cdot L_{\mathrm{v}}\cdot {\displaystyle \frac{\left(\right(q_{\text{af}}-q_{\text{f,sat}}\left)\right)}{\left(\right(r_{\mathrm{a}}+r_{\mathrm{s}}\left)\right)}}\cdot (\mathrm{1}-r).\end{array}$$

Here r_a is the aerodynamic resistance (Louis, 1979) and r_s is the stomatal resistance.

3 Evaluation of the modified WRF multi-layer green roof model

3.1 Surface energy balance and temperature

Performance of the original and modified WRF multi-layer green-roof models was assessed against measurements of latent (Q_E) and ground (Q_G) heat fluxes at the surface, surface temperature (T_s), and soil water content (θ), along with a residual estimate of sensible heat flux (Q_H).

Time series show small warm biases for both modelled surface temperatures (Fig. 2a and b). For latent heat flux (Q_E), the original and modified models show similar overall performance in summer, with identical RMSE and R² values (61.8 W m⁻² and R²=0.40, respectively). However, the original model has a slightly lower MAE and higher refined index of agreement (d) than the modified model (40.4 vs. 43.9 W m⁻²; d=0.67 vs. 0.64). In fall, the modified model shows a lower RMSE and higher R² than the original model (59.7 vs. 64.7 W m⁻²; R²=0.29 vs. 0.17), while the original model has a slightly lower MAE and slightly higher d (39.6 vs. 41.0 W m⁻²; d=0.60 vs. 0.59) (Table 1; Fig. 2e and f). For ground heat flux (Q_G), the largest model observation discrepancy occurs at night, when the observed Q_G signal becomes nearly flat and approaches zero, whereas both model versions simulate a stronger nighttime heat release. The experimental configuration can influence measured Q_G throughout the full diurnal cycle because the elevated modules leave an air space beneath the green roof system and the individual modules introduce exposed lateral boundaries. These effects appear to be strongest at night due to smaller energy balance forcing (i.e., Q^∗) at this time, but may also influence daytime measurements. We therefore additionally report daytime-only Q_G metrics (Q_G_daytime) in Table 1. Nevertheless, the air gap may also impact the daytime measurements (further discussed later in this section based on measurements in 2025).

Overall, for Q_G the modified model delivers a clear, season-independent improvement, with significant reductions in total error. The removal of bias evident in the original model simulations (lower RMSE and higher index of agreement; Table 1) is visible in both seasons (Fig. 2c and d). This indicates a substantive correction to substrate heat storage, which impacts modelled Q_E and critically, Q_H.

Figure 2 Comparison of the original and modified multi-layer green roof models against observations during 1 July to 31 August 2014 (Summer; left column) and 1 September to 31 October 2014 (Fall; right column). In each panel, the original model is shown with blue circles and the modified model with orange triangles. Panels (a) and (b) show surface temperature, T_s (°C); panels (c) and (d) show soil heat flux, Q_G (W m⁻²); and panels (e) and (f) show latent heat flux, Q_E (W m⁻²). The dashed red line indicates the 1:1 line, and R², RMSE, and MAE are reported in each panel.

In terms of soil water content (θ), the original model shows slightly better statistical agreement with the observations during summer, with a lower RMSE and higher refined index of agreement (d) than the modified model (Table 1). In fall, however, the modified model shows slightly better agreement, with a lower RMSE and a less negative d. Overall, the two model configurations produce similar seasonal evolution of soil water content, particularly during summer (Figs. 3d and 4d). The relatively good summer performance of the original model does not necessarily indicate that its simpler hydrological representation is more physically realistic, because compensating differences in atmospheric water loss and drainage reduce overall errors. This behaviour is examined further in Sect. 3.2 based on the seasonal water-balance analysis presented in Sect. S3 and Table S1 in the Supplement.

Figure 3 Observed and modified model daily (a) surface temperature T_s, (b) latent heat flux Q_E, (c) ground heat flux Q_G, and (d) volumetric soil water content θ for a period from 1 July to 31 August 2014. In panel (b), gray lines indicate observed latent heat flux during rainy days.

Figure 4 Same as Fig. 3, except for the original model.

Both models slightly overestimate daytime surface temperature maxima; at night, the modified model shows a stronger cold bias, while the original exhibits a smaller bias (Figs. 3a and 4a). The original model underestimates daytime latent heat flux and largely misses the immediate post-rainfall increase, whereas the modified model better matches both the timing and magnitude of the peaks. For direct visual comparison of the timing and magnitude of the original and modified model responses, the corresponding time series are overlaid in Fig. S4. In particular, the sharp increase in Q_E after rainfall is clearly seen in the observations but is missed by the original model. By adding evaporation from intercepted water on the foliage, the modified model can now reproduce this post-rainfall boost in latent heat flux. This addition allows the model to respond quickly to wet conditions and better match the observed evapotranspiration peaks (Figs. 3b and 4b). To illustrate the post-rainfall latent-heat-flux response more clearly, Fig. 5 shows an event spanning 15–17 July 2014 (DOY 196–198). Following precipitation on DOY 196, the leaf wetness sensor signal increases sharply, indicating rapid wetting of the sedum canopy surface. At the same time, the modified model produces a strong, short-lived increase in modeled latent heat flux (Q_E), consistent with enhanced evaporation immediately after rainfall (including evaporation of intercepted water on the canopy). As the leaf surfaces dry, both the leaf wetness signal and modelled Q_E decrease. The lysimeter-derived observed Q_E also increases rapidly following rainfall and shows close temporal agreement with the modeled response, with only a small difference in the timing of the peak. In fact, observed Q_E is derived from changes in lysimeter mass over each measurement interval, which can contribute to small differences in the apparent timing of the responses. In contrast, the original model shows a much more muted response and largely fails to capture the rapid post-rainfall peak.

Figure 5 Event-scale comparison of modeled latent heat flux (left axis, W m⁻²) with leaf wetness sensor signal at the sedum surface (right axis, mV) and precipitation (red bars, mm) for DOY 196–198 (15–17 July 2014). The leaf wetness sensor reports the raw voltage output (mV) from the probe; higher values indicate a wetter surface. Modeled Q_E is shown in dashed orange for (a) modified model and (b) original model. Observed Q_E is shown as the solid black line with markers.

The most pronounced differences between the model versions appear in the ground heat flux: the modified model damps the exaggerated daytime spikes seen in the original model and exhibits a more realistic nocturnal recovery, consistent with improved substrate thermal conductivity and heat storage – and hence better phase and amplitude (Figs. 3c and 4c). In addition to the long-period evaluation, Fig. 6 shows an event-scale comparison of observed and modeled ground heat flux (Q_G) for 11–18 October 2025 (DOY 284–291), a week with clear observations. During this period, the experimental green-roof module was placed directly on the roof deck (i.e., without an airspace beneath), providing lower-boundary conditions for vertical heat-flux measurements that are both more typical of actual green roofs and which match the model assumptions. Under these conditions, the observed Q_G exhibits a physically consistent diurnal cycle, and the modified model reproduces the main phase and amplitude of Q_G variability, including the daytime peaks and nighttime recovery (R²=0.642).

Figure 6 Observed and modeled (modified model) ground heat flux (left axis) and precipitation (right axis) for 11–18 October 2025 (DOY 284–291).

Table 2 Summary of one-at-a-time sensitivity experiments in the WRF multi-layer green roof scheme. Short names correspond to titles in Sect. 2 and labels used in Figs. 8 and 9.

The original model exaggerates the diurnal cycle of Q_G, diverting energy from the turbulent fluxes and depressing morning values as well as midday peaks in Q_E and Q_H (Fig. 7a). Figure 7b demonstrates that the modified model captures both the magnitude and timing of peak latent heat flux Q_E, yielding a realistic daytime partitioning of available energy. The modified model tempers the large positive morning Q_G bias in the original model and as a result better reproduces the timing and amplitude of the daytime peaks of Q_H and Q_E. Differences in Q^∗ are small between the models.

Figure 7 Diurnal mean energy-flux composites for the original model (a) and the modified model (b). Solid lines are model output; dotted lines are observations. Curves are hourly means from 1 July to 31 August 2014.

Figures 8 and 9 present a one-at-a-time sensitivity analysis, where each of the five modifications summarized in Table 2 and described in detail in Sect. 2.3 is switched off individually while the other updates remain active. These experiments are used to isolate the role of each added process and to show how the different modifications shape the modelled surface energy partitioning relative to the observations.

Figure 8 Seasonal distributions of observed and modeled daytime (06:48–20:02 LST) latent heat flux Q_E for 2014, shown separately for summer (a, c, e; left column) and fall (b, d, f; right column). Panels (a) and (b) multi-layer plant water-uptake removed; reverts to the old uptake treatment (–root uptake), panels (c) and (d) canopy interception storage and its evaporation removed (–intercepted water), and panels (e) and (f) explicit soil-surface evaporation removed (–soil evap). Each panel compares observations with the full modified configuration (all on), the corresponding switch-off experiment, and the original model (all off).

Figure 9 Seasonal distributions of daytime (06:48–20:02 LST) observed and modeled soil heat flux Q_G for 2014, comparing observations with the modified configuration (all on), –ksoil (new soil thermal-conductivity formulation removed), –heat storage (vegetation-modulated effective surface thermal conductivity term removed), and the original model (all off).

In Fig. 8, the six panels focus on moisture related processes and how they affect the total latent heat flux (Q_E) in both summer and fall (focus on daytime-only). When multi-layer root uptake is removed (–root uptake; Fig. 8a and b), plants rely primarily on moisture from the uppermost soil layer. As this layer dries during dry-down periods, transpiration becomes increasingly water-limited, leading to a modest reduction in mean daytime Q_E, with a small decrease in the median, particularly in summer, bringing median modeled values below both the observations and the fully modified configuration, and a narrowing of the upper distribution. Removing canopy interception storage and its subsequent evaporation (–intercepted water; Fig. 8c and d) suppresses short-lived evaporation pulses following rainfall events, lowering the summer maximum and the fall upper quartile of the daytime Q_E distribution and shifting the median to lower values in both summer and fall relative to the observations and the fully modified model. Eliminating explicit soil-surface evaporation (–soil evap; Fig. 8e and f) further decreases daytime Q_E by removing direct evaporation from the exposed soil. In the present simulations, vegetation cover is high (86 %), limiting the exposed soil fraction and therefore the magnitude of this effect. Under lower vegetation coverage, however, soil-surface evaporation would contribute more substantially to daytime Q_E, and the sensitivity to this process would be expected to increase.

Figure 9 highlights the thermal-related modifications, which produce the clearest and most systematic sensitivity in the soil heat flux Q_G. Turning off the revised soil thermal-conductivity formulation (–ksoil) causes the Q_G distribution to widen substantially and more extreme values occur, showing that the updated substrate thermal treatment is one of the key controls on soil heat exchange. Removing the vegetation-modulated effective surface thermal conductivity term (–heat storage) also increases the spread and magnitude of Q_G. This modification primarily affects the mean behavior of Q_G by regulating how efficiently incoming radiation is conducted into the substrate. Without vegetation modulation, the roof surface behaves more like bare soil, absorbing and storing a larger fraction of solar energy, which systematically enhances heat storage and leads to stronger daytime Q_G peaks. Together, these results indicate that while the revised substrate thermal conductivity governs the range and extremes of Q_G, the vegetation-dependent surface thermal conductivity plays a dominant role in controlling the mean partitioning of energy into heat storage, thereby reducing persistent biases in Q_G.

Figure 10 Seasonal totals of observed rainfall (line) and cumulative drainage (bars) for summer (DOY 182–243) and fall (DOY 244–304) 2014. Drainage totals are shown for observations and for the original and modified model configurations.

3.2 Water balance and hydrological processes

Figure 10 compares the seasonal cumulative drainage over the full observation periods in summer (DOY 182–243) and fall (DOY 244–304) 2014, together with the total observed rainfall during each season. In summer, both model configurations produce lower drainage than observed (45.90 mm), with the original model giving 43.31 mm and the modified model giving 32.54 mm. In fall, both models produce more drainage than observed (104.90 mm), reaching 131.21 mm in the original model and 125.10 mm in the modified model. Overall, the seasonal totals indicate an underestimation of drainage in summer and an overestimation in fall, while the modified configuration consistently reduces drainage relative to the original model. As an additional conservation diagnostic, the seasonal water-balance closure errors were 1.95 % and 1.31 % for the original model in summer and fall, respectively, compared with only 0.13 % and 0.16 % for the modified model, after accounting for precipitation, net atmospheric water loss, drainage, and changes in water storage (Sect. S3 and Table S1). These results also help explain why the original model can reproduce soil water content reasonably well despite its simpler hydrological representation. In summer, the original model produces lower atmospheric water loss but higher drainage than the modified model, so these differences partly compensate and result in a similar seasonal change in soil water storage. Therefore, the slightly better summer soil water content statistics of the original model do not necessarily indicate a more physically realistic representation of the water balance.

4 Discussion

4.1 Model strengths and residual weaknesses after targeted upgrades

The performance of the modified WRF multi-layer green roof scheme can be evaluated relative to other physically based green roof models. Simulations of extensive green roofs, such as those presented by Heusinger et al. (2018) with a modified version of the EnergyPlus green roof model “EcoRoof”, explicitly represent vegetation processes including canopy transpiration, soil evaporation, and evaporation of intercepted water, supported by a detailed canopy microclimate sub-model that resolves radiative and turbulent exchanges between the atmosphere, the canopy, and the surface. Simulations with this modified EcoRoof (ER-modified) model for London indicate that the enhanced canopy microclimate formulation improves surface temperature predictions relative to the original WRF multi-layer green-roof scheme, but with a tendency to overestimate latent and ground heat fluxes (Martinez Mendoza et al., 2026). Because ER-modified and WRF-MLGR v2.0 are evaluated against the same observational dataset, a direct quantitative comparison can be made. ER-modified demonstrates stronger agreement for surface temperature and ground heat flux, with a index of agreement (d) of 0.77 for T_s and 0.74 for Q_G, when evaluated over the full July–October 2014 period (i.e., combined summer and fall). In comparison, WRF-MLGR v2.0 yields seasonal values of 0.56 (summer) and 0.62 (fall) for T_s and 0.35 (summer) and 0.37 (fall) for Q_G in WRF-MLGR v2.0. In contrast, WRF-MLGR v2.0 exhibits improved skill in simulating latent heat flux, with d values of 0.64 (summer) and 0.59 (fall), exceeding the long-period performance of ER-modified (d=0.55).

It is relevant to additionally compare other green-roof models. A separate WRF green-roof scheme has also been developed for the single-layer urban canopy model (Yang et al., 2015). This scheme is different from the multi-layer WRF-BEP+BEM green-roof scheme of Zonato et al. (2021) evaluated and modified in the present study. The single-layer scheme does not document a dedicated interception storage and evaporation process; however, it was shown to reproduce surface temperature and turbulent energy fluxes within a reasonable range when evaluated against observations over a 6 d period. Another approach is that the TEB-GREENROOF model directly couples the vegetated substrate to the roof layers without a canopy air space. When evaluated, it showed underestimated substrate moisture, overestimated drainage, and larger than observed temperature amplitudes, reflecting a warm bias and discrepancies in moisture dynamics (de Munck et al., 2013). In the following, we assess the behavior of the modified WRF multi-layer green roof scheme and its strengths, weaknesses, and the impacts of the applied updates.

This study develops a modified version of the WRF multi-layer green-roof scheme and evaluates it against observations in summer and fall of 2014 in London, Canada. The goal of the updates is to address specific process-based limitations identified in the original scheme. In brief, the modified configuration adds (i) a revised non-linear moisture-dependent substrate thermal conductivity formulation (replacing the original linear relationship), (ii) a vegetation-modulated effective surface thermal conductivity term, (iii) explicit soil-surface evaporation, (iv) multi-layer root water uptake, and (v) a canopy interception module that represents interception storage, evaporation of intercepted water, and dew formation on the canopy. Below, we summarize which aspects improve model results most clearly, where key biases remain, and what these results suggest for future development and applications within WRF.

Across both seasons, the largest improvement is in ground heat flux (Q_G). In summer, Q_G mean absolute error decreases from 87.4 to 17.6 W m⁻² and RMSE decreases from 105.9 to 24.0 W m⁻² (Table 1). In fall, Q_G RMSE decreases from 94.2 to 24.0 W m⁻² (Table 1). The time series and diurnal plots explain this improvement: the original scheme produces Q_G that swings too strongly, with peaks that are too large during the day and also exaggerated at night, especially an unrealistically sharp morning spike. The modified scheme reduces these exaggerated peaks, so the daily cycle has a more realistic timing and magnitude (Figs. 2–4 and 7).

It is also important to interpret Q_G in the context of the measurement setup. During the 2014 deployment the design of the experimental green roof array provides an elevated platform to allow installation of lysimeters and drainage units below the modules, which also results in an air layer beneath the green roof array. This configuration allows heat to exchange through the bottom (and potentially sides) of the module, which can make the plate-based Q_G measurement harder to interpret – particularly at night, when signals are smaller and boundary effects can be more influential. For this reason, we emphasize a separate one-week period (11–18 October 2025), when the modules were placed directly on the roof deck and the physical setup more closely matched the model's lower-boundary assumption of direct contact between the green roof soil and the roof deck underneath. During this week, the observed Q_G shows a clearer diurnal cycle, and the modified model aligns closely with the measurements in both timing and magnitude (Fig. 6), providing strong evidence that the revised physics substantially improves the simulation of ground heat flux (Q_G) under conditions that best match the model configuration (and the configuration of most actual green roofs).

Reducing the unrealistic daytime Q_G peak changes how the available energy is divided between the different heat fluxes. Net radiation (Q^∗) is very similar in the original and modified simulations. Therefore, when Q_G becomes smaller in the modified model, more energy is available for Q_E and Q_H.

The main reason for the large change in Q_G is the improved soil thermal conductivity. Ground heat flux depends on both the temperature difference within the soil and the soil thermal conductivity. Therefore, Q_G can change a lot even if the modeled surface temperature does not change very much. In the original model, the soil thermal conductivity became unrealistically high in the morning, which caused too much heat to move into the soil. In the modified model, the revised soil thermal conductivity and especially the vegetation-dependent surface conductivity reduce this excessive heat transfer. As a result, more of the available energy is directed to Q_E and Q_H, particularly earlier in the day, shifting the timing of the diurnal cycle of these fluxes.

The hydrological modifications further influence this coupled energy balance. Soil evaporation provides an additional pathway for water loss from the green roof, canopy interception allows stored rainfall to evaporate rapidly after rain, and multilayer root uptake allows plants to access water from deeper soil layers as the upper soil becomes dry. These processes generally enhance Q_E when water is available. The resulting changes in Q_E, Q_G, surface temperature, and moisture conditions also affect Q_H. Therefore, the modified model produces a more realistic overall partitioning and diurnal evolution of the surface energy fluxes, even though changes in surface temperature remain relatively small.

The hydrological modifications also help explain the changes in latent heat flux (Q_E). Figures 3b and 4b show the diurnal cycles with a noticeable increase in Q_E right after rainfall. The original model often represents this response too weakly or misses it, which contributes to an overall underestimation of Q_E during wet periods. In contrast, the modified model captures these short-lived post-rain increases more consistently and provides a better overall prediction of Q_E variability and magnitude. This behavior is especially clear during 15–17 July 2014 (DOY 196–198). Following rainfall, the leaf wetness sensor signal rises sharply, indicating rapid wetting of the sedum canopy surface, and the modified model simultaneously produces a strong, short-lived spike in modeled Q_E. As the leaf surfaces dry, both the sensor signal and modeled Q_E decrease. The original model, however, produces a much weaker response and largely fails to capture this rapid post-rainfall peak (Fig. 5). Overall, this suggests that the modified scheme better represents fast evaporation pathways immediately after wetting, consistent with the inclusion of evaporation from intercepted water on the canopy. The most direct explanation is the added canopy interception module, which allows intercepted water stored on the canopy to evaporate quickly after rainfall (and also represents dew formation), rather than relying only on slower changes in substrate moisture and evaporation. This is consistent with the broader land-surface modeling literature, where interception loss is recognized as an efficient, rapid pathway for returning precipitation to the atmosphere and is often needed to reproduce these short post-rain evaporation responses (e.g., Cammalleri et al., 2010; Heusinger et al., 2018).

Another important improvement comes from the updated root water uptake. In the original scheme, transpiration draws water mainly from the thin top layer, so Q_E drops when that layer dries and more water ends up draining from the bottom of the soil layer. In the modified model, uptake can access moisture from multiple soil layers, which supports more realistic Q_E during dry-down and reduces excess drainage. We also added explicit evaporation from the soil surface. Together, these changes improve the overall energy partitioning. The drainage comparison is consistent with this behavior: the modified model produces less drainage than the original model (Fig. 10). However, drainage remains underestimated in summer and overestimated in fall in both model configurations.

Even with these improvements, some biases remain. First, daytime surface temperature (T_s) is still too warm in both model configurations (Figs. 3a and 4a). A key reason for this remaining T_s bias is that the current WRF-MLGR formulation calculates only a single bulk surface skin temperature for the green-roof system. Separate computation of the soil/substrate surface temperature, vegetation/canopy surface temperature, and canopy-air temperature is not available in the present formulation. As a result, T_s cannot be directly interpreted as either a soil-surface temperature or a vegetation-canopy temperature, but instead represents an aggregated surface temperature produced by the simplified bulk representation of the green roof. This limits the ability to diagnose whether the remaining T_s bias arises primarily from the substrate thermal state, vegetation thermal response, or both. Second, Q_E is still low and shows less variability than the observations (Fig. 8). The sensitivity tests help explain this: when canopy interception storage and its evaporation (–intercepted water) or multi-layer root uptake is removed (–root uptake), Q_E decreases further and more energy can shift into Q_H and Q_G (Fig. 8). This shows that the added moisture pathways are working as intended. However, the remaining underestimation of Q_E variability suggests that the model still misses some real-world controls, such as microclimatic gradients within the canopy air space and plant physiological responses (for example, a more detailed canopy microclimate). For example, the EnergyPlus green roof (EcoRoof) model represents vegetation as an explicit canopy, including within-canopy radiative exchange (e.g., shading effects) and canopy air exchange processes that regulate evapotranspiration, which may help capture some of the canopy-microclimate controls that are still simplified in the present WRF configuration (Sailor, 2008). Future development should therefore separate the vegetation canopy, canopy air, and soil/substrate surface in the model, instead of representing them with one bulk surface temperature. This would allow the model to calculate their heat exchange processes separately, including radiation, turbulent heat exchange, and heat conduction into the substrate. It would also make it possible to identify whether any remaining errors in T_s come mainly from the vegetation or the soil surface, and to improve the inaccurate model component.

In terms of substrate moisture, the original model better captured the observed temporal variability in volumetric water content (θ) during summer, whereas the modified model showed improved agreement in fall (Table 1). However, both model versions show a small bias in drainage, modestly underestimating summer totals and overestimating fall totals (Fig. 10). During the summer of 2014, conditions were relatively dry, with even the strongest convective rainfall events producing less than 14 mm. In contrast, fall 2014 saw powerful thunderstorms, with notable events exceeding 73 mm of rainfall in just 3 h and another event surpassing 45 mm in 2 h. This difference in seasonal rainfall intensity helps explain the models biases. During drier summer periods, the models may underestimate drainage since soils rarely reach near-saturation. Conversely, during fall's intense rain events, soil often approach saturation, prompting the models to produce excessive drainage. The modified model produces marginally lower drainage than the original, consistent with the added water-loss pathways through soil evaporation, transpiration via root uptake across all soil layers, and evaporation of intercepted water from the canopy.

4.2 Limitations and future work

Several limitations affect how these results should be interpreted and also highlight priorities for future work. The main observational limitation is ground heat flux (Q_G), which needs to be considered in the context of the experimental setup. During the 2014 deployment, the elevated green roof modules left an air space between the modules and the underlying roof, while the use of individual modules introduced exposed lateral boundaries. These conditions allow additional heat exchange through the bottom and side boundaries that is not represented in the model. The influence appears to be strongest at night, when observed Q_G becomes nearly flat and approaches zero, suggesting weak thermal coupling between the green roof modules and the underlying roof structure. However, the same boundary effects may also influence daytime Q_G. We therefore report daytime-only Q_G statistics (Q_G,day; Table 1) as an additional evaluation metric, while recognizing that daytime measurements are not completely unaffected by the experimental configuration. The October 2025 comparison, conducted under a configuration more representative of the modeled lower boundary, shows a clearer observed diurnal Q_G cycle, and the modified model reproduces the main timing and magnitude of the observed variability more closely (Fig. 6).

A related limitation is that the “observed” sensible heat flux (Q_H) is not measured directly. Instead, it is calculated as a residual: net radiation minus the sum of measured Q_E and measured Q_G. This means that any bias in Q_G, during the day and especially at night – will automatically appear as a bias in the residual Q_H, which makes nighttime Q_H comparisons harder to interpret. Importantly, this does not change one of the main conclusion of this study: the modified model dramatically reduces the large Q_G errors seen in the original scheme. Rather, the residual nature of Q_H mainly limits how confidently we can interpret the remaining nighttime Q_H differences as resulting from model physics.

From a model-development perspective, the upgrades in this study are focused on specific processes, but several simplifications remain. Even with the new interception module, the canopy representation is still basic. For example, the scheme does not fully represent a separate canopy microclimate (its own temperature and humidity near the leaves) or detailed radiation exchange within the canopy. It also simplifies how vegetation modifies near-surface aerodynamics, such as wind reduction, roughness, and the resistances that control heat and moisture transfer. Probably most important is separation of the soil surface and canopy radiation and energy balances together with addition of a clearer resistance/aerodynamic formulation. A likely contributor to the surface temperature bias is that the observed canopy surface (as seen by the infrared radiation thermometer) has substantial three-dimensional structure, with a mixture of sunlit and shaded elements during the day, and locally reduced sky-view factors at night that keep sheltered portions of the canopy warmer than exposed upper surfaces. However, the model represents the vegetation canopy as a flat, uniform layer, thus generating an apparent overheat during the day and too much nocturnal cooling relative to the measured values. These additions would likely help reduce the remaining surface-temperature biases (too warm during the day and too cool at night), improve the spread and magnitude of Q_E (Fig. 8), and further improve the daytime heat storage and nighttime release. Incorporating these processes consistent with field studies showing that green-roof energy partitioning depends strongly on soil moisture, vegetation state, and canopy wetness (e.g., Heusinger and Weber, 2017).

Also, this evaluation focuses on warm-season conditions (summer and fall 2014). For Canadian and broader climate applications, cold-season processes may also be important. Future development should include snow accumulation and melt, freezing and thawing of the substrate, and temperature-dependent hydraulic and thermal properties, because these processes can strongly modify wintertime heat exchange and spring moisture availability. In addition, the scheme could be extended to represent green roof uptake of trace gases – such as CO₂ exchange through photosynthesis (e.g., Mirebeau et al., 2025) and a simple deposition/uptake pathway, so the model can also be used to explore potential green-roof impacts on urban air quality.

Finally, since the ground-level vegetation scheme in WRF-BEP+BEM is based on this green roof model, these upgrades can be applied in the future to further improve the accuracy of urban vegetation simulations within that framework.

5 Conclusion

This study introduces targeted process-based modifications to the WRF multi-layer green-roof scheme and evaluates the updated model using field observations from 2014 for continuous summer and fall periods. To better match the model's lower-boundary setup and to interpret ground heat flux (Q_G) more clearly, we also analyze a one-week period (11–18 October 2025) when the green roof modules were placed directly on the roof deck.

The targeted process changes to the green roof model address key limitations in the original scheme. These updates include (i) a revised non-linear, moisture-dependent treatment of substrate thermal conductivity, (ii) a vegetation cover-modulated effective surface thermal conductivity, (iii) explicit soil-surface evaporation, (iv) multi-layer root water uptake, and (v) a canopy interception module that accounts for interception storage, evaporation of intercepted water, and dew formation on the canopy.

Across both seasons, the largest improvement is in ground heat flux (Q_G). In the original scheme, the diurnal Q_G cycle was too strong, with unrealistically large peaks during morning and midday and overly large variations at night. The modified scheme reduces these overestimations and better reproduces both the timing and the magnitude of the daily cycle, leading to much lower Q_G errors (Table 1). This improvement is also supported by the October 2025 test week: when the modules were placed directly on the roof deck, the observed Q_G shows a clearer diurnal pattern, and the modified model matches it closely.

Improving Q_G also ameliorates partitioning of available energy among the remaining surface fluxes. When Q_G is too large, the model directs too much energy into the ground and leaves too little for turbulent exchange, which overly reduces the morning and midday latent heat flux (Q_E) and sensible heat flux (Q_H). By mitigating this unrealistic Q_G draw, the modified scheme reproduces the timing and magnitude of the daytime Q_E and Q_H peaks more accurately. The modified model also represents Q_E more realistically: it captures brief increases in Q_E after rain more consistently, supported by evaporation of intercepted water (including dew-related canopy exchanges), added soil-surface evaporation, and multi-layer root uptake that allows transpiration to draw water from deeper layers. Together, these changes reduce drainage compared with the original model.

Despite this progress, some biases remain. Daytime surface temperature (T_s) is still too warm, nighttime T_s can be too cool, and daytime Q_E is still somewhat low and shows less variability than the observations.

Overall, these results show that the targeted physical updates greatly improve how the WRF MLGR scheme simulates energy and water exchanges, mainly by correcting the substrate thermal treatment as well as adding or modifying select hydrological processes. The modified scheme therefore provides a stronger basis for WRF-based urban climate applications.