# Finite element analysis of trees in the wind based on terrestrial laser scanning data

Toby Jackson, Alexander Shenkin, Axel Wellpott, Kim Calders, Niall Origo, Mathias Disney, Andrew Burt, Pasi Raumonen, Barry Gardiner, M. Herold, et al.

► To cite this version:

Toby Jackson, Alexander Shenkin, Axel Wellpott, Kim Calders, Niall Origo, et al. Finite element analysis of trees in the wind based on terrestrial laser scanning data. Agricultural and Forest Meteorology, 2019, 265, pp.137-144. 10.1016/j.agrformet.2018.11.014. hal-02089532

HAL Id: hal-02089532

https://hal.science/hal-02089532

Submitted on 22 Nov 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Finite element analysis of trees in the wind based on terrestrial laser scanning data

<figure><img src="image_3.png" /></figure>

T. Jackson<sup>a,*</sup>, A. Shenkin<sup>a</sup>, A. Wellpott<sup>b</sup>, K. Calders<sup>c</sup>, N. Origo<sup>d,e</sup>, M. Disney<sup>e,f</sup>, A. Burt<sup>e</sup>, P. Raumonen<sup>g</sup>, B. Gardiner<sup>h</sup>, M. Herold<sup>i</sup>, T. Fourcaud<sup>j</sup>, Y. Malhi<sup>a</sup>

## ARTICLE INFO

### Keywords:

TLS

Terrestrial laser scanning

Wind damage

Finite element analysis

Critical wind speed

Resonant frequency

## ABSTRACT

Wind damage is an important driver of forest structure and dynamics, but it is poorly understood in natural broadleaf forests. This paper presents a new approach in the study of wind damage: combining terrestrial laser scanning (TLS) data and finite element analysis. Recent advances in tree reconstruction from TLS data allowed us to accurately represent the 3D geometry of a tree in a mechanical simulation, without the need for arduous manual mapping or simplifying assumptions about tree shape. We used this simulation to predict the mechanical strains produced on the trunks of 21 trees in Wytham Woods, UK, and validated it using strain data measured on these same trees.

For a subset of five trees near the anemometer, the model predicted a five-minute time-series of strain with a mean cross-correlation coefficient of 0.71, when forced by the locally measured wind speed data. Additionally, the maximum strain associated with a 5 ms⁻¹ or 15 ms⁻¹ wind speed was well predicted by the model (N = 17, R² = 0.81 and R² = 0.79, respectively). We also predicted the critical wind speed at which the trees will break from both the field data and models and find a good overall agreement (N = 17, R² = 0.40). Finally, the model predicted the correct trend in the fundamental frequencies of the trees (N = 20, R² = 0.38) although there was a systematic underprediction, possibly due to the simplified treatment of material properties in the model. The current approach relies on local wind data, so must be combined with wind flow modelling to be applicable at the landscape-scale or over complex terrain. This approach is applicable at the plot level and could also be applied to open-grown trees, such as in cities or parks.

## 1. Introduction

Wind damage to forests is a significant but poorly understood driver of the terrestrial carbon cycle (Espírito-Santo et al., 2014). Trees respond to strong winds by leaning and swaying in a dynamic manner. This response has been characterized using static models (Hale et al., 2015), dynamic models with simplistic tree geometry (Kerzenmacher and Gardiner, 1998) and dynamic models which represent the effect of branches (James et al., 2014). The literature has mostly focused on

wind damage in conifer plantations due to their importance to the forestry industry (Dupont et al., 2015; Hale et al., 2015, 2012). In broadleaf forests the variety of tree sizes and shapes influence their response to wind forcing, meaning that finite element analysis is necessary to model the detailed effects of wind on trees (James et al., 2014; Rodriguez et al., 2008; Spatz and Theckes, 2013). Previous finite element analysis studies have been limited to small numbers of trees due to the difficulty of manually mapping tree geometry (Moore and Maguire, 2008; Sellier et al., 2006). However, terrestrial laser scanning

<footnote><sup>a</sup> Environmental Change Institute, School of Geography and the Environment, University of Oxford, OX1 3QY, UK</footnote>

<footnote><sup>b</sup> Facility for Airborne Atmospheric Measurements, Wharley End, Bedford, MK43 0AL, UK</footnote>

<footnote><sup>c</sup> CAVELab - Computational & Applied Vegetation Ecology, Ghent University, Belgium</footnote>

<footnote><sup>d</sup> Earth Observation, Climate and Optical Group, National Physical Laboratory, Hampton Road, Teddington, Middlesex, TW11 0LW, UK</footnote>

<footnote><sup>e</sup> Department of Geography, University College London, WC1E 6BT, UK</footnote>

<footnote><sup>f</sup> NERC National Centre for Earth Observation (NCEO), UK</footnote>

<footnote><sup>g</sup> Tampere University of Technology, Korkeakoulunkatu 10, 33720, Tampere, Finland</footnote>

<footnote><sup>h</sup> EFI Planted Forest Facility, 69 route d'Arcachon, 33612, Cestas, France</footnote>

<footnote><sup>i</sup> Laboratory of Geo-Information Science and Remote Sensing, Wageningen University, Wageningen, the Netherlands</footnote>

<footnote><sup>j</sup> AMAP, University of Montpellier, CIRAD, CNRS, INRA, IRD, Montpellier, France</footnote>

<footnote>* Corresponding author.</footnote>

<footnote>E-mail address: tobydjackson@gmail.com (T. Jackson).</footnote>
(TLS) data offers the potential for dramatic advances in this area. Recent developments in the extraction of accurate 3D tree models from TLS data can provide detailed geometrical information for a scanned tree in a matter of minutes (Raumonen et al., 2015; Calders et al., 2015; Gonzalez de et al., 2018).

We used 3D tree models derived from TLS data as the basis of a finite element analysis, applied using the engineering software Abaqus (Simulia Software Company, 2017). This approach allowed us to simulate the dynamics of broadleaf trees without the need to manually map their geometry. In addition to its application in natural forests, this approach could be useful in testing proposed interventions on urban amenity trees, or valuable veteran trees, in a model environment before implementation (Reiland et al., 2015).

This paper tests the hypothesis that finite element models developed using TLS data can be successfully used to simulate the behaviour of broadleaved trees under applied wind loading. We describe the method and discuss the assumptions involved. We also provide software to implement this technique. We test the model against field data collected in Wytham Woods, UK, for which we have both TLS data and field measurements of strain in the tree stem and wind speed. This validation takes three forms:

1 We predict the exact time-series of strain produced during a five-minute period based on high time-resolution wind data.

2 We calculate the maximum hourly strain produced on a tree given maximum hourly wind speeds of 5 ms⁻¹ or 15 ms⁻¹. We also estimate the critical wind speeds (CWS) at which trees will break by extrapolating from long term field data (Hale et al., 2012).

3 We predict the fundamental frequencies of the trees and compare these to the frequencies extracted from the field data using Fourier analysis.

Finally, we discuss the results of this study in the context of other work in the field.

## 2. Materials and methods

The workflow involved field data collection and analysis, TLS data collection and modelling work, followed by three separate tests of the model against the field data (Fig. 1).

### 2.1. Study site and field data

Wytham Woods (51°46'27.2" N 1°20'20.1" W) is a mature temperate deciduous woodland near Oxford, UK. The woodland contains approximately 950 trees per hectare and the dominant species are sycamore (Acer pseudoplatanus), ash (Fraxinus excelsior), beech (Fagus sylvatica) and English oak (Quercus robur). We measured the bending strains produced on 21 trees (eight ash, nine sycamore, four birch (Betula spp.) and one oak) ranging from 12 to 24 m tall. One of the

<figure><img src="image_4.png" /><figcaption>Fig. 1. Workflow diagram showing the main analysis and validation steps.</figcaption></figure>

sycamore trees had four main stems below 1.3 m, each of which are treated as separate trees in the forest census and were instrumented as such in this study, we therefore refer to them as separate trees in this paper. Seven of the trees were located near a canopy walkway structure, and the remaining 14 trees were located 230 m away. Average canopy height was 21–23 m in both locations.

Strain was measured using the method developed by Moore et al., (2005): pairs of strain transducers were fixed to the trees at approximately 1.3 m and perpendicular to each other. These transducers were then connected to Campbell Scientific (Logan, UT, USA) data loggers and measured at 4 Hz (two CR23xs and three CR1000s). The limited storage capacity of the CR23x data loggers was supplemented using a Raspberry Pi (Cambridge, UK).

In the field, strain gauge resistances drift with environmental factors. We partially corrected these in situ using a variable resistor within a wheat-stone bridge circuit. We further accounted for this drift by subtracting a 5000 point (20 min) running mode from the signal. This is based on the assumption that the tree will regularly pass through its equilibrium position. This same technique was used in Wellpott (2008) and Hale et al. (2012, 2015) and we found our results to be relatively insensitive to this window length (SI-1).

Wind speed was measured in winter using three cup anemometers (A100LK/5 M, Vector Instruments, Rhyl, Wales) mounted at 5 m, 10 m and 15 m heights above ground level on a canopy walkway next to the focal trees and supplemented, notably during the 5-minute high-resolution validation period, with a 3D sonic anemometer (Campbell Scientific CSAT3) mounted at 15 m. In summer, we used a Gill (Gosport, UK) Sonic-1 anemometer recording at 1 Hz mounted at 15 m. Local climate data are available through a meteorological station approximately 1 km away, just outside the forest, operated by the Environmental Change Network (ECN, www.ecn.ac.uk). This station records mean hourly wind speed as well as the maximum 5 s mean within that hour (5 s maximum gust speed). The wind data taken from anemometers within the canopy were considered unrepresentative, especially in summer, because of localized sheltering effects. We therefore used the ECN meteorological station's hourly maximum wind speeds for all critical wind speed analysis and select the maximum strain on each tree within each hour. A comparison of the different measures of wind speed is given in SI-2.

### 2.2. Development of quantitative structure models

TLS data were collected in Wytham Woods, Oxford, as part of a separate project (Calders et al., 2018, 2016). Trees were manually extracted from the point cloud and cylinder fitting algorithms, developed primarily for biomass estimation, were used to extract the tree shapes (Åkerblom, 2017; Burt, 2017). This method represents each tree as a series of cylinders – a quantitative structure model (QSM).

We adapted these QSMs by removing cylinders under 2 cm diameter, which are generally the least accurate in the QSM reconstruction due to limitations of TLS scanning resolution (see Fig. 2). We then longitudinally averaged neighbouring cylinders, replacing each pair of neighbouring cylinders with a single cylinder with the average length, radius and direction of the original pair (Åkerblom, 2017). This reduced the number of cylinders from approximately 4000 per tree to 500, without significantly affecting the total woody volume of the tree. This process resulted in a new set of QSMs comprised of cylinders with a higher length to radius ratio, allowing us to use beam theory to build a dynamic model of the tree without significantly altering the overall structure (Carrera et al., 2011). It is important to note that the Wytham Woods TLS data were collected in ideal conditions based on a 10 m grid using a Riegl (Horn, Austria) VZ-400 with 0.04° angular resolution (Wilkes et al., 2017). Most of the trees were scanned in winter, but those near the walkway were scanned in summer and the leaves removed through filtering (Vicari, 2017). Optimization of the QSM approach shows that there is very little variance between the QSMs in
<figure><img src="image_5.png" /><figcaption>**Fig. 2.** Left – QSM directly from the fitting process. Right – simplified QSM with longer cylinders and fewer small branches - this is our input for finite element analysis.</figcaption></figure>

these conditions (Burt, 2017).

## 2.3. Development of finite element models

To build a finite element model in Abaqus, each cylinder in the QSM was represented as a two node Euler-Bernoulli beam element (B31) and the orthotropic features of wood were neglected (Sellier et al., 2006). Since we intend this approach to be applicable wherever TLS data is available, we did not measure the material properties of the trees individually. Values of green wood density, $ρ_{gw}$, elasticity, $E_{gw}$, and modulus of rupture, $MOR_{gw}$, were taken from the literature (Lavers, 1983; Niklas and Spatz, 2010) (Table 1). Breaking strain was defined as $ε_{break} = MOR_{gw}/E_{gw}$ (Niklas and Spatz, 2010). The intra-specific variation in elasticity is between 14 and 20% (defined as standard deviation divided by mean value) in sycamore, ash and oak. No estimates of variation in green wood density were given. The inter-specific variation in wood density is 20% and elasticity 25% in Wytham Woods (defined as the range divided by the minimum value, Table 1). The values used in finite element analysis are given in Table 1.

Root anchorage was modelled as a universal joint with two rotational degrees of freedom and symmetric elasticity $E_{root} = 5$ MPa in each direction. Energy dissipation at the root-soil boundary was modelled by a dashpot with coefficient $C_s = 10$. Previous authors used tree pulling data to tune these model parameters, which gives accurate results in the case of a few trees (Sellier et al., 2006; Sellier and Fourcaud, 2009). However, we are attempting to generalize to many trees without intensive field measurements, so the root parameter values were chosen because the model is insensitive to these parameters in this range (Table 2). Roots are therefore included in the model primarily so that future studies can tune these parameters if the necessary data are available.

Applying gravity in a realistic way proved difficult, since the trees were scanned in their gravity deformed positions. Branches develop reaction wood and grow in irregular shapes to counteract the force of gravity (Scurfield, 1973), whereas our model assumes cylindrical branches with isotropic material properties. We tested three approaches

to modelling the effect of gravity and found a small but interesting difference between them (see SI-3). For all analysis in the main text of this paper we first applied an inverse gravity force to the model tree, extracted the nodal displacements and then restarted the analysis applying a realistic gravitational force to the displaced geometry. This is the most physically realistic approach since branches are not displaced before the analysis begins and the additional load due to self-weight of the displaced crown is accounted for.

## 2.4. Dynamic simulation

The dynamic simulation of a tree's response to wind forcing was conducted using the Abaqus Aqua module (Simulia Software Company, 2017). Each tree was modelled individually, meaning that the dynamic response does not include crown clashing which is thought to be significant in some forest structures (Rudnicki et al., 2003, 2001). In all simulations, we used a non-linear solution method to account for large displacements in the tree structure. The equation of motion takes the form

$$M\ddot{q} + D\dot{q} + Kq = G + F_D \quad (1)$$

where $q$, $\dot{q}$ and $\ddot{q}$ are the nodal displacement, velocity and acceleration vectors. M, D and K are the mass, damping and stiffness matrices respectively, and G and $F_D$ represent the gravitational and wind induced forces. The output variable, strain, is calculated internally by Abaqus and we use the strain at the approximate position of the strain gauges (1.3 m), so as to correspond with the field data.

While horizontal wind velocities were fully represented, the vertical profile of wind speed was limited to a power law fit (Eq. 2), using field data from three heights within the canopy (Figure SI-2).

$$v = v_0 \left( \frac{h}{h_{ref}} \right)^{\gamma} \quad (2)$$

where $v$ is the wind speed at a given height $h$, $v_0$ is the reference wind speed at reference height $h_{ref}$ and $\gamma$ is a variable parameter determined by the local wind profile (Figure SI-2). The wind drag force, $F_D$, per unit length on each cylinder of the QSM is given by

$$F_D = Amp \cdot \rho_{air} C_D r \Delta V \Delta V_{fn} \quad (3)$$

where Amp is a user defined amplitude that describes the horizontal variation in wind speed, $ρ_{air}$ is the density of air, $C_D$ is the drag factor, $r$ is the radius of the beam, and $\Delta V_{fn}$ is velocity difference between the wind and the beam. As in previous work (Ciftci et al., 2013; Moore and Maguire, 2008; Sellier et al., 2006; Sellier and Fourcaud, 2009), we assumed a constant drag factor, although in fact the drag factor varies with wind speed and with the size of the branch.

The most significant source of damping, that due to aerodynamic drag, was modelled implicitly in the Aqua module. Material damping was modelled using Rayleigh's hypothesis (Rayleigh, 1877).

$$D = \alpha_w M + \beta K \quad (4)$$

where $α_w M$ is the mass proportional damping and $βK$ is the stiffness proportional damping. Following Sellier et al. (2006) we assumed stiffness proportional damping would be insignificant in the high

**Table 1**

Material properties of green wood used in finite element analysis. Values in brackets are standard deviation and number of samples tested. Quercus rubra was used as a substitute for Quercus robur, since data for the latter was not available from the same testing method. This affects only a single tree in the frequency analysis validation.

<table><thead><tr><th>Latin Name</th><th>Common Name</th><th>Density (kgm<sup>-3</sup>)</th><th>Elasticity (MNm<sup>-2</sup>)</th><th>MOR (MNm<sup>-2</sup>)</th></tr></thead><tbody><tr><td>Acer pseudoplatanus</td><td>Sycamore</td><td>721</td><td>8400 (1340, 43)</td><td>66 (8.8, 43)</td></tr><tr><td>Fraxinus excelsior</td><td>Ash</td><td>801</td><td>9500 (1880, 204)</td><td>66 (10.9, 204)</td></tr><tr><td>Quercus robur</td><td>English Oak</td><td>865</td><td>10500 (1280, 40)</td><td>72 (7.9, 40)</td></tr><tr><td>Betula spp.</td><td>Birch</td><td>801</td><td>9900 (410, 6)</td><td>63 (5.1, 6)</td></tr></tbody></table>
<figure><img src="image_6.png" /><figcaption>**Fig. 3.** Right - TLS data showing trees near the canopy walkway. The coloured trees are the ones used in this analysis. The black arrow shows the prevailing wind direction and is located at the point of the wind measurement. Left - simulated strain (coloured lines) against field data (black lines, inverted for clarity). The time-lag (field strain - simulated strain) as well as the cross correlation co-efficient are given for each signal.</figcaption></figure>

<figure><img src="image_7.png" /></figure>

Reynolds number regime in question, and so set $\beta = 1 \times 10^{-3}$. Contact between branches was allowed but no attempt was made to accurately model the energy transfer of these contacts.

Multiple simulations were carried out to test the sensitivity of the model. All continuous parameters were tested at ± 20%. Non-continuous parameters, such as the simplification of the QSM, were tested at different levels.

## 2.5. Validation techniques

All the strain data used for validation were collected in winter, when the absence of leaves makes the modelling and interpretation simpler. Future work will parameterize the effects of leaves, but the change in Reynolds number and drag coefficient make the full treatment of wind flow over leaves highly complex (de Langre, 2008).

### 2.5.1. Time series validation

For the time-series validation, we use a five-minute subset of data from 14:00-14:05 on 28th March 2016, chosen for high wind speeds and high (10 Hz) wind measurement resolution. Two trees were removed due to problems with the field data, in the first case the balancing circuit failed and in the second the physical connection to the tree came undone, leaving five trees with reliable data. We therefore compared measured and modelled strain time-series for five trees close to the anemometers. We did not expect an exact correspondence for the 5-minute time-series, since the wind forcing applied to each model tree was identical and reflects the wind regime measured at the position of the canopy walkway. In the field, each individual tree experiences a slightly different wind forcing due to its position (Kerzenmacher and Gardiner, 1998). A previous finite element analysis study found a reduced correlation at 12 m from the wind measurements, although no correlation statistics were given (Sellier et al., 2008).

### 2.5.2. Hourly maximum strain validation

In this validation, we used the five trees from time series validation

as well as twelve additional trees situated approximately 230 m away. Two of these additional trees were excluded due to a poor connection between the tree and the strain gauge. We used hourly maximum wind speed measurements taken from the ECN meteorological station. For each hour we extracted the maximum strain value for every tree. All but four trees had over 500 h of data available, and for these four trees we had 165 h of data. We tested a robust maximum selection technique, using the mode of the maxima of subsections of the hours data, but found no increase in coefficient of determination nor reduction in root mean squared error. Details of the strain data processing and reference to the software and data repository to repeat this analysis are given in SI-1. We regressed the hourly maximum strain data against the squared hourly maximum wind speed and so predicted the strain at 5 and $15 \text{ ms}^{-1}$ wind speeds. We used this regression line to extrapolate the relationship up to the breaking strain, giving an estimate of the critical wind speed from the field data. In the simulation, we artificially increased the wind forcing until the lowest beam in the trunk reached its breaking strain.

### 2.5.3. Fundamental frequency validation

We used a Welch's power spectral density function (Krauss et al., 1994) to extract the fundamental frequency for each hour of field data during winter and take the mean. This analysis only requires one strain gauge per tree, so three of the four strain gauge failures had no effect here, meaning that we can carry out this analysis for 20 trees. Using Abaqus, we predicted all the undamped sway modes from the QSMs using a subspace method (Viberg, 1995). This can be performed separately from the dynamic analysis and is much faster, it is also available in the free version of Abaqus.

## 3. Results

### 3.1. Time-series validation

The modelled five-minute time-series of strain was similar to that
<figure><img src="image_8.png" /><figcaption>**Fig. 4.** Strain produced by wind speeds of 5 ms⁻¹ (A) and 15 ms⁻¹ (B) from field data and finite element analysis. C – critical wind speed estimates from field data and finite element analysis. Pearson's correlation coefficient and the concordance correlation coefficient (CCC) are given for all figures. Trees with an outer black circle are situated near the walkway, approximately 230 m from the other trees and were scanned in summer. The four trees that are independent stems of a large sycamore (branching below 1.3 m) tend to cluster together and are marked with a central black point (two stems overlap in panels A and B).</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Fig. 5.** Field measured fundamental frequency against modelled fundamental frequency for 20 trees (one tree was excluded due to errors in the simulation). The goodness of fit is highly sensitive to the two birch trees for which the fundamental frequency was overestimated. Removing these two trees increases the goodness of fit (R² to 0.66, CCC = 0.83).</figcaption></figure>

measured in the field (Fig. 3, mean cross-correlation = 0.71). We found the lag between measured and simulated strain was partially explained by the position of the trees in relation to the wind measurements. However, all modelled trees displayed some peaks that were not evident in the field data, and vice-versa. This was presumably due to differences in wind at the anemometer position and at the position of the tree. This time-series validation is similar to that given in Moore and Maguire (2008) and Sellier et al. (2008). It is useful to demonstrate that the modelled response resembles that of the actual trees, but only applicable at short distances from the anemometers due to the spatial variability of wind speeds.

## 3.2. Hourly maximum strain validation

This section examines the model's accuracy in predicting the strains produced by wind speeds of 5 and 15 ms⁻¹, and in extrapolating up to

the critical wind speed (Fig. 4). It is important to note that, in the case of CWS, both the x and y axes are estimates, since the CWS was estimated from field data by extrapolating beyond the available measurements. Therefore, this is not a straightforward validation but rather a demonstration that predictions from field and model approaches are correlated. We used hourly maximum wind speeds (highest 5 s mean per hour) measured at the nearby meteorological station, since this avoids the need for a multiplication factor (gust factor). A more detailed discussion of the data analysis choices relating to the wind data is given in SI-2.

The model predicted the strain produced at 5 ms⁻¹ with an R² = 0.81 and concordance correlation coefficient (CCC) of 0.91. At 15 ms⁻¹ wind speed the fit was slightly reduced R² = 0.79 and CCC = 0.89. CWS predictions from field data and the model displayed a worse agreement, R² = 0.40 and CCC = 0.66, presumably because the model behavior departed from a pure square law before reaching the CWS (due to streamlining), while the best fit line derived from the field data does not. We saw no departure from a square law in our field data, but we did not have measurements at high enough wind speeds to test this.

## 3.3. Fundamental frequency validation

We found that finite element analysis predicts the trend of fundamental frequency ($f_0$) correctly (R² = 0.38 and CCC = 0.65, N = 20). Using sample mean material properties, instead of the species-specific material properties which are often unavailable, only decreased the R² value by 0.05. Interestingly, there is a systematic under-prediction of $f_0$ across the 20 trees (Fig. 5). This is likely due either to our simplistic treatment of material properties, which are constant throughout the tree, or our arbitrary trimming of all branches with a diameter under 2 cm. As noted in the figure legend, the fit is radically improved by removing the two birch trees whose fundamental frequencies were over predicted.

The free version of Abaqus can perform this analysis for any tree model consisting of fewer than 1000 beams. As previously noted this method predicts a dense distribution of modes about the main frequencies (Rodriguez et al., 2012), representing the slightly different shapes or directions in which the tree can sway. In the field data, this takes the form of a wide resonance peak, instead of distinct modes.
<caption>**Table 2**</caption>

List of parameters used in simulations with sensitivity analysis results. Continuous parameters were increased by 20%, non-continuous parameters were varied as indicated in brackets (→ test value).

<table><thead><tr><th rowspan="2"></th><th rowspan="2">Parameter</th><th rowspan="2">Default<br>(→ test)</th><th colspan="3">Sensitivity – % increase for 20%<br>increase in parameter</th></tr><tr><th>Strain at<br>15 ms<sup>-1</sup></th><th>CWS</th><th>f<sub>0</sub></th></tr></thead><tbody><tr><td rowspan="2">QSM</td><td>Averaging iterations</td><td>1 (→ 2)</td><td>43</td><td>-14</td><td>-20</td></tr><tr><td>Beam type</td><td>B31<br>(→ B33H)</td><td>-7</td><td>3</td><td>1</td></tr><tr><td rowspan="3">Wood<br>properties</td><td>Green wood<br>density, ρ<sub>gw</sub></td><td>Species<br>specific</td><td>7</td><td>-2</td><td>-12</td></tr><tr><td>Elasticity, E<sub>gw</sub></td><td>Species<br>specific</td><td>-20</td><td>12</td><td>12</td></tr><tr><td>Damping factor, α<sub>W</sub></td><td>0.2</td><td>0</td><td>0</td><td>0</td></tr><tr><td rowspan="2">Root<br>model</td><td>Root elasticity,<br>E<sub>root</sub></td><td>5e6 GPa</td><td>0</td><td>0</td><td>0</td></tr><tr><td>Dashpot coefficient</td><td>10</td><td>0</td><td>0</td><td>0</td></tr><tr><td rowspan="2">Aero-<br>dynamics</td><td>Drag coefficient, C<sub>D</sub></td><td>0.82</td><td>23</td><td>-9</td><td>0</td></tr><tr><td>Wind decay, γ</td><td>0.05</td><td>-1.5</td><td>1</td><td>0</td></tr></tbody></table>

### 3.4. Sensitivity analysis

In this section, we examine the sensitivity of the model to the input parameters. We used three output variables as indicators of the models' response to a prescribed change in input parameters, the strain at a wind speed of 15 ms⁻¹, the critical wind speed estimate and the fundamental frequency (Table 2). The model was highly sensitive to the details of the QSM used as input. We found that the simplification step, used to reduce the number of cylinders in the original QSM for finite element analysis, has a large effect on the result. The second largest sensitivity was due to material properties, especially wood elasticity. Natural variation in wood elasticity, where data are available, is between 14% and 19%, which is comparable to the range of our sensitivity analysis. Unfortunately, green wood elasticity is not often measured and different measurement techniques are often not directly comparable (Ruel et al., 2010).

The beam type is a property of the finite element analysis and determines the solution method used to calculate the strains in each beam (Simulia Software Company, 2017). The choice of beam type is one of the smaller sources of sensitivity in the model. The drag factor, C<sub>D</sub>, is an experimentally determined parameter which determines how much force is transferred to an object for a given wind speed (see Eq. 3). We treated the drag factor as a constant, following (Ciftci et al., 2012; James, 2010; Moore and Maguire, 2008; Sellier and Fourcaud, 2009), whereas in fact it changes with the Reynolds number of the flow (de Langre, 2008). These changes are small in the low wind speed regime of our measurements and several projects are currently investigating the changes in the drag factor for trees at high wind speed. The model is insensitive to root parameters in this range, which are poorly constrained in the literature.

### 4. Discussion

Overall, finite element analysis using QSMs from TLS data is a useful method for modelling the behaviour of trees in the wind (Table 3). It can predict the time-series of strain produced in the trunk of a tree, if

the model is forced by wind speed data measured nearby. It can also predict the maximum strain for a given maximum wind speed, the critical wind speed at which the tree will break and the fundamental frequency of a tree. Since TLS data are becoming increasingly available for trees from a wide range of environments from tropical forests to city parks, this technique could now be used to model the dynamic behavior of a wide range of trees. This will allow us to revisit questions about the relationship between tree architecture and dynamics with realistic tree shapes instead of theory-based fractal or modelled trees. In particular, we could explore the relative importance of variations in material properties and tree architecture (Sellier and Fourcaud, 2009) and the differences between open grown trees and those in forests (James et al., 2014; Kane et al., 2014).

The characteristics of the high time-resolution response of a tree to wind are important, since it is at these scales that dynamical amplification and damping processes occur, which alter the likelihood of wind damage (Ciftci et al., 2013; James et al., 2006; Spatz et al., 2007). However, wind speeds are highly spatially variable, so an array of anemometers is needed to represent the wind forcing. The variability of wind speeds can broadly be described by coherent gusts of approximately one tree height in the cross-stream direction and significantly longer in the stream-wise direction (Finnigan, 2000). In the current

<caption>**Table 3**</caption>

Summary of validation results.

<table><thead><tr><th>Validation type</th><th># Trees</th><th>Fit statistics</th></tr></thead><tbody><tr><td>Time-series validation</td><td>5</td><td>Mean x-corr = 0.71</td></tr><tr><td rowspan="3">Hourly maxima</td><td>5 ms<sup>-1</sup></td><td>R<sup>2</sup>=0.81, CCC=0.91</td></tr><tr><td>15 ms<sup>-1</sup></td><td>R<sup>2</sup>=0.79, CCC=0.89</td></tr><tr><td>CWS</td><td>R<sup>2</sup>=0.40, CCC=0.66</td></tr><tr><td>Fundamental frequency</td><td>20</td><td>R<sup>2</sup>=0.38, CCC=0.65</td></tr></tbody></table>
model, the wind regime at the position of each tree must be provided
and we assume these wind regimes are identical. In the high-time re-
solution validation this led to some gusts hitting the trees which were
not measured by the anemometer and vice versa.

The static response of the QSMs to an applied wind forcing was
similar to the hourly maxima measured in the field. This is in agreement
with recent work which demonstrated that, for *Pinus sylvestris*, a large
part of the tree response was due to simple static effects (Schindler and
Mohr, 2018). Within the range of wind speeds measured in this study
the correlation between model and field data exceeded R² = 0.79,
CCC = 0.89. As the wind speed increased this correlation decreased,
possibly because the model accounts for non-linear effects due to the
displacement of the crown while the extrapolation from field data did
not. The critical wind speeds were systematically overestimated by the
simulation. It is important to note that there are likely to be further non-
linear effects at high-wind speeds that are not accounted for in either
the field data or the model. These are associated with changes in the
drag factor, streamlining of the tree crown and possible plastic de-
formation of the wood. However, we have no data with which to
quantify these effects and assume they are small in comparison with the
main effect of increasing wind speed.

The fundamental frequency is an important parameter since a lower
fundamental frequency means that energy is more easily passed from
the wind to the tree and the tree is at a higher risk of damage (Schütz
et al., 2006). The model predicted the fundamental frequency of trees
with R² = 0.38, CCC = 0.65, but displayed a clear negative offset. This
offset is likely due to either our simplistic treatment of material prop-
erties or the simplification of the QSM. We do not possess independent
validations of the accuracy of the QSM reconstructions and therefore
cannot separate errors due the QSMs and those due to the treatment of
material properties. A review of the material properties of trees
worldwide (Niklas and Spatz, 2010) found a clear relationship between
wood density and elasticity at the species level, but the potential for
wide variation within species, especially in those selected for rapid
growth, and within individual trees remains. Our sensitivity analysis
shows that wood elasticity plays a larger role than wood density in
driving tree dynamics, but it is also harder to measure (Ruel et al.,
2010).

Previous authors increased the accuracy of the finite element ana-
lysis by adjusting the material properties in order to fit the model to
static field tests performed on the tree (Sellier et al., 2006). This method
could prove highly useful for small numbers of urban or valuable ve-
teran trees where proposed treatments, such as removing a large branch
or adding supports, can be tested in a model environment before being
carried out. Another possible addition to the model would be a more
accurate representation of leaves and twigs, which must influence the
motion of the tree. This is currently limited by the difficulty of re-
trieving small twigs and leaves from TLS data, but a parameterization
could possibly be developed using summer (leaf-on) data in combina-
tion with the winter (leaf-off) data.

5. Conclusion

Previously, finite element analysis of trees has been confined to a
small number of relatively short or even idealized trees. We demon-
strate that combining TLS data and finite element analysis allows us to
model the detailed dynamics of realistic and complex broadleaf trees
without manually mapping their 3D architecture. We show that this
approach can model the strains produced in a tree for a given time-
series wind forcing. Finite element analysis can also predict the effects
of a static wind forcing, which was found to follow a similar trend to the
hourly maximum strain - wind speed relationship. It can also predict the
correct trend of fundamental frequencies for a set of trees given their
geometry, although there is a systematic offset.

Several uncertainties remain, in particular relating to the accuracy
of the TLS data in representing 3D tree architecture, the mechanical

properties of roots, and the change in wind-tree interaction at high
wind speeds. This approach is applicable in natural forest plots and
would also be valuable in predicting the dynamics of park or city trees,
and how they will respond to a proposed intervention.

Data availability

All the field data collected for this study is available online (https://
doi.org/10.5285/533d87d3-48c1-4c6e-9f2f-fda273ab45bc)as are the
data processing scripts (https://github.com/TobyDJackson/
WindAndTrees_Wytham). The matlab scripts that convert QSMs to
Abaqus input files are also available online (https://doi.org/http://
doi.org/10.5281/zenodo.894543). We used the commercial finite ele-
ment analysis software Abaqus (Simulia Software Company, 2017), but
a student version can be freely obtained on which a more limited
analysis can be run – this is clearly explained in the software readme
files.

Acknowledgements

We would like to thank the Oxford Department of Engineering,
especially Mike Stockford, for the use of their workshop and access to
Abaqus. Also, Clint Davies-Taylor for his expert advice on Abaqus
scripting. This research was supported by a NERC studentship to TJ on
the DTP Environmental Research (NE/L0021612/1). The TLS data
collection was part-funded through the Metrology for Earth
Observation and Climate project (MetEOC-2), grant number ENV55
within the European Metrology Research Programme (EMRP). The
EMRP is jointly funded by the EMRP participating countries within
EURAMET and the European Union. MD was supported in part by NERC
Standard GrantsNE/N00373X/1 and NE/P011780/1, and via the NERC
National Centre for Earth Observation for capital purchase of lidar
equipment. AS and YM are supported by NERC GrantNE/P012337/1 (to
YM). YM is also supported by the Frank Jackson Foundation.

Part of this study has been funded in the frame of the 3DForMod
project, ERA-NET FACCE ERA-GAS (ANR-17-EGAS-0002-01). FACCE
ERA-GAS has received funding from the European Union's Horizon
2020 research and innovation programme under grant agreement No
696356.

**Appendix A. Supplementary data**

Supplementary material related to this article can be found, in the
online version, at doi:https://doi.org/10.1016/j.agrformet.2018.11.
014.

References

Åkerblom, M., 2017. Inversetampere/Treeqsm: Initial Release. https://doi.org/10.5281/zenodo.844626.

Burt, A.P., 2017. New 3D Measurements of Forest Structure. University College London.

Calders, K., Newnham, G., Burt, A., 2015. Nondestructive estimates of above-ground biomass using terrestrial laser scanning. Methods Ecol. Evol. 6, 198–208.

Calders, K., Burt, A., Origo, N., Disney, M., Nightingale, J., Raumonen, P., Lewis, P., 2016. Large-area virtual forests from terrestrial laser scanning data. in: International Geoscience and Remote Sensing Symposium (IGARSS) 1765–1767. https://doi.org/10.1109/IGARSS.2016.7729452.

Calders, K., Origo, N., Burt, A., Disney, M., Nightingale, J., Raumonen, P., Åkerblom, M., Malhi, Y., Lewis, P., 2018. Realistic forest stand reconstruction from terrestrial LiDAR for radiative transfer modelling. Remote Sens. (Basel) 10, 933. https://doi.org/10.3390/rs10060933.

Carrera, E., Giunta, G., Petrolo, M., 2011. The Euler-Bernoulli and timoshenko theories. Beam Structures. John Wiley & Sons, Ltd, pp. 9–22. https://doi.org/10.1002/9781119978565.ch2.

Ciftci, C., Brena, S., Kane, B., Arwade, S., 2012. Effects of crown architecture on dynamic behavior of trees: understanding the relationship between branches and stem. 2012. Int. Soc. Arboric. Annu. Conf. Trade Show 2012.

Ciftci, C., Brena, S.F., Kane, B., Arwade, S.R., 2013. The effect of crown architecture on dynamic amplification factor of an open-grown sugar maple (Acer saccharum L.). Trees - Struct. Funct. 27, 1175–1189. https://doi.org/10.1007/s00468-013-0867-z.

de Langre, E., 2008. Effects of wind on plants. Annu. Rev. Fluid Mech. 40, 141–168.
https://doi.org/10.1146/annurev.fluid.40.111406.102135.

Dupont, S., Pivato, D., Brunet, Y., 2015. Wind damage propagation in forests. Agric. For. Meteorol 214-215, 243-251. https://doi.org/10.1016/j.agrformet.2015.07.010.

Espírito-Santo, F.D.B., Gloor, M., Keller, M., Malhi, Y., Saatchi, S., Nelson, B., Junior, R.C.O., Pereira, C., Lloyd, J., Frolking, S., Palace, M., Shimabukuro, Y.E., Duarte, V., Mendoza, A.M., López-González, G., Baker, T.R., Feldpausch, T.R., Brienen, R.J.W., Asner, G.P., Boyd, D.S., Phillips, O.L., 2014. Size and frequency of natural forest disturbances and the Amazon forest carbon balance. Nat. Commun. 5. https://doi.org/10.1038/ncomms4434.

Finnigan, J., 2000. Turbulence in plant canopies. Annu. Rev. Fluid Mech. 32, 519-571. https://doi.org/10.1146/annurev.fluid.32.1.519.

Gonzalez de J. Tanago, Lau, A., Bartholomeus, H., Herold, M., Avitabile, V., Raumonen, P., Martius, C., Goodman, R., Disney, M., Manuri, S., Burt, A., Calders, K., 2018. Estimation of above-ground biomass of large tropical trees with Terrestrial LiDAR. Methods Ecol. Evol. 9, 223-234. https://doi.org/10.1111/2041-210X.12904.

Hale, S.E., Gardiner, B.A., Wellpott, A., Nicoll, B.C., Achim, A., 2012. Wind loading of trees: influence of tree size and competition. Eur. J. For. Res. 131, 203-217. https://doi.org/10.1007/s10342-010-0448-2.

Hale, S.A., Gardiner, B., Peace, A., Nicoll, B., Taylor, P., Pizzirani, S., 2015. Comparison and validation of three versions of a forest wind risk model. Environ. Model. Softw. 68, 27-41. https://doi.org/10.1016/j.envsoft.2015.01.016.

James, K.R., 2010. A dynamic structural analysis of trees subject to wind loading. Melb. Sch. L. Environ. Doctor 277.

James, K.R., Haritos, N., Ades, P.K., 2006. Mechanical stability of trees under dynamics loads. Am. J. Bot. 93 (10), 1522-1530. https://doi.org/10.3732/ajb.93.10.1501.

James, K.R., Dahle, G.A., Grabosky, J., Kane, B., Detter, A., 2014. Tree biomechanics literature review: dynamics. Arboric. Urban For. 40 (1), 1-15. https://doi.org/10.13140/RG.2.1.1089.1765.

Kane, B., Modarres-Sadeghi, Y., James, K.R., Reiland, M., 2014. Effects of crown structure on the sway characteristics of large decurrent trees. Trees - Struct. Funct. 28, 151-159. https://doi.org/10.1007/s00468-013-0938-1.

Kerzenmacher, T., Gardiner, B., 1998. A mathematical model to describe the dynamic response of a spruce tree to the wind. Trees - Struct. Funct. 12, 385-394. https://doi.org/10.1007/s004680050165.

Krauss, T., Shure, L., Little, J., 1994. Signal Processing Toolbox for Use With MATLAB®: User's Guide.

Lavers, G.M., 1983. The strength properties of timbers. For. Prod. Res. Lab. Bull. 60, 50.

Moore, J.R., Maguire, D.A., 2008. Simulating the dynamic behavior of Douglas-fir trees under applied loads by the finite element method. Tree Physiol. 28, 75-83. https://doi.org/10.1093/treephys/28.1.75.

Moore, J.R., Gardiner, B.A., Blackburn, G.R.A., Brickman, A., Maguire, D.A., 2005. An inexpensive instrument to measure the dynamic response of standing trees to wind loading. Agric. For. Meteorol. 132, 78-83. https://doi.org/10.1016/j.agrformet.2005.07.007.

Niklas, K.J., Spatz, H.C., 2010. Worldwide correlations of mechanical properties and green wood density. Am. J. Bot. 97, 1587-1594. https://doi.org/10.3732/ajb.1000150.

Raumonen, P., Casella, E., Calders, K., Murphy, S., Åkerbloma, M., Kaasalainen, M., 2015. Massive-scale tree modelling from TLS data. ISPRS Ann. Photogramm. Remote Sens. Spat. Inf. Sci. II-3/W4, 189-196. https://doi.org/10.5194/isprsannals-II-3-W4-189-2015.

Rayleigh, R., 1877. On progressive waves. Proc. London Math. Soc. s1-9, 21-26. https://doi.org/10.1112/plms/s1-9.1.21.

Reiland, M., Kane, B., Modarres-Sadeghi, Y., Ryan, H.D.P., 2015. The effect of cables and leaves on the dynamic properties of red oak (Quercus rubra) with co-dominant stems. Urban For. Urban Green. 14, 844-850. https://doi.org/10.1016/j.ufug.2015.08.010.

Rodriguez, M., De Langre, E., Moulia, B., 2008. A scaling law for the effects of architecture and allometry on tree vibration modes suggests a biological tuning to modal compartmentalization. Am. J. Bot. 95, 1523-1537. https://doi.org/10.3732/ajb.0800161.

Rodriguez, M., Ploquin, S., Moulia, B., de Langre, E., 2012. The multimodal dynamics of a walnut tree: experiments and models. J. Appl. Mech. 79, 44505. https://doi.org/10.1115/1.4005553.

Rudnicki, M., Silins, U., Lieffers, V.J., Josi, G., 2001. Measure of simultaneous tree sways and estimation of crown interactions among a group of trees. Trees 15, 83-90. https://doi.org/10.1007/s004680000080.

Rudnicki, M., Lieffers, V.J., Silins, U., 2003. Stand structure governs the crown collisions of lodgepole pine. Can. J. For. Res. 33, 1238-1244. https://doi.org/10.1139/x03-055.

Ruel, J.C., Achim, A., Herrera, R.E., Cloutier, A., 2010. Relating mechanical strength at the stem level to values obtained from defect-free wood samples. Trees - Struct. Funct. 24, 1127-1135. https://doi.org/10.1007/s00468-010-0485-y.

Schindler, D., Mohr, M., 2018. Non-oscillatory response to wind loading dominates movement of Scots pine trees. Agric. For. Meteorol 250-251, 209-216. https://doi.org/10.1016/j.agrformet.2017.12.258.

Schütz, J.P., Götz, M., Schmid, W., Mandallaz, D., 2006. Vulnerability of spruce (Picea abies) and beech (Fagus sylvatica) forest stands to storms and consequences for silviculture. Eur. J. For. Res. 125 (3), 291-302. https://doi.org/10.1007/s10342-006-0111-0.

Scurfield, G., 1973. Reaction wood: its structure and function. Science 80 (179), 647-655.

Sellier, D., Fourcaud, T., 2009. Crown structure and wood properties: influence on tree sway and response to high winds. Am. J. Bot. 96, 885-896. https://doi.org/10.3732/ajb.0800226.

Sellier, D., Fourcaud, T., Lac, P., 2006. A finite element model for investigating effects of aerial architecture on tree oscillations. Tree Physiol. 26, 799-806. https://doi.org/10.1093/treephys/26.6.799.

Sellier, D., Brunet, Y., Fourcaud, T., 2008. A numerical model of tree aerodynamic response to a turbulent airflow. Forestry 81, 279-297. https://doi.org/10.1093/forestry/cpn024.

Simulia Software Company, 2017. Abaqus Theory Guide Abacus 6.4. 2014. URL. (Accessed 7.7.17). http://50.16.225.63/v6.14/books/stm/default.htm.

Spatz, H., Theckes, B., 2013. Plant science oscillation damping in trees. Plant Sci. 207, 66-71. https://doi.org/10.1016/j.plantsci.2013.02.015.

Spatz, H.C., Brüchert, F., Pfisterer, J., 2007. Multiple resonance damping or how do trees escape dangerously large oscillations? Am. J. Bot. 94, 1603-1611. https://doi.org/10.3732/ajb.94.10.1603.

Viberg, M., 1995. Subspace-based methods for the identification of linear time-invariant systems. Automatica (Oxf) 31, 1835-1851. https://doi.org/10.1016/0005-1098(95)00107-5.

Vicari, M.B., 2017. Leaf Separation in TLS Data. https://doi.org/10.5281/ZENODO.802233.

Wellpott, a., 2008. The stability of continuous cover forests. Sch. Geosci. PhD 160.

Wilkes, P., Lau, A., Disney, M., Calders, K., Burt, A., Gonzalez de Tanago, J., Bartholomeus, H., Brede, B., Herold, M., 2017. Data acquisition considerations for Terrestrial Laser Scanning of forest plots. Remote Sens. Environ. 196, 140-153. https://doi.org/10.1016/j.rse.2017.04.030.