Modelling tree growth taking into account carbon source and sink limitations

Increasing CO2 concentrations are strongly controlled by the behaviour of undisturbed forests, which are believed to be a major current sink of atmospheric CO2. There are many models which predict forest responses to environmental changes but they are almost exclusively carbon source (i.e. photosynthesis) driven. Here we present a model for an individual tree that takes into account also the intrinsic limits of meristems and cellular growth rates, as well as control mechanisms within the tree that influence its diameter and height growth over time. This new framework is built on process-based understanding combined with differential equations solved by the Runge-Kutta-Fehlberg (RKF45) numerical method. It was successfully tested for stands of beech trees in two different sites representing part of a long-term forest yield experiment in Germany. This model provides new insights into tree growth and limits to tree height, and addresses limitations of previous models with respect to sink-limited growth. Author Summary Greenhouse gas emissions, in particular of CO2, have emerged as one of the most important global concerns, and it is therefore important to understand the behaviour of forests as they absorb and store a very large quantity of carbon. Most models treat forests as boxes with growth only driven by photosynthesis, while their actual growth depends also on many other important processes such as the maximal rate at which individual cells can grow, the influences of temperature and soil moisture on these cells, and the control that the tree has on itself through endogenous signalling pathways. Therefore, and with inspiration from process-based understanding of the biological functioning of trees, we have developed a model which takes into account these different factors. We first use this knowledge and additional basic assumptions to derive a system of several equations which, when solved, enable us to predict the height and the radius of an individual tree at a given time, provided that we have enough information about its initial state and its surroundings. We use the Runge-Kutta-Fehlberg mathematical method to obtain a numerical solution and thus predict the development of the height and radius of an individual tree over time under specified conditions.


1
Forests are an important component of the global carbon cycle and are currently a 2 major sink of atmospheric CO 2 ( [1]; [2]). Being able to predict the future response of 3 forests is therefore of great interest. Many models have been used to address this issue, 4 but they are almost exclusively carbon source-driven, with plants at any particular 5 location treated as a box, or boxes, of carbon mainly driven by photosynthesis 6 (e.g. [3]; [4]; [5]). However, it is likely that many other factors, such as the intrinsic 7 limits of meristems and cellular growth rates, as well as control mechanisms within the 8 tree, have large influences on forest responses ( [6]; [7]). A few research groups have 9 addressed the issue of sink-limited growth in a modelling context. The potential for 10 sink-limited growth to affect carbon storage and treeline position was addressed by [8] 11 using a global vegetation model. However, their approach was highly empirical and only 12 addressed temperature limitations. A more mechanistic approach was presented by [9], 13 in which a tree-level carbon balance model was constructed with both carbon source 14 and sink parameterisations. The sink parameterisations were based on thermal-time and 15 included different ontogenetic effects between tissue types and xylogenetic processes for 16 secondary growth. However, while this paper makes a significant contribution, the 17 various parameterisations were very simply incorporated, with no effect of moisture on 18 wood growth, fixed durations for xylem enlargement, and no overall tree growth across 19 years. [10] also presented a tree-level model parameterisation that addressed the effect 20 of growth processes independently of photosynthesis, and in this case looked particularly 21 at soil moisture effects. However, they did not explicitly treat meristem growth, but 22 instead used modified allocation coefficients depending on temperature and soil water. 23 Here we present a new framework for addressing tree growth responses to environmental 24 change, building on knowledge of tree physiology to develop approaches for predicting 25 the development of an individual tree, and thereby enabling a better understanding of 26 forest responses to environmental change than purely source-driven models can achieve, 27 as well as addressing the limitations of previous sink-limited approaches. 28 In this paper we suggest that this objective can be realised using differential models, 29 that is to say evolution models using differential equations. We propose here a first 30 differential model taking into account control mechanisms and the intrinsic limits of 31 meristems and cellular growth within trees together with the carbon balance. The paper 32 is presented as follows. First we make several hypotheses concerning physiological 33 considerations, then we give a description of the model. We show how the model is used 34 to predict tree height and stem volume in two sets of three different stands of beech 35 trees. Finally, we discuss the results and the behaviour of the model in several 36 situations. 37 2 Model and framework 38

39
To describe the development of the tree, we suppose that it has two types of meristems, 40 apical and lateral, and that the apical meristem increases the height through sustaining 41 primary growth, while the lateral meristem (i.e. the vascular cambium) increases the 42 radius through sustaining secondary growth. We recognise that trees will usually have 43 2/20 many separate apical meristems distributed across many branches, but for convenience 44 we treat these together as one apical meristem; we also ignore apical root growth for 45 now. We further assume that the stem can be represented by a cylinder and that the 46 crown (i.e. the branches and foliage) occupies a cylindrical volume and has dimensions 47 that are proportional to the stem dimensions, i.e. the crown depth is proportional to the 48 height of the tree stem and its radius to the radius of the tree stem. In this way we need 49 only obtain information about the radius r and the height h of the stem to describe the 50 growth of the tree. Future developments will include introducing other crown and stem 51 shapes, as well as root meristems. Here our objective is to keep the model as simple as 52 possible in order to explore the realism of its fundamental assumptions.

53
In order to derive the dynamics of r and h we need three types of equations:

54
• constitutive equations -related to the structure of the tree;

55
• carbon balance equations -related to the balance of carbon between the tree and 56 its environment, and the balance within the tree; and 57 • control equations -representation(s) of the intrinsic controls that the tree has on 58 itself. subsequently enlarge and increase in mass ( [11]). We suppose that for each type of 62 meristem, growth is proportional to its volume, that is to say that the mass of carbon 63 allocated to a meristematic region will depend proportionally on the volume of 64 meristem, which determines the rate of new cell production. We also suppose that a 65 tree controls the activities of the meristems and therefore the relative demand for 66 carbon between the meristematic regions: it can favour either height or diameter growth 67 in this way, depending on environmental signals. Moreover, we suppose that this is the 68 only control the tree has on its growth. With this hypothesis we soon get: where α 2 is the ratio of activities between the apical and lateral meristems, i.e. when 70 the activity of the lateral meristem produces 1 g of carbon growth, α 2 g of carbon 71 growth occurs from the apical meristem. It is the control parameter.

72
Note that we suppose here that the volume of the apical meristem is proportional to 73 the top surface area of the (cylindrical) stem and that the volume of the lateral 74 meristem is proportional to the lateral surface area of the stem. This hypothesis is used 75 as it is more realistic than assuming that the volume of the lateral and apical meristems 76 are proportional to the volume of the tree, i.e. that the thickness of the lateral meristem 77 area increases proportionally with the radius and the thickness of apical meristem 78 increases proportionally with height. Although even under this alternative hypothesis 79 we would get the same type of equation with α 2 being the proportion of growth per unit 80 of volume of apical meristem relative to proportional lateral growth, only adding a 81 numerical coefficient in front of α 2 .

82
Knowing the definition of the volume and integrating we get the following 83 constitutive equations: (2) We recognize Equation (2) as the usual allometric relationship where γ is a constant. 85 As it is true for any time t, it can be expressed as h 0 /r 2α2 0 , with h 0 and r 0 the height 86 and radius at some time origin t = 0.

87
Equation (3) is simply the definition of the volume of the stem plus canopy 88 (containing the branches and leaves) in relation to its mass, M c . g 1 is the 89 proportionality ratio that takes into account the density of the stem, the density of the 90 branches and foliage, their contribution to the mass of the tree, and a geometric factor 91 such that g 1 πr 2 h is the mass of the tree. g 1 can therefore be found using the 92 proportionality coefficients we assumed between crown radius and stem radius and 93 crown height and stem height, and their relative densities. 94 2.3 Control mechanism of the tree 95 We want now to take into account the intrinsic control that the tree has on itself, 96 parameterised in our model by α 2 , the ratio between the activities of the apical and 97 lateral meristems. For this we need to consider what happens physiologically (for 98 convenience we use teleological terminology): the tree uses internal controls to be able 99 to adapt its shape and physiology to the surrounding environment, for instance other 100 competing trees (e.g. [12]), through phytohormonal signalling networks ( [13]; [14]).

101
Many behaviours could be taken into account, but here we focus on the behaviour of a 102 tree competing for light by trying to grow taller faster than its competitors by 103 increasing α 2 . The modelled tree is assumed to detect the presence of surrounding trees 104 by sensing the ratio of downwelling red radiation (i.e. wavelengths between 655 and 665 105 nm) to downwelling far-red radiation (i.e. wavelengths between 725 and 735 nm). A low 106 ratio signals potentially more neighbouring trees and therefore a potential threat of 107 shading. In that case the tree reacts by increasing the activity of the apical meristem 108 relative to the lateral meristem, while a higher ratio potentially means no threat and so 109 the tree reacts by balancing the activities of the two meristems as it needs to compete 110 less for light ( [15]). We hypothesise that this reaction can be modelled as: Where α 0 is the highest limit for α 2 beyond which the tree breaks due to mechanical 112 failure, and R a is the ratio of received downwelling red:far-red ratio. α 3 is a scaling 113 parameter such that α 0 exp(−α 3 ) is the value of α 2 when the tree is unshaded (i.e.

114
when it detects no other tree around).

115
This approach is based on the findings of several studies (e.g. [15]; [16]). In 116 particular, [16] observed the laboratory behaviour of two similar young trees grown for 117 21 days with the same intensity of photosynthetically active radiation (PAR), but 118 different red:far-red ratios. They experimentally obtained a relationship which, 119 extended by limited development for short trees, is coherent with the relationship we 120 give in Equation (4).

121
A potential limitation to our approach is that we may not have access to direct 122 measurements of R a . For a lone tree this would not be a real problem: the allometric 123 relationship would be constant in our model and the red:far-red ratio maximal. For 124 individual trees growing in a stand in a forest, however, the situation is different as the 125 trees shade each other. In order to overcome this problem, we can try to obtain R a from 126 more accessible variables. Intuitively, the ratio should depend on the density of trees 127 and the area they shade horizontally. This ratio should also allow a finite maximum 128 when the density tends to 0, that is to say there are no neighbouring trees. Therefore, 129 we could propose the following dependency: With d the density of trees (ha −1 ), r n the average radius of the neighbouring trees 131 (m), and a, c 1 , and c 2 are three constants such that log(a + c1 c2 ) is the maximal value of 132 the ratio (i.e. when there are no trees around). This formula is inspired by several 133 studies (e.g. [12]; [17]), where red:far-red ratios were measured together with the density 134 of the stands, and which suggest that for young trees the ratio has a logarithmic 135 increase with the inverse of the density. However, it is likely that this is our most 136 limiting relationship as we are here trying to model the internal control that the tree 137 has on itself, which is probably the most complicated phenomenon to take into account, 138 and likely the least-well understood component of our model (cf. [14]; [18]). It is also 139 easy to see that the r n variable makes this formulation challenging to use: even with an 140 average value, finding the solution of h = γr 2α2 becomes difficult and in practice would 141 increase the complexity of the differential equations. Also, it should be taken into 142 account that this formula becomes very sensitive to any small variation of r due to 143 measurement error or the consequences of a particular year's conditions. This is the 144 reason why, to keep the model simple, we consider that when we do not have access to 145 measurements of the red:far-red ratio, α 2 is constant for each stand, as it would be for a 146 lone tree. When the value of this constant will be needed, it will be estimated as a 147 typical value for the species we consider from a stand with similar density.

149
The variation of carbon in the tree over time can be written as the sum of the carbon 150 used for volume growth and the carbon stored: where M c is the carbon used for volume growth, S c = Sg 1 πhr 2 the carbon stored, 152 and S is defined as carbon stored per unit volume of stem. This distinction has a 153 physical interpretation and therefore in this model the density of the tree can then be 154 seen as a standard density related to the carbon used for growth plus an increase of 155 density due to the storage. 156 We also know that the variation of carbon is given by: where A is the carbon assimilation due to gross photosynthesis minus a foliage 158 respiration part proportional to photosynthesis, while C is a constant grouping a part of 159 the respiration and a litter component proportional to the volume, together with the 160 geometrical factor g 1 π (such that C 1 /g 1 π is the respiration and litter component 161 constant). We assume here that there is no other litter component, that is to say that 162 litter is completely proportional to the volume of the stem. We note that this also 163 means that the litter is proportional to the volume of the crown with our previous 164 assumptions. Finally, C 2 is the same type of constant but for the storage. A is assumed 165 to depend on PAR (Q1), the concentration of CO 2 in the atmosphere (C g ), temperature 166 (T ), water potential (ψ), and the dimensions of the tree (h and r).

167
Therefore: The system is so far incomplete as there needs to be a relationship that controls 169 which part of the carbon is allocated to growth and which part is allocated to storage. 170 Also, the volume growth of the tree is sustained by meristematic cell division and 171 5/20 enlargement, with a maximal rate determined by intrinsic physiological limits or 172 environmental controls. 173 We assume that the tree grows in volume as much as it can, as long as it keeps 174 enough storage to avoid being in danger, such as for repairing damage and surviving 175 poor growing seasons. Therefore, if there is enough carbon assimilation then the rate of 176 growth is equal to the meristem-sustained growth-rate limit and the storage can 177 increase using the difference of these two. Also, if the carbon assimilation rate becomes 178 too low but there is enough storage, then the storage delivers carbon to the 179 meristematic region to maintain a rate of growth equal to its maximal sustainable rate. 180 Then density decreases as the carbon of the storage is used to maintain the volume 181 growth. This is likely to occur for instance when buds appear while photosynthesis is 182 not yet high enough to support the maximal growth potential due to meristem

192
If there is enough light or storage, i.e. then:

196
We assume the following dependencies of photosynthesis on water potential, 197 atmospheric CO 2 concentration, PAR, temperature, radius, and height: 198 1. Photosynthesis increases with CO 2 concentration and PAR exponentially and 199 tends to saturate: 6/20 2. As a first approximation photosynthesis of each leaf decreases linearly with a 201 reduction in water potential. We suppose also that it decreases linearly with 202 height: Where: With ψ soil the water potential in the soil and ψ lim the limit water potential below 205 which any photosynthesis cannot be performed due to a lack of turgor pressure in 206 the leaves and increased probability of (xylem) cavitation. R h is the hydraulic 207 resistance per unit length.

208
So, integrating across the full tree, we get: Note that β is a numerical coefficient that only depends on the ratio between the 210 crown depth and the height. If we suppose that the crown starts at around half 211 the height of the tree then β = 3/4. Also, the dependency of both height and 212 water potential is given above. 213 3. We assume a dependency with temperature given by an asymmetric parabola: With T i being 21 • C if T > T opt and 25 • C otherwise. We set here T opt = 19 o C.  4. We suppose that photosynthesis is proportional to the surface area of the crown, 217 and so proportional to r 2 .

218
Overall we get: 2.6 Maximal meristem-sustained growth 220 Intuitively, as meristem-sustained growth occurs at a cellular level, we assume that the 221 limit to meristem-sustained growth per volume increases with temperature up to a 222 certain temperature and increases with water potential (cf. [20]). Therefore, we assume 223 the following dependencies: 2. The meristem-sustained growth limit decreases with height as:

7/20
With ψ soil the water potential in the soil and ψ lim2 the limit water potential below 228 which meristem-sustained growth cannot occur. Physically this gives us a 229 maximal height h 2 due to the limit on meristematic growth. 230 We define therefore R max0 the reference maximal meristem-sustained growth per 231 volume independant of water potential and temperature by: 3. As stated previously, we suppose that the volume of meristems is expressed as:

234
V me = (t 1 πr 2 + 2πrht 2 ), where t 1 is the effective thickness of the apical meristem 235 and t 2 is the effective thickness of the lateral meristem. Using the previous elements we finally get: depending on whether carbon supply from photosynthesis plus storage is sufficient to 240 meet demand for growth (Eqns. 23 and 24) or not (Eqns. 25 and 26). The most accessible are probably Cr and Qr. These are usually well known and can 245 be obtained by looking at the curve of photosynthesis activity as a function of 246 atmospheric CO 2 mixing ratio (in ppmv) or as a function of PAR, respectively.

247
A max can be obtained by direct measurements of the carbon exchange at a given 248 height, ambient CO 2 concentration, temperature, PAR, and compensating for 249 respiration as in [21]. Acting similarly and measuring litter with litter traps gives access 250 to C. These values can also be deduced from net primary production and gross primary 251 production and the typical litter flux rates. Note that here we use data coming from 252 forest stands, although in the model we suppose that we have an individual tree without 253 shading because data for lone standing trees are seldom available. This means that we 254 may underestimate the growth of isolated trees, but at the same time the values of 255 A max and C obtained enable to take into account the shading effect by neighbours of a 256 forest with the same density, which is a hidden parameter of this model (i.e. not 257 explicitly taken into account). Also, depending on the typical time length of the study, 258 we choose a model timestep and therefore we use averaged values of the parameters on 259 8/20 this timestep. In this analysis, we chose one year as the timestep although smaller 260 timesteps could be considered.

261
Finding h 1 is equivalent to finding the minimal water potential that can sustain 262 photosynthesis. h 2 is its equivalent with respect to growth. Combined, these two 263 parameters take into account all water potential-dependant limits in trees such as 264 cavitation in the xylem conduits, minimal turgor potential to maintain the shape of the 265 leaf and enable stomatal opening, etc. While h 1 (or equivalently ψ max1 ) has been 266 extensively studied (e.g. [22]; [23]; [24]), and is believed to be equal around 90 m for 267 redwood trees ( [23]), h 2 seems to be less known. However, knowing h 1 and the other 268 parameters as well as the effective limit height for the considered species, we can deduce 269 h 2 relatively easily. 270 Q 10 is a physiological parameter that can be found by studying the activity of 271 meristem cells in a controlled environment. We can find R max0 , at least approximately, 272 by looking at the tree rings of many trees across several geographical locations and 273 assuming that the largest ones correspond to the maximal increment dr achievable in a 274 year, which is directly linked to R max0 , the maximal meristem-sustained growth rate 275 under perfect temperatures and water potentials. Then we can deduce the value of 276 R max0 , or at least its order of magnitude.

Symbol
V alue Quantity g 1 477 kg.m −3 Density of the equivalent cylinder modelling the tree C/(g 1 π) 72 kg.m −3 .y −1 Respiration and litter coefficient Reference maximal corresponds to a maximal meristem-sustained radius expansion of 1.5 cm.y −1 lateral growth per area independent of water potential and temperature t 2 /t 1  3 Results

280
The model described above is ordinary, first order, non-linear, and has no widely known 281 explicit analytical solution in the case where 2α 2 is not an integer. Therefore, it was 282 solved numerically by the Runge-Kutta-Fehlberg (RK45) method ( [25]). Coherence 283 both with the analytical solution in simplified cases and with other solving methods (i.e. 284 the trapezoidal rule and numerical differentiation formula) was tested.

285
As some parameters can depend, a priori, on the considered species, we derived the 286 parameters explained previously for beech trees from the data provided in [21] and [26]. 287 Beech trees were chosen because of their importance in European forests, especially in 288 Germany where they represent more than 15% of the forest and in France where they 289 comprises 15% of the non-conifer forest. Therefore, a great number of studies have 290 focussed on beech trees, enabling easier comparisons with our current study.

291
As the computation of the parameters is based on independent measurements that 292 have sometimes a non-negligible margin of error, we allowed the parameters thus 293 obtained to vary with a 20% margin to account for this error, and to allow a small 294 adaptation of the model (as it is only a simplification of reality) by fitting them with 295 reference measurements of a beech tree stand in Fabrikschleichach, Germany. As R max0 296 was more coarsely estimated we allowed it a variation of 50% in order to maintain the 297 right order of magnitude. Also, as no precise data on the soil water potentials over the 298 tested period A.D. 1870 to 2010 were available for these stands, we allowed a 50% 299 variation of h 2 , the height where meristem-sustained growth cannot occur due to cell 300 turgor limits. We obtained then a final vector of parameters for our model using a 301 non-linear least-square method (trust-region-reflective algorithm: 302 http://uk.mathworks.com/help/optim/ug/equation-solving-algorithms.html). 303 We then used our model to predict the time evolution of height and volume of an 304 average tree in several stands of beech trees. It should be emphasised that there is a 305 substantial difference between fitting and prediction. Fitting a model to some data gives 306 a partial understanding of the data but does not usually enable prediction for another 307 tree or even for future data points as the fitting is done on a restricted set of 308 measurements. Prediction, on the other hand, is much harder to achieve as it supposes 309 the use of previously-derived parameters and environmental information (e.g. air 310 temperature, soil moisture, CO 2 concentration, etc.) to predict the measurements.

311
Usually most models are fitted rather than predicting as it is much harder to predict 312 anything without any fitting, although when it works it gives much more information: a 313 better understanding and a reliable way to deduce future measurements before they 314 occur even for other trees (e.g. height or volume prediction). 315 We first considered three stands of beech trees that are part of a long-term  D. 1990. For each stand we considered the averaged height of the 100 trees 326 with the largest stem diameter at breast height. The first stand was the reference stand 327 used to partially fit the parameters (within the 20% margin).

328
As expected, our partially fitted model agrees well with the measurements from this 329 stand (Fig. 1 a; R 2 = 0.975). Predictions with the model were then performed for the 330 two other stands (Fig. 1 b,c). Figure 1. Size of the modelled tree with age at the first site. Continuous curves are computed with the model while circles correspond to experimental measurements. a) represents the results of the partial fit on the reference stand. b) and c) represent the results of the prediction with the parameters obtained by partial fitting.

11/20
The results give a very good prediction for the third stand (R 2 > 0.94), and a small, 332 although non-negligible, error for the third one (R 2 = 0.905). However, this could be at 333 least partially explained by the change in the density of trees during the 120 yr due to 334 different rates of thinning in the different stands (i.e. from 3,500 tree ha −1 to 200 335 tree ha −1 in the first stand, 6,400 tree ha −1 to 200 tree ha −1 in the second, and 2,400 336 tree ha −1 , to 200 tree ha −1 in the third). Therefore, the red:far-red radiation ratio 337 changes with time in the real forest while it was assumed constant in our model.

338
To test the generality of the model, it was also used to predict the height of 339 dominant trees in three other stands at a different site in the Spessart, Hain (50.0 • N; 340 9.3 • E), a low wooded mountain range in central Germany, located approximately 100 341 km from the first site (Fig. 2). As data were missing concerning the water potential 342 during the growth period from A.D. 1870 to compare with the first stand, we allowed 343 the model to adapt its value of h 2 within 20% by keeping the parameters derived 344 previously and fitting only this one on the first of the three new stands. Then the set of 345 parameters were used to predict the heights of the three stands according to the model. 346 Again, to avoid the potential error induced by the different variations of the red:far-red 347 ratios within the stands, we restrict our analysis to data when the density was lower 348 than 1,000 tree ha −1 . As before, the predictions are in strong accordance with the 349 observations for all three stands, with a nearly exact correspondence with the 350 measurements (R 2 = 0.98, R 2 = 0.91, R 2 = 0.97 respectively), except for the second 351 stand. 352 We should note that the correspondence is very close but not exact, and that the 353 data we have might limit the accuracy: firstly, from our equations we can see that there 354 is a propagation of error from the initial conditions. An error on the initial height of 1 355 m could imply an error of the predicted height of up to 2.5 m after 50 years even if the 356 model were perfect. Also, no data were available for the difference of water potentials 357 between the stands on a same site, while it seems logical that for a certain density of 358 trees, the stands with higher density will have less water available per tree. Therefore, 359 this could induce a small error in h 2 and h 1 that could be corrected if we knew the 360 average water potential during the growth period in each stand.

361
Finally, mortality is not addressed by the model, although the model does produce 362 cessation of growth under stressful circumstances. Therefore tree mortality due to 363 critical events (e.g. disturbance, pathogens, etc.) other than a limit on growth due to 364 the external conditions, can create discontinuities in the data which are not captured by 365 the model. Although the continual use of the 100 trees with largest DBHs tends to 366 average this discontinuity (as the probability that a large change in those 100 trees 367 occurs in one year is low), this might have created another limit to the accuracy. For 368 instance, we note that we have observed in the 2010 measurements, i.e. after 188 years, 369 the death of several of the 100 trees with the largest DBHs, which seems to have caused 370 a decline of the average height of the 100 living trees with the largest DBHs as taller 371 trees were replaced by smaller trees in this group, which is a result of size-related 372 mortality dynamics (cf. [27]).

373
If we wanted to address further this question, especially if we wanted to simulate 374 many individual trees that interact together in a forest, we could multiply the height h 375 by a random process M (t) that would be equal to 1 when the tree is alive and 0 when 376 the tree is dead. Then the probability p(M(t)=0|M(s<t)=1) of transition between a 377 living tree and a dead tree could be derived using experimental measurements and even 378 be a function of the size of the tree. So far, though, this doesn't seem to be needed in 379 the present application.

380
So far, the only things we assumed to be known for predicting the height are the 381 external variables (i.e. temperature, water potential, PAR, and CO 2 concentration), the 382 species of the tree, and the allometric relationships for each stand. No other knowledge 383 Figure 2. Plot of the size of the tree with age at the second site. Continuous curves are computed with the model and start when the tree density is lower than 1,000 tree ha −1 , while circles correspond to experimental measurements. Parameters for beech trees are obtained from the first site and h 2 was allowed to vary up to 20% to compensate for the lack of data concerning water potential from A.D. 1870 on the two sites. a) to c) represent respectively the results for stands 1 to 3.
13/20 about the trees was used to perform the predictions. Measurements in themselves were 384 used only to compare our predictions with reality. 385 We could avoid assuming known allometric relationships for each stand and try to 386 estimate them as described in Section 2.3 ("Control mechanisms of the tree"). In that 387 case we would still get good results with R 2 around 0.8 − 0.9. The small difference of 388 accuracy between this case and the previous one is mainly due to the model becoming 389 more sensitive to measurement error when we estimate the allometric relationship, and 390 also to the estimated allometric relationship being a constant as explained in Section 391 2.3, whereas the allometric relationship in reality changes with changes in tree density 392 and therefore time. Lastly, our model for α 2 is probably too simple as we also explained 393 in Section 2.3. Nevertheless, we can note that the accuracy of the results is still very 394 good considering that we try here to predict the heights of three stands without 395 knowing anything but the species, the CO 2 concentration, the temperature, the soil 396 moisture, the PAR, the number of trees per ha, and the initial tree dimensions. This 397 last point also underscores an additional motivation to build such a relatively simple 398 model: collecting data about the external variables can be difficult. Therefore, data can 399 be limiting and in order to be practical the model should stay as simple as it can 400 relative to the required data, provided that it can still give accurate results. 401 We chose to predict height as there is usually a relatively small error of measurement 402 in height compared to other variables. It should be noted that the same procedure is 403 possible taking instead the stem volume as the quantity to predict. However, the results 404 would be less precise because they will be susceptible to cumulative measurement errors 405 of the radius and the stem height, among other reasons, but they would still be accurate 406 (the highest error is around 25% after 75 years) and are presented in Fig. 3. Also, we 407 chose to predict height instead of radius or diameter as the measurements performed on 408 the radius were the DBH which might be different from the equivalent cylindrical 409 diameter, especially for the 100 trees with the largest DBH, and this would have added 410 an artificial error of measurement.  Firstly, it is known that forest stand growth dynamics have significantly changed 415 since 1870, as shown in [28], and many experiments have been conducted to measure the 416 impact of climate and CO 2 change on forests. Climate and CO 2 change also shows an 417 impact on this individual tree model, as the model depends on external conditions such 418 as temperature and soil moisture, as well as CO 2 concentration. For instance, using the 419 model with a constant temperature CO 2 concentration (both from 1870) rather than 420 the actual temperatures CO 2 concentrations gives biased results and a slower growth of 421 trees (Fig 4). This seems coherent with both the idea that climate change had an 422 impact on forest stand growth and that present forest stands are growing more rapidly 423 than comparable stands before.

424
Also, if we consider a small tree in a forest, shaded by the canopy above, it would 425 receive only a small amount of light. Therefore the model would predict a very limited 426 maximum height and the tree would remain near this height with very small growth 427 until a big tree nearby falls. Then the light received would be higher and enough for the 428 tree to grow taller very quickly. If the tree were very close to the dead big tree then the 429 aperture is large enough and the tree would reach the canopy. If not, and if the tree is 430 still significantly partially shaded compared to the rest of the trees, then its limit height 431 would be lower than the forest height and it would reduce its growth progressively until 432 it were close to this height. This behaviour would seem to be qualitatively in accordance 433 with reality [29], but due to a lack of data we have not been able test this quantitatively. 434 Figure 3. Plot of the stem volume with age at the first site. Continuous curves are computed with the model and start when the tree density is lower than 1,000 tree ha −1 , while circles correspond to experimental measurements. a) and b) represent the results of the prediction with the parameters obtained by partial fitting. c) represents the partial fit on the reference stand. Our model starts from physiological considerations about trees and several hypotheses 436 to derive differential equations that are solved over time to obtain the evolution of a 437 tree and predict its size. Although the model is fairly simple and far from taking every 438 possible parameter of the tree into account, it seems to obtain very good results and 439 predictions which agree closely with observations from different, independent, stands of 440 trees with different surrounding conditions.

441
Interestingly, our approach gives rise to two limiting tree heights: h 1 , the limit 442 height above which photosynthesis cannot occur, and h 2 , the height above which growth 443 is limited by the capacity of the meristematic regions themselves to grow at low water 444 potential. These two heights represent different processes: h 1 includes the factors that 445 limit photosynthesis such as minimal turgor potential and possible xylem cavitation (i.e. 446 source limitations), whereas h 2 takes into account intrinsic meristem factors such as the 447 physiological limits of cell division and growth limits such as cell wall extension when 448 the water potential is too low (i.e. sink limitations).

449
While h 1 (and, more generally, source processes) is often considered as the only 450 limiting height to explain the maximal height of trees (e.g. [23]), under some external 451 conditions of light and temperature it could be that h 2 (i.e. sink activity) becomes 452 limiting (cf. [22]). Although in this model, water potential is still the critical factor 453 which limits the height of trees, the amount of light received and the temperature 454 during the growth period also play roles that could produce one or other limit and 455 influence its value. Hence, our approach could give a new and more complete 456 understanding of the limits of tree height, and growth more generally, through a 457 balanced consideration of both source and sink limitations.

458
As a first dynamic model for an individual tree using such considerations about 459 control mechanisms and the maximal meristem-sustained limited rate of growth, 460 improvements could probably be made by introducing additional parameters, especially 461 taking the roots into account explicitly, and improving the function of the red:far-red 462 ratio. Nevertheless, the model already gives very good results while being relatively 463 simple, allowing its wider use such as through allowing many trees to interact with each 464 other via shading, water potential, etc., and obtaining a model for the whole forest. So 465 far the partial fitting of the parameters within bounds of 20% takes less than 120 s on a 466 desk-top computer, and computing height over 100 yr less than 0.1 s, which suggests 467 that its use as a growth simulator within a forest model could be achieved with today's 468 computer resources.

469
The relative simplicity of the model and its extremely encouraging results should 470 stimulate further use of similar differential models in the future. Starting from 471 physiological considerations to derive differential equations and an individual tree model 472 as a cornerstone for simulating a whole forest may considerably improve our 473 understanding of forest behaviour, including, but not limited to, the prediction of the 474 overall reaction of forests to increasing CO 2 concentrations in order to address the 475 future terrestrial carbon sink.

476
To conclude, we have proposed a differential model for modelling tree growth over 477 time under external conditions such as temperature, soil moisture, and CO 2 478 concentration. This model takes into account not only photosynthesis and the carbon 479 balance but also meristem behaviour and cellular growth limits. We established a 480 procedure to parametrise the model with measurable quantities and reference 481 measurements. This model seems not only to fit data very well but also to give accurate 482 predictions for tree height and tree volume.