On the resolution requirements for accurately representing interactions between plant canopy structure and function in three-dimensional leaf-resolving models

While functional–structural plant models (FSPMs) have been proposed as a tool for better analysing and predicting interactions between plant structure and function, it is still unclear as to what spatial resolution is required to adequately resolve such interactions. Shadows cast by neighbouring leaves in a plant canopy create extremely large spatial gradients in absorbed radiation at the sub-leaf scale, which are usually not fully resolved in ‘leaf-resolving’ plant models. This failure to resolve sharp radiative gradients can propagate to other depend-ent biophysical models, and result in dramatic overprediction of whole-plant and -canopy fluxes with errors significantly higher than that of a statistical ‘big leaf’ or turbid medium model. Under-resolving radiative gradients creates a diffusive effect in the probability distribution of absorbed radiation, and smears out the effect of canopy structure, effectively undermining the original goal of a leaf-resolving model. Errors in whole-canopy fluxes of photosynthesis increased approximately linearly with increasing LAI , projected area fraction G , and decreased logarithmically as the fraction of incoming diffuse radiation was increased. When only one discrete element per leaf was used, errors in whole-canopy net CO 2 flux could be in excess of 100 %. Errors due to sub-leaf resolution decreased exponentially as the number of elements per leaf was increased. These results prompt closer consideration of the impact of sub-leaf resolution on model errors, which is likely to prompt an increase in resolution relative to current common practice.


IN TROD UCTION
The complex spatial structure of plant systems often creates high heterogeneity in biophysical processes across a wide range of spatial scales at any instant in time. The interception of radiative energy from the sun, which is the primary driver of these biophysical processes, exhibits particularly high spatial heterogeneity. Unlike most continuum-based processes, incident radiation propagates in approximately straight lines with minimal diffusion over length scales characteristic of plant systems. This can create extremely sharp spatial gradients in absorbed radiation, moving from full sun to full shade over the scale of single organs. When attempting to model nearly any physical process, sharp gradients in space or time create difficulties in schemes for numerical representation (Press et al. 2007). This is primarily because directly representing these sharp gradients usually requires an extremely fine discretization in time or space, and failure to do so typically results in an artificial numerical diffusion effect that reduces the gradient magnitude and can have dramatic effects on the resulting model performance. A classical example is in numerical simulation of turbulent fluid flow at high Reynolds number, where the range of relevant turbulent length scales can be as wide as six orders of magnitude or more depending on the Reynolds number (Pope 2000). If the numerical grid is not fine enough to fully represent this entire range of scales (and associated gradients), dissipation of turbulent energy is not correctly represented, potentially resulting in unphysical behaviour of the model (Komen et al. 2017).
An alternative and less computationally demanding approach to modelling processes with small characteristic length scales is to represent their effects statistically. For interception of radiative energy by a plant canopy, a statistical approach assuming the canopy consists of a quasi-homogeneous 'turbid' medium is commonly used. Rather than directly representing every sunfleck and shadow in the canopy, most models assume that interception of radiation by a medium of vegetation is a stochastic process following an exponential probability distribution, which can be readily derived from first principles (e.g. Bailey et al. 2020). Once the probability distribution describing radiation interception is known, it is relatively straightforward to integrate the distribution across the canopy to yield the cumulative absorption, and fractions of sunlit and shaded surface area (DePury and Farquhar 1997;Wang and Leuning 1998). However, the challenge with this approach is that one must be able to assume that the interception probability distribution is homogeneous over some spatial scale much larger than that of individual leaves. Thus, for canopies with high vegetative heterogeneity, assuming a constant probability distribution across the entire canopy can yield considerable errors (Ponce de León and Bailey 2019; Bailey et al. 2020).
Functional-structural plant models (FSPMs) have been more recently proposed as a general class of models that attempt to better represent interactions between plant structure and function by resolving heterogeneity at the leaf scale (Vos et al. 2010;Prusinkiewicz and Runions 2012;Evers et al. 2018). It is often presupposed that such a detailed spatial representation is necessary or at least advantageous in accurately representing plant-environmental interactions, particularly aspects related to interactions between vegetative structure and radiation transport. Most FSPMs include detailed ray-tracing-based models for radiation transfer that directly represent the fate of photon 'beams' as they propagate through a geometrically resolved, three-dimensional (3D) canopy. Thus, they can more directly include the effects of scattering, transmission, leaf angle, leaf area, etc. than simpler, statistically based models that aggregate groups of vegetative elements into homogeneous volumes (e.g. Sinclair et al. 1976;Wang and Jarvis 1990;Cescatti 1997;Monsi and Saeki 2005;Bailey et al. 2014). This detailed representation of plant geometry and radiation interception can be coupled with other physiological models, which enable direct investigation of the impacts of leaf topology, phyllotaxy, and within-plant diffusive and convective transfer, among others.
Recent work by Kent and Bailey (2021) has suggested that, while leaf-resolving plant models may be able to accurately predict integrated radiation attenuation with relative ease, resulting leaf-averaged fluxes provided as inputs to coupled non-linear leaf-level biophysical models can create large errors in whole-canopy fluxes. Because of the non-linear response of many biophysical processes to absorbed radiation flux, Kent and Bailey (2021) demonstrated that failure to resolve sub-leafscale heterogeneity in radiation interception in a leaf-resolving model could create errors in canopy-level fluxes that were on average larger than predictions using a multilayer big-leaf (Norman 1993;Pyles et al. 2000;Baldocchi and Wilson 2001;Clark et al. 2011) or two-big-leaf model (DePury andFarquhar 1997;Wang and Leuning 1998).
The prior results of Kent and Bailey (2021) reveal a fundamental, yet under studied knowledge gap in 3D plant modelling, which is that the resolution necessary to adequately resolve biophysical processes is not fully understood, nor is its full impact on the representation of plant structure and function interactions. The goal of this work was to more thoroughly investigate and understand how under-resolving radiative gradients in 3D plant simulations impacts the ability to represent interactions between plant structure and function. It is not the intention of this work to precisely determine errors corresponding to any previous or future 3D plant model simulation, as clearly such errors will be case-specific. Rather, the aim is to increase awareness of the existence of such errors and provide a basic overview of their behaviour for select example cases.

Simulation experiment
Radiation interception and photosynthetic production was simulated using Helios software version 1.2.8 (Bailey 2019) for a range of homogeneous canopies to explore the effects of sub-leaf resolution on wholecanopy fluxes. A homogeneous canopy geometry was generated in the model by randomly sampling leaf positions from a uniform distribution within a horizontal slab (Fig. 1). Leaf angles were generated by randomly sampling from one of four of the canonical leaf angle distributions introduced by Wit (1965): spherical, uniform, planophile and erectophile. The number of leaves in the canopy was chosen to achieve one of four LAI values: 0.5, 1.0, 2.0 and 3.0. The horizontal extent of the canopy represented in the simulations was finite and equal to 10 × 10 m 2 , but a periodic boundary condition was applied in the horizontal to effectively achieve an infinite canopy. The canopy height was set to 1.0 m, and leaves were squares with one-sided area of 0.05 × 0.05 m 2 . In order to test the impact of sub-leaf resolution, leaves were represented geometrically using a uniform grid of squares with resolutions of N sub equal to 1, 3 2 , 5 2 , 10 2 and 15 2 elements per leaf.
Helios simulates the absorbed radiation flux for every leaf element (or sub-element) in the model domain using a reverse raytracing approach that is described and evaluated in detail in Bailey (2018). A collimated solar radiation source was placed above the canopy. For simplicity and to consolidate analysis, only a single sun direction was considered, which is a solar zenith of θ s = 0. Leaves and the ground were also considered to be black for simplicity. Separate simulations were performed with varying diffuse radiation fraction: f diff of 0, 0.1 and 0.2. Diffuse radiation was isotropic, which means that conditions were much more diffuse than on a typical clear-sky day in which most diffuse radiation is concentrated near the direction of the sun (Harrison and Coombes 1988). The collimated radiation source had a flux in the photosynthetically active band of (1200 μmol m −2 s −1 )·(1 − f diff ) normal to the sun direction, while the diffuse radiation flux on a horizontal surface above the canopy was (1200 μmol m −2 s −1 )·(f diff ). The simulations used 150 direct and 500 diffuse rays per element.
The calculated radiative flux value for each leaf sub-element in the canopy was used as input for a model of photosynthesis. Helios includes an implementation of the biochemical photosynthesis model of Farquhar et al. (1980), which also generally requires coupling with a model for stomatal conductance and leaf temperature. In order to avoid convoluting effects of these other models, leaf-level photosynthesis in this study was represented using a simple rectangular hyperbolic function of light where A is net CO 2 assimilation flux on a per unit leaf area basis (μmol CO 2 m −2 s −1 ), Q is the absorbed photosynthetic photon flux density (μmol photons m −2 s −1 ), A max is the asymptotic maximum in the assimilation flux at infinite Q, α is a shape parameter describing the response of A to varying Q and R d is the dark respiration rate. Although Eq. (1) is an over simplified model of photosynthetic response to light, it provides a functionally similar response as the model of Farquhar et al. (1980) at constant temperature, ambient CO 2 concentration and photosynthetic capacity. Assumed model parameters are given in Table 1, and a plot showing the resulting light response of A is shown in Fig. 2. A visualization of the resulting distribution of A across a partially shaded leaf is shown in Fig. 3. A factorial simulation experiment design was used to run all possible combinations of LAI, leaf angle distribution, f diff and N sub cases, which resulted in 240 total simulations. For the purpose of presenting results, a default set of parameters was chosen in some instances, which are tabulated in Table 1.

Analysis of model outputs
Model outputs consisted of fluxes of absorbed photosynthetically active radiation Q and net CO 2 flux A for every leaf or sub-leaf element for each canopy case. These fluxes were then used to calculate discrete probability distribution functions (PDFs) across the entire canopy using 25 bins.
Integrated whole-canopy fluxes of absorbed radiation Q or net CO 2 A were calculated as where Q c and A c are, respectively, the photosynthetic photon flux density and net CO 2 on a per unit ground area basis, N is the number of leaf elements in the simulated canopy, a ,i is the one-sided surface area of the ith leaf element, Q i is the absorbed photosynthetic photon flux density of the ith leaf element, A i is the net CO 2 flux of the ith leaf element and a g is the ground surface area occupied by the simulated canopy. The error in integrated whole-canopy fluxes Q c or A c due to leaf element resolution was represented using the percentage difference between the simulated Q c or A c value and the corresponding value of Q c or A c for the particular case with N sub = 225: where ε Q (N sub ) and ε A (N sub ) are, respectively, the resolution error in Q c or A c at a resolution of N sub , and Q c (225) and A c (225) are, respectively, the reference Q c or A c value at a resolution of N sub = 225. By convention, we take errors as positive when fluxes are overpredicted. It will be later shown that further increasing the resolution above N sub = 225 yields a negligible reduction in resolution error, and thus this case is taken to be the reference. Therefore, by definition ε Q (225) = ε A (225) = 0. The level of heterogeneity in the PDF (i.e. deviation from a uniform distribution) was quantified by calculating the entropy of the distribution where S is the distribution entropy, N bins is the number of discrete bins in the PDF and P(Q i ) is the probability density based on the value of Q in the ith discrete bin. The normalization condition for the probability density is which means that for a uniform distribution P(Q i ) = 1 and S = 0. As the leaf angle distribution tends more towards horizontal leaves, the fraction of leaf area projected in the direction of the sun increases (since the sun direction is vertical). An integral measure of leaf area projection is the Ross-Nilson G-function (Ross 1981), which can be interpreted intuitively as the fraction of leaf surface area projected in the direction of the sun. Since only a vertical sun direction is considered, here the G-function is the fraction of leaf area projected on a horizontal plane and is denoted by G(0). A description of the calculation of G(0) for the four different leaf angle distributions considered is given in the Appendix. Figure 4 characterizes the impact of LAI, leaf angle distribution and diffuse radiation flux on probability distributions of absorbed photosynthetically active radiation flux Q and net CO 2 flux A using the finest sub-leaf resolution (N sub = 225). The impact of LAI is largely manifested in terms of a change in partitioning of sunlit and shaded leaf area, which determines the weighting of the leftmost distribution bin containing Q = 0 relative to all other bins. Increasing LAI increases  Note that a complex leaf shape is used for illustration purposes, and is different from the simplified rectangular geometry shown in Fig. 1.

Characterization of probability distributions across cases
the fraction of shaded leaf area, which in turn increases the probability within the bin containing Q = 0. This variation in the distribution of Q translates into a similar effect in the distribution of A, which is that increasing LAI increases the probability of A being near the value of −R d and scales the distribution downward for discrete bins of A above the first bin. Variation in the leaf angle distribution has an impact on both the relative partitioning of sunlit to shaded leaf area, as well as the entire shape of the distribution. As is to be expected, leaf angle distributions that increase the fraction of leaf area projected in the direction of the sun (in this case the vertical direction) has a similar effect as increasing LAI. The planophile distribution, which has leaf angles tending towards the horizontal, has both a relatively large fraction of shaded leaf area (evidenced by the high probability in the bin containing Q = 0) and a relatively high probability that Q is near the maximal value of 1200 μmol m −2 s −1 . In contrast, the erectophile distribution, which has leaf angles tending towards the vertical, has relatively low probabilities that Q is near 0 or 1200 μmol m −2 s −1 . The spherical distribution yields a probability distribution of Q that is uniform above the bin containing Q = 0, and the uniform distribution, which has leaf angles tending towards horizontal, yields a probability distribution of Q that falls roughly between that of the spherical and planophile distributions.
The probability distributions of A for planophile and uniform leaf angle distributions are similar in shape to the probability distributions of Q, where there is a relatively high probability that A is at the extreme ends of the distribution. The spherical leaf angle distribution yields a highly non-uniform probability distribution of A, although it yields a uniform Q distribution for Q > 0. This is due to the non-linear response of A to Q, as illustrated in Fig. 2. The probability distribution of A is closest to uniform for the erectophile leaf angle distribution, although there is still substantial non-uniformity.
The addition of diffuse radiation substantially changes the shape of the probability distributions of both Q and A. As is to be expected, diffuse radiation reduces the probability of a fully shaded leaf, and for the cases of f diff of 0.1 and 0.2, the most likely Q value is greater than zero. These changes in the distribution of Q due to diffuse radiation translate into even more pronounced differences in the distribution of A. Although the addition of diffuse radiation only shifted the lower peak in the distribution of Q a small amount, non-linearity in the light response of A creates a peak in the distribution of A at intermediate A values (e.g. 25-50 % of A max ) when diffuse radiation is added.

Impacts of sub-leaf resolution on probability distributions
The number of sub-elements per leaf that was used had a strong impact on the probability distribution of Q, as illustrated in Fig. 5A using the default parameter set. When only 1 element per leaf was used, the distribution of Q was nearly uniform, although the 'true' distribution (assumed to be that of N sub = 225) was highly non-uniform. Increasing sub-leaf resolution increased peaks in the distribution on the extreme ends of the distribution. The peak at the left end of the distribution near Q = 0 corresponds to shaded leaf area, which is largely eliminated for N sub = 1 due to averaging across the entire leaf. The peak at the right end of the distribution is created by the uniform leaf angle distribution in which leaves tend towards horizontal orientations. This leaf angle effect is largely eliminated when only 1 element per leaf is used (N sub = 1). There is little visible difference in the distribution of Q between N sub = 100 and N sub = 225, indicating that convergence occurred somewhere around N sub = 100. While the number of sub-elements per leaf had a significant impact on the probability distribution of A (Fig. 5B), these impacts were manifested somewhat differently than for the distribution of Q. Towards the left end of the distribution of A, the impact of sub-leaf resolution was similar as was observed in the distribution of Q, in that using only 1 element per leaf nearly eliminated the peak in the distribution near Q = 0. Increasing N sub revealed a large peak in the distribution of A that was similar to that of the distribution of Q. However, across the rest of the distribution of A, the effect of sub-leaf resolution had a much weaker effect. Even with N sub = 1 , there was a strong peak at the right end of the distribution of A, which only exhibited a moderate decrease as N sub was increased. These results can be attributed to the fact that the light response of A for low values of Q is approximately linear with a steep slope, but shows little variation with Q at high values (Fig. 2). This means that a flattening of the distribution of Q at high values of A will translate into relatively smaller changes in the distribution of A.

Impacts of sub-leaf resolution on errors in whole-canopy fluxes
Sub-leaf resolution had a negligible impact on errors in whole-canopy absorbed radiation flux Q c (Fig. 6A, C and E). The relative error was typically on the order of 10 −3 and varied randomly with N sub . Although the probability distribution of Q was shown to be highly sensitive to sub-leaf resolution (Fig. 5), the total amount of energy absorbed by the canopy as indicated by Q c is relatively insensitive to N sub because Q c does not reflect how that energy is distributed throughout the canopy. The error in whole-canopy CO 2 flux A c , however, was highly sensitive to sub-leaf resolution (Fig. 6B, D and F). Regardless of canopy configuration, the decline in the error ε A with increasing N sub was nearly exponential as indicated by the nearly linear decline in the error on a logarithmic scale. The error at N sub = 1, ε A (1), was highly sensitive to canopy configuration (Fig. 7). The error increased linearly with increasing LAI, increased linearly as the leaf angle distribution tended more towards horizontal (i.e. increasing G and effective LAI), and decreased logarithmically with increasing diffuse radiation fraction f diff . These trends in the error ε A remained relatively consistent as N sub was increased above 1 (Fig. 6B, D and F).

Relationship between radiation entropy and errors in whole-canopy fluxes
The entropy of the probability distribution of light interception Q can be used as a metric for explaining much of the variation in the wholecanopy error in net CO 2 flux ε A (1). Figure 8 plots the variation in the whole-canopy error ε A (1) as a function of the radiation distribution entropy S calculated from the highest resolution case (N sub = 225) for cases with varying LAI and G(0). Results illustrate that there is an approximately linear relationship between S and ε A (1) as LAI and G(0) are varied. As the product of LAI G(0) is reduced, this linear relationship becomes weaker. This reduction in the product of LAI G(0) corresponds to a reduction in shaded leaf area, which corresponds to the large peak in the probability distribution of Q near Q = 0. When diffuse radiation is added (Fig. 8B), the slope of the linear relationship between S and ε A (1) decreases (i.e. lower maximum error), and there is slightly more scatter in the relationship. As with reducing the product of LAI G(0), adding diffuse radiation also decreases the fraction of shaded leaf area.
Variation in radiation entropy as sub-leaf resolution was increased is shown in Fig. 9. As intuitively expected based on the probability distributions shown in Fig. 5, entropy decreased as sub-leaf resolution was increased. In most instances, this decrease was nearly exponential, as indicated by the roughly linear decline in S with log N sub . It was observed that under-resolved simulations caused the probability density of Q to tend towards uniform, which corresponds to an increase in  entropy. Sufficient sub-leaf resolution is needed to fully generate the reduction in distribution entropy that is created by canopy structure.
The erectophile and spherical distributions created relatively weak dependence of radiation entropy on N sub . Close examination of the erectophile probability distribution of Q (Fig. 4C) reveals that the erectophile and spherical leaf angle distributions result in a Q probability distribution that decreases monotonically with increasing Q, whereas the other distributions are closer to U-shaped. The U-shaped distribution caused large reductions in the probability of Q at both ends of the distribution when N sub was decreased, resulting in large increases in distribution entropy. For the erectophile and spherical distributions, which only had a large peak in probability near Q = 0, the increase in entropy with decreasing N sub was smaller.

DISCUSSION
Although it is often implicitly assumed that leaf-scale resolution is adequate for representing biophysical processes in 3D plant models, the results above demonstrate the critical importance of sub-leaf resolution in representing non-linear processes. While the accuracy of wholecanopy radiation interception prediction is essentially independent of sub-leaf-scale resolution (Fig. 6A, C and E), it plays an important role in accurate prediction of how the radiation flux is distributed across the canopy (Fig. 5). When the computed leaf-level radiation flux is input to other non-linear sub-models, substantial errors can result (e.g. >100 %) if sub-leaf resolution is insufficient.
As was illustrated in the previous work of Kent and Bailey (2021), an under-resolved leaf-resolving simulation can have much larger errors in whole-canopy fluxes than a statistical 'big leaf ' model. The present study reveals more clearly why this is the case. Inability to resolve subleaf-scale heterogeneity due to shadows reduces the amount of shaded leaf area, which effectively leads to a homogenization of the probability distribution of absorbed radiation. When coupled with a non-linear model, this can significantly affect prediction of whole-canopy fluxes. Statistical or 'turbid media' models do not seek to resolve these shadows explicitly, but rather represent them statistically by dividing leaf area into sunlit and shaded fractions (e.g. DePury and Farquhar 1997; Wang and Leuning 1998;Sinoquet et al. 2001). While the traditional sunlit-shaded statistical approach can have issues depending on the shape of the coupled non-linear model light response (Kent and Bailey 2021), they can perform better than a leaf-resolving model in a homogeneous canopy if sub-leaf resolution is insufficient.
While there have been hundreds of studies using 3D plant models that couple a model for radiation interception with nonlinear sub-models (e.g. photosynthesis, stomatal conductance, energy balance; see Vos et al. 2010;Evers et al. 2018, for reviews), few if any explicitly discuss the choice of sub-leaf resolution. Often, the number of elements per leaf used is either not specified by the authors or chosen to be 1 (e.g. Allen et al. 2005; Costes   Sarlikioti et al. 2011;Barillot et al. 2019). When multiple elements per leaf are used, the number of elements is generally chosen based on desire to represent the perimeter shape or curvature of the leaf, and is usually on the order of 10 elements per leaf (e.g. Kahlen et al. 2008;Song et al. 2013;Schmidt and Kahlen 2018). The results presented herein suggest a general need in this field to more rigorously consider and evaluate the choice of subleaf resolution and its impact on predicted fluxes, and ultimately growth and development.
In recent years, much effort has been devoted towards development of methods for accurate measurement of plant canopy structure, in order to provide robust geometrical inputs for 3D models. These include manual electronic measurement of the position and orientation of all leaves in the plant (Smith and Curtis 1995;Sinoquet et al. 1998), photogrammetry-based measurement of whole-canopy structure (Shlyakhter et al. 2001;Phattaralerphong and Sinoquet 2005;Pound et al. 2016) and LiDAR-based reconstruction approaches (Delagrange and Rochon 2011;Côté et al. 2012;Bailey and Ochoa 2018). However, such efforts may be wasted if the 3D simulations are underresolved. Results indicated that insufficient resolution creates a diffusive effect that tends the probability distribution of radiative flux towards uniform (i.e. artificially increases entropy). Decreases in entropy are a direct result of canopy structure, and thus insufficient resolution smears out the effects of canopy structure in the resulting probability distribution. While the impacts on whole-canopy radiation interception are negligible, insufficient resolution strongly impacts the distribution of intercepted radiation, which can create substantial errors when coupled to other non-linear biophysical models (Figs 6 and 7).
Using a resolution of hundreds of elements per leaf introduces a substantial challenge related to computational expense. Whole-plant and -canopy simulations at the leaf scale are already quite computationally expensive, and often limit domain sizes to one or a few plants depending on their size. Using O(100) elements per leaf effectively increases computational expense by 100 times, thus severely limiting the possible range of scales that can be feasibly simulated, and potentially putting accurate simulation of even a single large plant (e.g. tree) out of reach for many models. Helios has the capability to represent O(10 8 ) elements in a single simulation based on current GPU technology, which corresponds to roughly O(10 6 ) leaves with 100 elements per leaf and potentially a dozen or more large trees. More efficient models and/or improved computational hardware may still be needed to accurately simulate plant-scale heterogeneity at the field scale for some applications. The approximately exponential decrease in error with increasing sub-leaf resolution shown here and by Kent and Bailey (2021) fortunately means that even modest increases in sub-leaf scale resolution short of O(100) elements per leaf can substantially reduce errors. It is also acknowledged that sub-leaf-scale resolution is one of many considerations when using leaf-resolving plant models, and depending on application, may be outweighed by other sources of error (e.g. uncertainty in photosynthesis model parameters or environmental inputs).
These results related to resolution dependence are not unique to plant modelling, and indeed such attention towards the impact of resolution is typical when considering whether to use a continuum (statistical) or fully resolved model. For example, performing 'turbulence-resolving' simulations at a grid resolution insufficient to fully resolve all turbulent scales can create errors much larger than using a simple coarse resolution simulation that represents smallscale turbulence statistically (e.g. Komen et al. 2017). In the field of particulate dispersion, employing a Lagrangian modelling viewpoint that tracks each individual particle can introduce non-physical behaviour with errors much larger than an equivalent Eulerian averaged model that represents particulates statistically if model resolution is insufficient (e.g. Bailey 2017). Many other examples can be found in other fields, such as in molecular dynamics simulation and smooth particle hydrodynamics model applications.
The results of this study are clearly limited to only a few canopy configuration and model parameter combinations. While this means that quantitative error values for a specific case likely cannot be extrapolated to any particular simulation outside of this study, the trends revealed by the simulation experiment provide valuable information that can help qualitatively guide resolution choice in future work. Canopy configurations that reduce radiation entropy such as high effective LAI (LAI G), high horizontal homogeneity or low scattering are likely to be more sensitive to discretization errors. Although not explicitly included in this study, Kent and Bailey (2021) also showed that penumbral effects, radiation scattering and the shape of the photosynthetic response to light (governed by α) can reduce discretization errors, which this study has revealed is due to their smoothing effect on radiative gradients and associated increases in the entropy of the distribution of Q or A. The effect of leaf size was also investigated by Kent and Bailey (2021), which had a relatively small effect on averaging errors. This is likely because the leaf grid size decreased proportionately with leaf size and thus shadow size.

CON CLUSIONS
Accurate simulation of whole-canopy fluxes with a non-linear response to radiation using a leaf-resolving model requires that the model resolves shadows at the sub-leaf scale. This requires a relatively fine discretization of fluxes on individual leaves. Using only one to a few discrete elements to represent leaf-scale fluxes had a negligible impact on whole-canopy fluxes of absorbed radiation, but could create errors in whole-canopy photosynthesis in excess of 100 % for the limited set of simplified example cases considered. Errors due to sub-leaf-scale resolution decreased exponentially as the number of discrete sub-leaf elements was increased, such that convergence was achieved for >O(100) elements per leaf. Canopy configurations that increase radiative gradients and thus decrease radiation entropy are prone to higher discretization errors, and any process that increases radiation entropy is likely to decrease errors (e.g. penumbral effects, scattering, diffuse radiation, leaf flutter). Future work should more closely consider the role of sub-leaf resolution on model errors, which is likely to prompt an increase in resolution relative to current practice. While it may not be possible to produce generalizable a priori relationships for resolution errors given the wide range of potential models, assumptions and canopy architectures, resolution sensitivity analysis is a feasible means for assessing the chosen resolution for a given simulation. As is common practice in many model-based fields, resolution can be progressively refined for a select number of cases to verify that further refinement has a negligible impact on results.

A PPENDI X : C A LCUL ATION OF THE G -FUN CTION
The G-function is the fraction of leaf area projected in a given direction (usually the direction of the sun), and is defined mathematically as G( n s ) = 1 2π 2π 0 π/2 0 g L ( n L )| n s · n L |dθ L dφ L , where G( n s ) is the G-function in the direction of the sun, g L ( n L ) is the probability that the leaf normal vector is oriented in the direction n L, n s is a unit vector pointing in the direction of the sun, θ L is the zenithal component of the leaf normal vector n L and φ L is the azimuthal component of the leaf normal vector. In this work, only azimuthally symmetric leaf angle distributions were considered, and the sun was always directly overhead. In this case, G is independent of φ L and corresponds only to a solar zenith of θ s = 0. This value of the G-function is denoted as G(0) and can be defined mathematically through simplification of Eq. (6) G(0) = π/2 0 g L (θ L ) cos θ L dθ L .
Values of g L and G(0) for the four different leaf angle distributions considered in this study are tabulated in Table 2.

SO URCE S OF FUNDIN G
This research was supported by the USDA NIFA (Hatch project no. 1013396), and the Almond Board of California (project no. PREC1).