Stronger wind, smaller tree: Testing tree growth plasticity through a modeling approach

Plants exhibit plasticity in response to various external conditions, characterized by changes in physiological and morphological features. Although being non-negligible, compared to the other environmental factors, the effect of wind on plant growth is less extensively studied, either experimentally or computationally. This study aims to propose a modeling approach that can simulate the impact of wind on plant growth, which brings a biomechanical feedback to growth and biomass distribution into a functional–structural plant model (FSPM). Tree reaction to the wind is simulated based on the hypothesis that plants tend to fit in the environment best. This is interpreted as an optimization problem of finding the best growth-regulation sink parameter giving the maximal plant fitness (usually seed weight, but expressed as plant biomass and size). To test this hypothesis in silico, a functional–structural plant model, which simulates both the primary and secondary growth of stems, is coupled with a biomechanical model which computes forces, moments of forces, and breakage location in stems caused by both wind and self-weight increment during plant growth. The Non-dominated Sorting Genetic Algorithm II (NSGA-II) is adopted to maximize the multi-objective function (stem biomass and tree height) by determining the key parameter value controlling the biomass allocation to the secondary growth. The digital trees show considerable phenotypic plasticity under different wind speeds, whose behavior, as an emergent property, is in accordance with experimental results from works of literature: the height and leaf area of individual trees decreased with wind speed, and the diameter at the breast height (DBH) increased at low-speed wind but declined at higher-speed wind. Stronger wind results in a smaller tree. Such response of trees to the wind is realistically simulated, giving a deeper understanding of tree behavior. The result shows that the challenging task of modeling plant plasticity may be solved by optimizing the plant fitness function. Adding a biomechanical model enriches FSPMs and opens a wider application of plant models.

1 Introduction

Environmental conditions, such as light, temperature, and humidity, influence the physiological processes and structural development of trees. Wind, being almost ubiquitous in nature, affects tree growth, reproduction, and even survival (Moulia et al., 2006; de Langre, 2008; Lopez et al., 2011). The “compensating mechanism” whereby advantageous modifications of a tree’s phenotype can emerge to adapt to the wind conditions has already been discovered (Whitehead, 1962; Whitehead & Luti, 1962). The term thigmomorphogenesis has been coined, which refers to a phenomenon that trees demonstrate plasticity in reaction to the external mechanical stimulus (Jaffe, 1973; Jaffe & Forbes, 1993; Telewski & Pruyn, 1998). Considerable support to this hypothesis has been provided by experimental observations: when trees were exposed to wind forces or other biomechanical perturbations, total dry weight, in particular stem dry weight, decreased, and biomass allocation among organs altered (Bonnesoeur et al., 2016). A reduction in both stem height and leaf area, but an increase in diameter (Whitehead & Luti, 1962; Jaffe & Forbes, 1993; de Langre, 2008), has been observed. Allometry of axes, e.g., the ratio between length and radius, is also affected. Wind sways induce strain force on trees, and it is demonstrated that plants can sense strains they are subjected to (Coutand & Moulia, 2000; Coutand et al., 2014). In some circumstances, the secondary growth (increase in the stem diameter) is more sensitive to changes in environmental conditions than the primary growth (increase in the stem length) (Collet et al., 2001).

How does a plant adjust its behavior accordingly? A plant simulation system can be a helpful way to support or test biological hypotheses (Trewavas, 2017; Calvo et al., 2019). The functional-structure plant model (FSPM) is a key tool for biomechanical computation as it provides detailed information on plant structure and takes into account the interaction between plant architecture and physiological processes during plant growth (Fourcaud et al., 2008). Typically, FSPMs include two components: a structural model and a functional model. The one-directional or bidirectional interaction between these two sub-models is necessarily a feature of FSPMs (Perttunen et al., 1996; de Reffye et al., 2020). The bending of a tree stem has been under study when a very early coffee tree structure was simulated (de Reffye, 1976). The mechanics of orientation of a plant stem due to the changes in loading and stresses have been studied with biomechanical models (Fournier et al., 1994; Almeras et al., 2002; Taylor-Hell et al., 2005). Considering the dynamic growth of trees, an incremental biomechanical model has been introduced to calculate the deformation of a single growing stem (Fourcaud et al., 2003; Fourcaud & Lac, 2003), simulating the influence of wind and self-weight. The modeling work has been done at the scale of the stem (Coutand et al., 2011; Guillon et al., 2012) or even the whole plant (Fourcaud et al., 2003). However, the simulation of tree thigmomorphogenesis under a wind environment is hardly seen. The feedback between mechanics and plant growth has rarely been taken into account in growth models at the architectural scale (Eloy et al., 2017).

This study aims to explore tree growth plasticity under a wind environment in a mechanistic way through a modeling approach. To achieve this goal, firstly, a FSPM GreenLab (de Reffye et al., 2020) is coupled with a biomechanical model (Ancelin et al., 2004b) to simulate the effect of wind on biomechanical variables inside a full tree structure; secondly, a sensitivity analysis is performed to evaluate the effect of key parameters controlling the secondary growth; thirdly, tree plasticity is simulated through parameter optimization for different wind speeds searching the best parameter values that benefit tree growth. Thigmomorphogenesis is modeled as an emergent result of optimization in biomass reallocation.

The validation of this theoretical result can be made by experiments with different wind environments since parameters for secondary growth can be estimated inversely from measured data (Letort et al., 2008). However, as the sampling in such experiments is very costly, an in silico experiment remains an appealing tool to test hypotheses and inspire new thought. In a former study, the optimization of wood yield with the constraint of stability has been conducted (Qi et al., 2009), which also attempts to find the best parameter for secondary growth. The current work differs from it in several aspects: first, the wind environment is considered, while the former is limited to gravity. Secondly, this work couples the biomechanical and growth models, so that the incremental growth of trees is considered, while the former evaluates the mechanical properties of a straight adult tree. Thirdly, the biomechanical properties for all stems, including branches, are computed, while the previous work was limited to the main stem. The difference also lies in the constraint of the biomechanical stability of the tree.

2 Material and methods

2.1 System overview

The components of the whole system are shown in Figure 1. The plant model simulates step by step the tree structure and growth driven by a source–sink balance, giving the topological and geometrical tree structure with detailed information of the size and weight of individual organs (leaves, internodes, etc.). The wind environment model provides the wind profile and corresponding drag force. As this work is dedicated to wind, other environmental factors like light or temperature are not considered here, regarded as non-limiting. The biomechanical model computes the average force, moment, and new position of each internode under wind and gravity, for each growth step. The Factor of Safety Model checks whether the critical point of breakage is reached for every stem. The optimization model simulates plant plasticity by searching parameter values to maximize the objective function of tree growth while keeping the main stem unbroken.

FIGURE 1

Figure 1 Schematic diagram of the modeling system.

For a plant growth model, it is necessary to clarify at which level of detail the model is built. Here, the temporal scale corresponds to a growth cycle (GC), i.e., the time it takes to create a new growth unit (GU, Figure 2A), which is a year for temperate trees. The spatial scale corresponds to a phytomer, composed of an internode and its axillary leaves. A GU is composed of one or more phytomers of the same feature.

FIGURE 2

Figure 2 Illustration of tree development, primary and secondary growth, as well as the simulation of wind and gravity effect during a growth cycle. (A) plant top with two physiological ages (shown with different textures) at GC n - 1; (B) organogenesis, primary growth (blank rectangles) and secondary growth (new rings on existing rectangles) at GC n; (C) permanent deformation caused by increment in gravity; (D) temporary deformation caused by wind.

The integration of a dynamic plant growth process for plants and biomechanical model is as follows. At each GC, the organogenesis and organ growth take place, which creates new organs and increases the existing stem diameter (Figure 2B). A biomechanical model is called to compute force and moment as well as the new stem position caused by the increment in self-weight (Figure 2C). The deformation caused by self-weight is regarded as permanent. Given wind load, the wind-induced force leads to new destabilizations, and the biomechanical model is called again until a new equilibrium status, as well as the displacement of each internode, is reached (Figure 2D). However, after the withdrawal of wind, the plant recovers its original position as one can see in nature. This is achieved by running again the biomechanical module after the removal of wind force. This step is necessary as branch breakage induced by wind may have taken place, which influences self-weight distribution. In the next growth cycle, the plant development continues based on the deformed or broken structure of the last cycle.

2.2 Biomechanical model

Stem deformation is calculated at each internode k using the Incremental Transfer Matrix Method (ITMM) (Ancelin et al., 2004b), as in Eq. 1.

$\begin{array}{ll}\{\begin{array}{l}{D}_{k+1}=T_{11}{D}_k+T_{12}(C_k-S_k)+M_{1}F+M_a\\ S_{k+1}=T_{21}{D}_k+T_{22}(C_k+S_k)+M_{2}F\end{array}& (1)\end{array}$

where D is the displacement vector, which includes the positive translations of the cross-section center (m) and the positive rotations of the cross-section normal (radians) (see components in Supplement A), S is the stress and moment vector; C is the vector of the sum of all concentrated forces and moments at the branching node, with the incremental weight of the branch computed by the plant model; F is the vector of linear forces, caused by both the incremental self-weight and wind load (Eq. 9); and T and M are transformation matrices (see Eqs. S3-S6 and Eqs. S7-S12 in Supplement A). M_a denotes the maturation strains (MS) at the periphery of stems, associated with the formation of reaction wood (see Eq. S13 in Supplement A). The link between the biomechanical model and plant model lies in the concentrated forces and moments at the branching node (C), and the self-weight of internode (F), as the FSPM gives the 3D structure and the concentrated force of each sub-structure.

During the dynamic growth of plants, in each growth cycle, tree secondary growth takes place by adding a new layer to the existing stems. In the biomechanical model, a tree stem is thus regarded as a sequence of multilayer three-dimensional beam elements (Ancelin et al., 2004b). The biomechanical calculation consists of two steps. Firstly, compute the new incremental loads and the moments (caused by self-weight) on each internode. This is done using the idea of sub-structure (Kang et al., 2008), starting from the twigs of the highest branching order until the trunk (main stem). In ITMM, only the increment of weight must be considered in a growth step. Leaves, fruits, and flowers are considered as concentrated forces and moments for each internode. Dead or pruned branches are removed from the insertion point on the parent branch. Secondly, compute the new incremental displacement based on the last equilibrium, until the top of the main stem. Once this process is completed, the result is validated by examining the biomechanical boundary condition (BC) at the free top of the tree (Ancelin et al., 2004a), so that there is hardly any force or moment at each tip. If the BCs are not met, the last two steps are repeated.

When the wind is imposed, the displacement of the internode will take place, which changes the gravity center. With this increment of load, at each growth cycle, ITMM needs to be run multiple times to reach a new equilibrium, as if the plant swings before becoming stable. Similarly, when the wind calms down, iterated computation is called again. As a result, at each growth cycle, the equilibrium status of the tree with and without wind is updated. The iteration number can be set respectively for the main stem and branches.

2.3 Incremental weight in plant structure

The incremental weight of individual organs during plant growth is computed using the GreenLab model. GreenLab is a generic FSPM that has been applied for beech tree (Fagus sylvatica L.) (Letort et al., 2008), poplar tree (Salicaceae) (Yang et al., 2011), and Chinese pine (Pinus tabulaeformis Carr.) (Guo et al., 2012). For trees, a model with a finer temporal scale with stochastic behavior has been proposed (Kang et al., 2018), but here a deterministic model based on the growth cycle is adopted, where internodes created during the same GU have the same size. We recall here mainly components for secondary growth. See Supplement B for a more complete presentation, especially the primary growth.

For biomechanical computation, the output needed from GreenLab include the length and diameter of individual internodes; the weight, position, and gravity center of the aggregated branch; and the initial geometrical shape of the tree structure for computing the gravity center. The last is computed from the topological structure, the size of organs, and the insertion and phyllotaxy angle, which is not detailed here. The poplar (Populus) tree is selected for this study. Parameter settings used in this study are shown in Table 1.

TABLE 1

Table 1 Parameter values of poplar tree.

2.3.1 Weight of organs and substructure

The plant growth is based on source–sink regulation. Based on the pipe model, the demand for the secondary growth ( $D_{\text{sec}}^{\text{layer}}$) is proportional to the number of leaves in the last cycle $N_p^b(n-1)$ , using Eq. 2:

$\begin{array}{ll}{D}_{\text{sec}}^{\text{layer}}(n)=S_{\text{layer}}{\displaystyle \sum _{p=1}^{P_m}{N}_p^b}(n-1)& (2)\end{array}$

where S_layer denotes the sink strength of a new layer contributed by an individual leaf. It controls the proportion of produced biomass allocated to the secondary growth at each GC. Subscript p is the physiological age (PA), from 1 to maximum physiological age P_m. Superscript b denotes the leaf blade. Symbol n denotes growth cycle n. In GreenLab, new rings of wood formed around the stem periphery. Each ring on a stem is contributed by numerous layers. Rings are modeled in two steps. Firstly, calculate the amount of biomass allocated to the secondary growth at GC n at the whole-plant scale, denoted by $Q_{\text{sec}}^{\text{layer}}(n)$, as in Eq. 3:

$\begin{array}{ll}{Q}_{\text{sec}}^{\text{layer}}(n)=D_{\text{sec}}^{\text{layer}}(n)\frac{Q(n)}{D(n)}& (3)\end{array}$

Secondly, $Q_{\text{sec}}^{\text{layer}}(n)$ is allocated to each internode with two possible modes. Let parameter λ describe the way of biomass allocation associated with the secondary growth for each phytomer; the incremental biomass at GC n for an internode of PA p that showed GC j is described in Eq. 4:

$\begin{array}{ll}{q}_{\text{sec},p}^{\text{layer}}(j,n)=\left(\frac{1-\lambda }{D_1(n)}+\frac{\lambda ·s_p^a\left(j,n\right)}{D_2\left(n\right)}\right)·l_p\left(j\right)·P_p^{\text{r}}·Q_{\text{sec}}^{\text{layer}}(n)& (4)\end{array}$

See Eq. S22 for demand for primary growth D₁(n) and secondary growth D₂(n), as well as other parameters. If λ = 0, it means that at each cycle the biomass for the layer is uniformly allocated to internodes of the same PA; tapering can appear because older internodes have more layers. If λ = 1, the allocation is also proportional to the total number of living leaves above the internode; since upper internodes have fewer leaves above, thus tapering is more obvious compared to the case of λ = 0. The weight of an individual internode that is born in cycle j is thus the result of primary and secondary growth:

$\begin{array}{ll}{q}_p^i(j,n)=q_{\text{pri},p}^i(j)+{\displaystyle \sum _{i=j+1}^n{q}_{\sec ,p}^{layer}}(j,i)& (5)\end{array}$

Two key parameters that control the biomass allocation to the secondary growth are (1) parameter S_layer controlling the amount of biomass allocated to the secondary growth at the whole-plant level (Eq. 2) and (2) parameter λ controlling the allocation of biomass for secondary growth among all internodes (Eq. 4) and, accordingly, stem shape, as shown in Supplementary Figure S1.

2.3.2 Weight of plant and substructures

The weight of a branch (sub-structure) of PA p is the total weight of its bearing axis (stem) and all axillary axes (Kang et al., 2008). With the sub-structure method, the summing up starts from the branch of the highest PA (no sub-branch) until the main stem of PA 1. Similarly, the incremental weight of the branch can be obtained.

The total stem weight of the structure of the plant is the sum of the weights of the main stem and all branch stems, which is one of the object functions for optimization:

$\begin{array}{ll}{Q}_{\text{stem}}={\displaystyle \sum _{p=1}^{P_m}}{\displaystyle \sum _{j=1}^n{N}_p^i}(j)q_p^{\text{i}}(j,n)& (6)\end{array}$

The number of internodes that are created at cycle j, , is updated in each cycle considering both organogenesis and possible stem breakage.

2.4 Size of organs

Considering that wood density is dependent on the radial position and height of the internode, parameter ${\rho }_p^{\text{layer}}(j,i)$ is set to describe the density of the layer in GC i for internode of PA p created at GC j. Then, at cycle i, the section area of layer added to such internodes is computed as

$\begin{array}{ll}{s}_{\text{sec},p}^{\text{layer}}(j,i)=\frac{q_{\text{sec},p}^{\text{layer}}(j,i)}{{\rho }_p^{\text{layer}}(j,i)*l_p(j)}& (7)\end{array}$

The radius of internode is finally the result of layer accumulation:

$\begin{array}{ll}{r}_p^{\text{i}}(j,n)=\sqrt{s_p^{\text{i}}(j,n)/\pi }=\sqrt{[s_p^{\text{pith}}(j)+{\displaystyle \sum _{\text{i}=\text{j}+1}^{\text{n}}{s}_{\text{sec},p}^{\text{layer}}}(j,i)]/\pi }& (8)\end{array}$

where l_p(j) and $s_p^{\text{pith}}(j)$ are the length and section area of pith from primary growth (see Supplement B).

2.5Wind-induced force

According to literature (Mayhead, 1973; Niklas & Spatz, 2000; Ancelin et al., 2004a; Diener et al., 2010), the wind-induced drag force Fw by wind on the kth element at the height z aboveground can be described as in Eq. 9:

$\begin{array}{ll}\begin{array}{l}F{\text{w}}_k(z)=\frac{1}{2}{\rho }_{air}{C}_d{A}_{k}u(z{)}^2\\ A_k=2l_k(z)r_k(z)\left|\sin \chi \right|\end{array}& (9)\end{array}$

where ρ_air is the air density (1.226 kg·m^-3 when the air temperature is 15°C). C_d (dimensionless) is the drag coefficient, assumed to be 0.25 according to Koizumi et al. (2010). Differing from Sellier et al. (2008) where the exposed area is estimated from the biomass ratio, here A_k (m²) is the frontal area of the kth internode, i.e., the projected area of the internode against the wind direction. It is computed with the length, diameter of the internode, and angle to the vertical direction of the internode, where l and r are the length and radius of the internode, respectively; X is the angle between this internode with the wind direction. u(z)is the wind speed at vertical height z which is given by a logarithmic wind profile.

Wind profile refers to the vertical distribution of wind speed. There are different wind profiles in the literature describing the horizontal mean wind speed according to the height aboveground (Niklas & Spatz, 2000). Some researchers (Irvine et al., 1997; Ancelin et al., 2004a) distinguished the wind profiles in forests located inside the stand and at the stand edge, including above and below the canopy. For trees standing in an open terrain, as assumed in this paper, a logarithmic wind profile is used (Ancelin et al., 2004b), which has provided the best fit for wind speeds empirically recorded within and just above the canopy (Niklas & Spatz, 2000).

$\begin{array}{ll}u(z)=\text{u}(h_0)\frac{\text{ln}(z/z_o)}{\text{ln}(h_0/z_o)}& (10)\end{array}$

where u(z) (m·s^-1) is the wind speed at vertical height z; h₀ is the standard reference height of 10 m in Meteorology (Davenport, 1962); z₀ is the roughness length at stand edge (m). Ratio z_o/h₀ is set to 0.06 (Peltola & Kellomäki, 1993; Ancelin et al., 2004a).

2.6 Factor of safety

It is assumed that stem breakage takes place when its maximum longitudinal stress exceeds the Modulus of Rupture (MOR) of wood fibers (Jones, 1983; Gardiner et al., 2000). According to Jones (1983); Gardiner et al. (2000); Grotta et al. (2005), and Fjeld (2012), it is found that individual trees differed significantly in their MOR (range: 37–148 MN/m²). The maximum absolute value of longitudinal stress is found at the periphery of stem cross sections and depends on the bending moment (M, N·m ) of the stem and the third power of the stem diameter (Jones, 1983; Gardiner et al., 2000), as in Eq. 11:

$\begin{array}{ll}{\text{Stress}}_k=\frac{32M_k}{\text{π}{d}_k^3}& (11)\end{array}$

Differing from Gardiner et al. (2000) who assumed that the stress is constant and thus calculated the stress only at breast height, here for each internode this stress is computed. The bending moment M_k caused by the wind results from Eq. 1 (a component of S_k). d_k is the diameter of the kth stem internode as calculated in Eq. 8, and Stress_k (Pa) is the longitudinal stress of internode k.

The critical wind speed (CWS) is defined as the minimum wind speed causing the breakage of the trunk, with the given parameter set. Higher CWS means that the tree is more wind-resistant.

The wind resistance can also be indicated by the top stem deflection of the trunk, defined as the angle between the trunk’s top and vertical direction.

2.7 Objective function

As mentioned above, the adaptation of trees to varying circumstances is modeled as an optimization problem. The optimization objective of the acclimatization is set to be the wood fresh weight Q_stem (total internode weight, Eq. 6) and the height of tree H_tree. The reason for adding tree height as one of the objectives is that reaching a certain height is beneficial for light interception. Since in stronger winds the tree is more prone to breakage, thus it is expected that the tree’s reaction to the wind is to increase the radial growth. However, over-allocation to secondary growth will diminish the photosynthesis production and tree growth; thus, there is a balance between primary and secondary growth. The aim of optimization is to find the best balance. We choose the sink strength parameter for the secondary growth S_layer (Eq. 2) as the control variable.

Stem breakage under loads includes both trunk breakage and branch breakage (Lopez et al., 2011). Breakage on the trunk is likely to be fatal, while tree survival is often possible after branch breakage. Thus, here trunk breakage is considered a constraint condition.

In summary, this optimization problem with constraints is given by Eq. 12:

$\begin{array}{ll}\{\begin{array}{c}Max H_{\text{tree}}(S_{\text{layer}})\\ Max Q_{\text{stem}}(S_{\text{layer}})\\ subject to Stres{s}_{\text{trunk}}<MOR\end{array}& (12)\end{array}$

The value of MOR is 45 MPa (Table 1). The multi-objective optimization algorithm of the Non-determined Sorting Genetic Algorithm II (NSGA-II) (Deb et al., 2000) from PAGMO (a C + + scientific library for massively parallel optimization) (Biscani & Izzo, 2020) is adopted. NSGA-II is one of the most popular multi-objective genetic algorithms. It has the advantages of fast-running speed and good convergence of the solution set, making it a benchmark of the performance of other multi-objective optimization algorithms. The population size is set as 100, the evolutionary algebra as 100, the crossover probability as 0.7, and the mutation probability as 0.4.

3 Results

3.1 The equilibrium of tree under wind load

As mentioned above, the biomechanical calculation (Eq. 1) consists of two steps and requires multiple computational iterations to reach equilibrium. Under wind, the stem deformation stabilizes after several iterations or simply breaks under strong wind (Figure 3). At a wind speed of 13 m/s, the tree sways during about 10 iterations before the stabilization. When the wind speed is 15 m/s, trunk breakage takes place after the first iteration. With more iteration steps, the computational time increases linearly. With 10 iterations, at age 8 of the plant which means 8 years, it takes about 20 s for running the whole model in a PC of Core™ i7-8550U CPU 1.80 GHz, 16 GB RAM. Let N_p be the iteration number for the substructure of PA p. According to the result, in the following computations, the iteration number of ITMM is set to 10 since the displacement error becomes negligible.

FIGURE 3

Figure 3 The stem deflection angle of the trunk under wind speeds of 13 m/s and 15 m/s. When wind speed is 13 m/s, it takes about 10 iterations before the stabilization. When wind speed is 15 m/s, trunk breakage takes place after the first iteration.

3.2 Numerical results

3.2.1 Force and moment

To visualize how the system works, the tree structure of GC 8 is simulated with parameter set S_layer = 4, λ = 0.1, under wind speed 15 m/s. The result is shown for the first iteration before trunk breakage. The displacement, force (Figures 4A-C, the first row), moment (Figures 4D-F, the second row), and stress_k (Figures 4G-I, the third row) of each stem under a wind environment (which blows in the x-direction from left to right) are shown through tree shape and color, for the x-, y-, and z-components, respectively. The force is a vector. When the direction of the force is inconsistent with the direction of the coordinate axis, the force takes a negative value.

FIGURE 4

Figure 4 Computed force (A-C), moment (D-F) and stress_k. (G-I) of x-, y-, z- component for the plant under wind speed 15 m/s, shown with different colors. S_layer = 4, λ = 0.1.

3.2.2 Effect of wind speed

To understand how wind speed can influence tree breakage, the tree shapes are shown in Figure 5. Here, S_layer and λ are set to 6 and 0.1, respectively. With increasing wind speed, more branches are lost because of breakage. The effect is evident on the top of the trunk.

FIGURE 5

Figure 5 The effect of wind speed on tree breakage: 13m/s (A), 15m/s (B), 18m/s (C), respectively. Result is shown under calm conditions. The remaining trunk biomass is 170.78, 159.85 and 120.58 kg, respectively. The circles indicate the points where the trunk breakage takes place. The plant age is 8, S_layer = 6, λ = 0.1.

3.3 Parameter sensitivity

3.3.1 Parameter S_layer

Before parameter optimization, it is useful to understand the model behavior using a sensitivity analysis. According to the GreenLab model, the biomass allocation to the secondary growth, regulated by parameter S_layer (Eqs. 2, 3), has twofold effects. On the one hand, the increase in biomass allocation to secondary growth increases the stem biomass. On the other hand, it can inhibit primary growth, thus decreasing the total leaf area and biomass production. In Figure 6, the stem biomass and top stem deflection for a trunk with different S_layer values are simulated. As expected, the peak value of the biomass of tree and trunk exists, which is independent of the wind. With a lower S_layer value, stem breakage is more evident; a sudden drop for S_layer<8 means trunk breakage. Supplementary Figure S2 shows 3D tree shapes with different S_layer values, where broken branches are removed. With a higher S_layer value, the tree structure is more complete but smaller, because of the smaller leaf area for photosynthesis.

FIGURE 6

Figure 6 Influence of the sink strength for the secondary growth (S_layer) on the wood biomass Q_stem and top step deflection for the trunk, under wind speed 15 m/s with λ = 0.4.

The critical value of wind speed (CWS) gradually increases with S_layer (Figure 7). It is because when more biomass was allocated to the secondary growth, the tree is easier to resist the breakage induced by wind. This shows from another point of view than the step deflection angle that with a bigger S_layer, the tree is more resistant to wind.

FIGURE 7

Figure 7 The critical wind speed (CWS) increases with S_layer values with λ = 0.4.

3.3.2 Parameter λ

Differing from the conclusion that λ should be equal to 1 to obtain the maximal wood biomass (Qi et al., 2009), here the influence of λ is more complex, as shown in Figure 8. At a lower S_layer value, e.g., S_layer = 0.5, the trunk is thin and the diameter of the tree top is too slim to support itself. Even very small wind can cause a breakage due to self-weight. The position with the maximum value of the breakage value also moves downward along the trunk. At a higher S_layer value, e.g., S_layer = 4, the critical wind bore by the trunk becomes stronger and becomes indifferent to the λ value. Table 2 shows the force and moment given MOR = 45 MN/m². Supplementary Figure S3 shows the 3D tree shapes under different λ, where the largest stress value on the trunk Stress_Trunk is marked with a black circle. Only the trunk and the primary branch structure are shown. λ does not change the height of the tree as it matters only the allocation of biomass for secondary growth, but it influences the stem shape and accordingly breakage. According to the numerical results, λ has little effect on optimization results; thus, it is not considered as a control parameter for tree reaction (set to λ = 0.1).

FIGURE 8

Figure 8 Effect of λ on the critical wind speed (CWS) with different S_layer values.

TABLE 2

Table 2 The force and moment given MOR 45 MN/m².

3.4 Plant plasticity simulation

3.4.1 Pareto optimal frontier

Tree plasticity is simulated by searching the best parameter value of S_layer to maximize the tree biomass and tree height while keeping the trunk unbroken, given the wind speed. Because of multi-objective optimization, the result is no more a single point but a point set, called Pareto set or Pareto frontier (Figure 9). In general, when the wind speed is smaller than 16 m/s, the distribution of points in the Pareto boundary coincides with each other, which means the reaction tree is close to the relatively small wind. When the wind speed is high, as seen at the wind speed of 16 m/s, the points in the Pareto optimal front are separated from the previous ones, which means a different reaction. In both situations, according to the Pareto sets, the optimal tree can be a taller tree with smaller biomass, or a shorter tree with bigger biomass. In general, the tree height shows an obvious downward trend with wind speed. The total weight does not change significantly when the wind speed is low but has wider distribution with the increase in wind speed.

FIGURE 9

Figure 9 Pareto optimal front for different wind speeds.

3.4.2 Optimal parameter values

Each dot in the Pareto frontier corresponds to an optimal parameter, so the box diagram is used to show the dispersion degree of the optimal S_layer for each wind speed (Supplementary Figure S4). When the wind speed is small, the box is small, which means the parameters are more concentrated. When the wind speed increases to 16 m/s or bigger, the maximum deviation of the “outlier” from the box gradually increases, indicating that the tree tends to allocate more biomass to the secondary growth to reinforce its structure under the wind. With the increase in wind speed, the median value of each box has a gradual upward trend.

Table 3 shows the optimal parameter values under different wind speeds according to the box diagram (Supplementary Figure S4) by taking representative values (the maximum, minimum, median, and the ratios of 10%, 25%, 75%, and 90%) of the box. In general, the trend is that the parameter value increased with wind speed.

TABLE 3

Table 3 Optimized parameter values for different wind speeds.

3.4.3 Morphological result

To better analyze the morphological changes of trees reacting to different wind conditions, three points (representing three strategies) at the two ends and the middle of the Pareto frontier are selected. For each wind speed, the corresponding S_layer value of the same strategy is shown in Table 4. The strategy I (0% H, 100% Q) in Table 4 represents the point corresponding to the minimal tree height and the maximal biomass in the Pareto optimal front. This point can be considered as the tree focuses on pursuing weight goals. In this case, when the wind speed is small (less than 12 m/s), the change of S_layer is not obvious, so it can be considered that it is not sensitive to the stimulation of the breeze. Strategy III (100% H, 0% Q) corresponds to the point with the largest tree height and the smallest biomass in the Pareto optimal frontier. It can be considered that in this case, trees focus on pursuing tree height but not the tree biomass. In this case, the S_layer value is taken in a small range. Strategy II (50% H, 50% Q) corresponds to a middle point in the Pareto front representing this set of optimal solutions. Its selection method is as follows: (1) average the optimal solution set corresponding to a Pareto front; (2) select an optimal solution belonging to the Pareto set that is closest to the average value to represent the whole set of optimal solutions. This strategy is a compromise between tree height and weight. The S_layer value gradually increases in strategy II.

TABLE 4

Table 4 Corresponding S_layer value with different wind speeds under strategies I, II, and III.

Through the representative optimal result of strategy II, the plasticity in tree morphological features can be seen. With the increase in wind speed, the diameter at breast height (DBH) increased slightly with wind speed but declined significantly at a strong wind speed; the height and leaf area of the tree progressively declined, which is consistent with the experiments in the literatures (Whitehead, 1962; Telewski & Pruyn, 1998). The partitioning of biomass to the secondary growth increased, especially those for the trunk. The morphology of the tree significantly varies with the gradual acceleration of wind speed. The tree structure corresponding to strategies I, II, and III changes with the wind speed, as shown in Figure 10. The corresponding optimal 3D shapes of simulated poplar trees for strategy II under different wind environments are shown in Figure 11.

FIGURE 10

Figure 10 Simulated tree morphological structures reacting to wind with strategies I, II, and III. (A) Plant biomass (kg); (B) plant height (m); (C) trunk diameter at breast height (cm); (D) tree leaf area (m²); (E) biomass partitioning to secondary growth (%); and (F) ratio for trunk biomass to the sub-branch biomass (%).

FIGURE 11

Figure 11 Simulated optimal tree plasticity for strategy II under different wind speeds. (A) 0.0 m/s; (B) 8 m/s; (C) 16 m/s; (D) 20 m/s (E) 22 m/s. Images are created under calm conditions.

4 Discussion

The interest in this work was initially inspired by the phenomena of thigmomorphogenesis (Hamant, 2013), which is the plant response to biomechanical signals. Such knowledge is applied in a greenhouse by posing artificial wind to have stronger plants. Under lower light, the plant tends to be thinner (with lower biomass) yet keeps tall to catch the light (with modified allometry, i.e., the ratio between length and diameter for the internode). This is a usually unexpected phenotype as plants will be prone to breakage. Blowing the plant with wind can compensate for plants grown in a greenhouse, where light is less intensive than in an open field. An integrated modeling system with light and wind can support quantitative evaluation. However, most modeling studies on plant plasticity concentrate on the response of plants to other environmental conditions like temperature (Liu et al., 2016), light (Kang et al., 2012), or irrigation (Wu et al., 2012), without mechanical constraint. Combining the effect of both light and wind is feasible as light environment simulation has become very common in the FSPM community (Wang et al., 2012). The balance between light and wind can be simulated: to have a better light interception, the leaves need to be more dispersed, while the biomechanical constraint may require a more compact distribution of leaves. An example can be found in Eloy et al. (2017) who explored tree fitness under both light and biomechanical signals. Light competition influences plant height through both the whole tree biomass and architecture, as simulated in Cournède et al. (2008) as well as the allometries; in some cases, stronger light competition leads to higher individual plants (Chen et al., 2010). The current paper deals with trees grown in the open field where the light condition is uniform, and it is not expressed explicitly in the environmental variable E(n) (Eq. S15). Its competition level is expressed with the ground projection area S_p (Eq. S15) as described in Cournède et al. (2008). For the indoor environment as in the greenhouse, the light effect can be considered in the integrated variable E(n) as done in Fan et al. (2015). The sophisticated interplay between light and wind can be better understood by introducing a realistic light component.

Tree reaction to a biomechanical stimulus can be split into two complementary processes: (1) the reallocation of biomass within the whole plant and stem allometry changing and (2) the differentiation of wood cells (e.g., formation of reaction wood) and active biomechanical control of stem shape due to wood fiber maturation strains (Moulia et al., 2006). This paper concentrates on the first complementary mechanism by finding the parameter controlling the reallocation of biomass to secondary growth. It may be called a goal-seeking approach (Takahara & Mesarovic, 2003) in optimizing dry matter distribution with mechanical constraints: instead of simulating directly the impact of wind on plant growth, it searches for suitable growth parameters that lead to optimal biomass and height at given wind load. A similar idea of simulating plant acclimation through optimization can be found in Eloy et al. (2017), with a conceptual model. In the current work, the reaction of the tree is reflected by a single parameter S_layer. Even with a constant S_layer, the allocation changes because of tree organogenesis (Figure 10E). The results fit the situation where trees receive stable wind. When trees grow in varying wind environments, a variable S_layer is expected. Such variable value has been found on maple trees (De Reffye et al., 2017).

The numerical results show that the response of the plant to the wind is mild for a breeze. The response takes place gradually according to Figure 10. It is applicable for the situation where plants receive frequent wind loads for a long period. However, when the tree starts to respond to the wind load is not studied here. How plants react instantly to critical wind load that causes top stem deflection is uninvestigated. Due to the limit of experiment conditions, most experiments were done for low wind speeds. The in silico experiment gives simulated results under strong wind that is hard to achieve and needs further verification. For trees, through measurement and model calibration, all sink–source parameters can be estimated, which is the basis of the in silico experiment. In this work, the model parameters are based on both literature and previous study. The role of wind is not limited to a biomechanical stimulus; it also influences the photosynthesis process by modifying leaf boundary layer conductance, which in turn influences leaf surface CO₂ concentration and humidity. As wind speed increases, the increasing transpiration rate can lead to a decrease in the water content of leaves (Yabuki, 2004). The simulations of such physiological processes are expected to bring more realistic results.

5 Conclusion

The compensating mechanism of plants in reacting to the wind environment is simulated using a system coupling a FSPM with a biomechanical model. A smaller tree (with smaller plant height and leaf area) with higher wind resistance is obtained at higher wind speed, as an emergent property of the system. Secondary growth increases to avoid trunk breakage. This work allows investigating the compensating mechanism of biomass computationally and increases the potential of FSPM for a wider application.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding authors.

Author contributions

MK, TF, and PDR planned and designed the research. HW performed the experiments and analyzed the data. MK, HW, JH, XW, XF and TF wrote and revised the manuscript. All authors contributed to the article and approved the submitted version.

Funding

This work is supported by the Major S&T project (Innovation 2030) of China (2021ZD0113704), the National Natural Science Foundation of China (62076239), the CAS-NSTDA Joint Research Program (GJHZ2076), the Natural Science Startup Foundation of Chongqing Technology and Business University (2056019) and the Science and Technology Research Program of Chongqing Municipal Education Commission (KJQN201900833).

Acknowledgments

The authors thank LIAMA for its support within the cPlant project framework.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2022.971690/full#supplementary-material

References

Almeras, T., Gril, J., Costes, E. (2002). Bending of apricot tree branches under the weight of axillary growth: test of a mechanical model with experimental data. Trees 16, 5–15. doi: 10.1007/s00468-001-0139-1

CrossRef Full Text | Google Scholar

Ancelin, P., Courbaud, B. I. T., Fourcaud, T. (2004a). Development of an individual tree-based mechanical model to predict wind damage within forest stands. For. Ecol. Manage. 203, 101–121. doi: 10.1016/j.foreco.2004.07.067

CrossRef Full Text | Google Scholar

Ancelin, P., Fourcaud, T., Lac, P. (2004b). Modelling the biomechanical behaviour of growing trees at the forest stand scale. part I: Development of an incremental transfer matrix method and application to simplified tree structures. Ann. For. Sci. 61 263–275. doi: 10.1051/forest:2004019

CrossRef Full Text | Google Scholar

Biscani, F., Izzo, D. (2020). A parallel global multiobjective framework for optimization: pagmo. J. Open Sour. Softw. 5, 2338. doi: 10.21105/joss.02338

CrossRef Full Text | Google Scholar

Bonnesoeur, V., Constant, T., Moulia, B., Fournier, M. (2016). Forest trees filter chronic wind-signals to acclimate to high winds. New Phytol. 210, 850–860. doi: 10.1111/nph.13836

PubMed Abstract | CrossRef Full Text | Google Scholar

Calvo, P., Gagliano, M., Souza, G. M., Trewavas, A. (2019). Plants are intelligent, here’s how. Ann. Bot. 125, 11–28. doi: 10.1093/aob/mcz155

CrossRef Full Text | Google Scholar

Chen, S., Zhang, J., Jia, P., Xu, J., Wang, G., Xiao, S. (2010). Effects of size variation and spatial structure on plastic response of plant height to light competition. Chin. Sci. Bull. 55, 1135–1141. doi: 10.1007/s11434-010-0107-5

CrossRef Full Text | Google Scholar

Collet, C., Lanter, O., Pardos, M. (2001). Effects of canopy opening on height and diameter growth in naturally regenerated beech seedlings. Ann. For. Sci. 58, 127–134. doi: 10.1051/forest:2001112

CrossRef Full Text | Google Scholar

Cournède, P.-H., Mathieu, A., Houllier, F., Barthélémy, D., de Reffye, P. (2008). Computing competition for light in the GREENLAB model of plant growth: A contribution to the study of the effects of density on resource acquisition and architectural development. Ann. Bot. 101, 1207–1219. doi: 10.1093/aob/mcm272

PubMed Abstract | CrossRef Full Text | Google Scholar

Coutand, C., Mathias, J. D., Jeronimidis, G., Destrebecq, J. F. O. (2011). TWIG: A model to simulate the gravitropic response of a tree axis in the frame of elasticity and viscoelasticity, at intra-annual time scale. J. Theor. Biol. 273, 115–129. doi: 10.1016/j.jtbi.2010.12.027

PubMed Abstract | CrossRef Full Text | Google Scholar

Coutand, C., Moulia, B. (2000). Biomechanical study of the effect of a controlled bending on tomato stem elongation: local strain sensing and spatial integration of the signal. J. Exp. Bot. 352, 1825–1842. doi: 10.1093/jexbot/51.352.1825

CrossRef Full Text | Google Scholar

Coutand, C., Pot, G., Badel, E. (2014). Mechanosensing is involved in the regulation of autostress levels in tension wood. Trees 28, 687–697. doi: 10.1007/s00468-014-0981-6

CrossRef Full Text | Google Scholar

Davenport, A. G. (1962). The spectrum of horizontal gustiness near the ground in high winds. J. R. Meteorol. Soc 88, 197–198. doi: 10.1002/qj.49708837618

CrossRef Full Text | Google Scholar

Deb, K., Agrawal, S., Pratap, A., Meyarivan, T. (2000). “A fast elitist non-dominated sorting genetic algorithm for multi-objective optimization: NSGA-II,” in Parallel problem solving from nature PPSN VI. Eds. Schoenauer, M., Deb, K., Rudolph, G., Yao, X., Lutton, E., Merelo, J. J., Schwefel, H.-P. (Berlin, Heidelberg: Springer Berlin Heidelberg), 849–858.

Google Scholar

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

CrossRef Full Text | Google Scholar

de Reffye, P. (1976). Modélisation et simulation de la verse du caféier, a l'aide de la théorie de la resistance des matériaux. Café Cacao Thé 20, 251–272. Available at: https://agritrop.cirad.fr/454749/

Google Scholar

de Reffye, P., Hu, B.-G., Kang, M., Letort, V., Jaeger, M. (2020). Two decades of research with the GreenLab model in agronomy. Ann. Bot. 127, 281–295. doi: 10.1093/aob/mcaa172

CrossRef Full Text | Google Scholar

De Reffye, P., Jaeger, M., Barthélémy, D., Houllier, F. (2017). Architecture des plantes et production végétale : Les apports de la modélisation mathématique (Versailles: Quae Edition).

Google Scholar

Diener, J., Rodriguez, M., Baboud, L., Reveret, L. (2010). Wind projection basis for real-time animation of trees. Comput. Graphics Forum 28, 533–540. doi: 10.1111/j.1467-8659.2009.01393.x

CrossRef Full Text | Google Scholar

Eloy, C., Fournier, M., Lacointe, A., Moulia, B. (2017). Wind loads and competition for light sculpt trees into self-similar structures. Nat. Commun. 8, 1014. doi: 10.1038/s41467-017-00995-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Fan, X.-R., Kang, M.-Z., Heuvelink, E., de Reffye, P., Hu, B.-G. (2015). A knowledge-and-data-driven modeling approach for simulating plant growth: A case study on tomato growth. Ecol. Model. 312, 363–373. doi: 10.1016/j.ecolmodel.2015.06.006

CrossRef Full Text | Google Scholar

Fjeld, L. (2012). Modeling MOR in Norway spruce (Picea abies (L) karst.) structural lumber with stands and tree characteristics. Holz Roh Werkst 66, 219–227.

Google Scholar

Fourcaud, T., Blaise, F., Lac, P., Cast E Ra, P., de Reffye, P. (2003). Numerical modelling of shape regulation and growth stresses in trees. II. implementation in the AMAPpara software and simulation of tree growth. Trees 17, 31–39. doi: 10.1007/s00468-002-0203-5

CrossRef Full Text | Google Scholar

Fourcaud, T., Lac, P. (2003). Numerical modelling of shape regulation and growth stresses in trees i. an incremental static finite element formulation. Trees 17, 23–30. doi: 10.1007/s00468-002-0202-6

CrossRef Full Text | Google Scholar

Fourcaud, T., Zhang, X., Stokes, A., Lambers, H., Körner, C. (2008). Plant growth modelling and applications: The increasing importance of plant architecture in growth models. Ann. Bot. 101, 1053–1063. doi: 10.1093/aob/mcn050

PubMed Abstract | CrossRef Full Text | Google Scholar

Fournier, M., Bailleres, H., Chanson, B. (1994). Tree biomechanics: Growth, cumulative prestresses, and reorientations. Biomimetics 2, 229–251. Available at: https://agris.fao.org/agris-search/search.do?recordID=US9506632

Google Scholar

Gardiner, B., Peltola, H., Kellomäki, S. (2000). Comparison of two models for predicting the critical wind speeds required to damage coniferous trees. Ecol. Model. 129, 1–23. doi: 10.1016/S0304-3800(00)00220-9

CrossRef Full Text | Google Scholar

Grotta, A. T., Leichti, R. J., Gartner, B. L., Johnson, G. R. (2005). Effect of growth ring orientation and placement of earlywood and latewood on moe and mor of very-small clear Douglas-fir beams. Wood Fiber. 37, 207–212.

Google Scholar

Guillon, T., Dumont, Y., Fourcaud, T. (2012). A new mathematical framework for modelling the biomechanics of growing trees with rod theory. Math. Comp. Model. 55, 2061–2077. doi: 10.1016/j.mcm.2011.12.024

CrossRef Full Text | Google Scholar

Guo, H., Lei, X., Courn E De, P.-H., Letort, V. (2012). Characterization of the effects of inter-tree competition on source–sink balance in Chinese pine trees with the GreenLab model. Trees 26, 1057–1067. doi: 10.1007/s00468-012-0683-x

CrossRef Full Text | Google Scholar

Hamant, O. (2013). Widespread mechanosensing controls the structure behind the architecture in plants. Curr. Opin. Plant Biol. 16, 654–660. doi: 10.1016/j.pbi.2013.06.006

PubMed Abstract | CrossRef Full Text | Google Scholar

Irvine, M. R., Gardiner, B. A., Hill, M. K. (1997). The evolution of turbulence across a forest edge. Boundary-Layer Meteorol. 84, 467–496. doi: 10.1023/A:1000453031036

CrossRef Full Text | Google Scholar

Jaffe, M. J. (1973). Thigmomorphogenesis: the response of plant growth and development to mechanical stimulation. Planta 114, 143–157. doi: 10.1007/BF00387472

PubMed Abstract | CrossRef Full Text | Google Scholar

Jaffe, M. J., Forbes, S. (1993). Thigmomorphogenesis: the effect of mechanical perturbation on plants. Plant Growth Regul. 12, 313–324. doi: 10.1007/BF00027213

PubMed Abstract | CrossRef Full Text | Google Scholar

Jones, H. G. (1983). Plants and microclimate. a quantitative approach to environmental plant physiology (Cambridge: Cambridge University Press).

Google Scholar

Kang, M. Z., Evers, J. B., Vos, J., De Reffye, P. (2008). The derivation of sink functions of wheat organs using the GREENLAB model. Ann. Bot. 101, 1099–1108. doi: 10.1093/aob/mcm212

PubMed Abstract | CrossRef Full Text | Google Scholar

Kang, M. Z., Fan, X.-R., Hua, J., Wang, H. Y., Wang, X. J., Wang, F.-Y. (2018). Managing traditional solar greenhouse with CPSS: A just-for-fit philosophy. IEEE T. Cybernet. 48, 3371–3380. doi: 10.1109/TCYB.2018.2858264

CrossRef Full Text | Google Scholar

Kang, M. Z., Heuvelink, E., Carvalho, S. M. P., de Reffye, P. (2012). A virtual plant that responds to the environment like a real one: the case for chrysanthemum. New Phytol. 195, 384–395. doi: 10.1111/j.1469-8137.2012.04177.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Koizumi, A., J-i, M., Sawata, K., Sasaki, Y., Hirai, T. (2010). Evaluation of drag coefficients of poplar-tree crowns by a field test method. J. Wood Sci. 56, 189–193. doi: 10.1007/s10086-009-1091-8

CrossRef Full Text | Google Scholar

Kord, B., Kialashaki, A., Kord, B. (2010). The within-tree variation in wood density and shrinkage, and their relationship in populus euramericana. Turk. J. Agric. For. 34, 121–126. doi: 10.3906/tar-0903-14

CrossRef Full Text | Google Scholar

Letort, V., Cournède, P. H., Mathieu, A., de Reffye, P., Constant, T. (2008). Parametric identification of a functional-structural tree growth model and application to beech trees fagus sylvatica. Funct. Plant Biol. 35, 951–963. doi: 10.1071/FP08065

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, B., Asseng, S., Liu, L., Tang, L., Cao, W., Zhu, Y. (2016). Testing the responses of four wheat crop models to heat stress at anthesis and grain filling. Global Change Biol. 22, 1890–1903. doi: 10.1111/gcb.13212

CrossRef Full Text | Google Scholar

Lopez, D., Michelin, S. E. B., de Langre, E. (2011). Flow-induced pruning of branched systems and brittle reconfiguration. J. Theor. Biol. 284, 117–124. doi: 10.1016/j.jtbi.2011.06.027

PubMed Abstract | CrossRef Full Text | Google Scholar

Mayhead, G. J. (1973). Some drag coefficients for British forest trees derived from wind tunnel studies. Agric. Meteorol. 12, 123–130. doi: 10.1016/0002-1571(73)90013-7

CrossRef Full Text | Google Scholar

Moulia, B., Coutand, C., Lenne, C. (2006). Posture control and skeletal mechanical acclimation in terrestrial plants: Implications for mechanical modeling of plant architecture. Am. J. Bot. 93, 1477–1489. doi: 10.3732/ajb.93.10.1477

PubMed Abstract | CrossRef Full Text | Google Scholar

Niklas, K. J., Spatz, H.-C. (2000). Wind-induced stresses in cherry trees: evidence against the hypothesis of constant stress levels. Trees 14, 230–237. doi: 10.1007/s004680050008

CrossRef Full Text | Google Scholar

Peltola, H., Kellomäki, S. (1993). A mechanistic model for calculating windthrow and stem breakage of scots pines at stand age. Silva Fenn. 27, 5504. doi: 10.14214/sf.a15665

CrossRef Full Text | Google Scholar

Perttunen, J., Sievänen, R., Nikinmaa, E., Salminen, H., Saarenmaa, H., Väkevä, J. (1996). LIGNUM: A tree model based on simple structural units. Ann. Bot. 77, 87–98. doi: 10.1006/anbo.1996.0011

CrossRef Full Text | Google Scholar

Qi, R., Letort, V., Kang, M., Cournede, P.-H., de Reffye, P., Fourcaud, T. (2009). Application of the GreenLab model to simulate and optimize wood production and tree stability: a theoretical study. Silva Fenn. 43, 465–487. doi: 10.14214/sf.201

CrossRef Full Text | Google Scholar

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

CrossRef Full Text | Google Scholar

Takahara, Y., Mesarovic, M. (2003). “Goal seeking system,” in Organization structure. Eds. Takahara, Y., Mesarovic, M. (Boston, MA: Springer US), 15–35.

Google Scholar

Taylor-Hell, J., Lawrence, D., Harder, D., Faramarz Samavati, D. (2005). Biomechanics in botanical trees (Master's thesis, The University of Calgary).

Google Scholar

Telewski, F. W., Pruyn, M. L. (1998). Thigmomorphogenesis: a dose response to flexing in ulmus americana seedlings. Tree Physiol. 18, 65–68. doi: 10.1093/treephys/18.1.65

PubMed Abstract | CrossRef Full Text | Google Scholar

Trewavas, A. (2017). The foundations of plant intelligence. Interface Focus 7, 20160098. doi: 10.1098/rsfs.2016.0098

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, H. Y., Kang, M. Z., Hua, J. (2012). “Simulating plant plasticity under light environment: A source-sink approach,” in 2012 IEEE 4th International Symposium on Plant Growth Modeling, Simulation, Visualization and Applications. 431–438. doi: 10.1109/PMA.2012.6524869

CrossRef Full Text | Google Scholar

Whitehead, F. H. (1962). Experimental studies of the effect of wind on plant growth and anatomy. II. helianthus annuus. New Phytol. 61, 59–62. doi: 10.1111/j.1469-8137.1962.tb06274.x

CrossRef Full Text | Google Scholar

Whitehead, F. H., Luti, R. (1962). Experimental studies of the effect of wind on plant growth and anatomy. i. zea mays. New Phytol. 61, 56–58. doi: 10.1111/j.1469-8137.1962.tb06273.x

CrossRef Full Text | Google Scholar

Wu, L., Le Dimet, F.-X., de Reffye, P., Hug, B.-G., Cournede, P.-H., Kang, M. Z. (2012). An optimal control methodology for plant growth-case study of a water supply problem of sunflower. Math. Comput. Simul. 82, 909–923. doi: 10.1016/j.matcom.2011.12.007

CrossRef Full Text | Google Scholar

Yabuki, K. (2004). Photosynthetic rate and dynamic environment (Dordrecht: Springer).

Google Scholar

Yang, G., Cao, W. Q., Liu, D. X., Kang, M. Z., Letort, V., Xing, M. J., et al. (2011). Growth and development simulation based on functional-structural model GreenLab for poplar (Salicaceae). Afr. J. Agric. Res. 6, 3071–3077. doi: 10.5897/AJAR10.855

CrossRef Full Text | Google Scholar