Application of the Functional-Structural Tree Model LIGNUM to Growth Simulation of Short-Rotation Eastern Cottonwood

The functional-structural tree growth model LIGNUM was developed as a general research tool that can be applied to several tree species. The growth simulation of short-rotation eastern cottonwood (Populus deltoides Bartr. ex Marsh.) inherits the basic LIGNUM modeling concepts including modular tree structure, L-system-based description of structural development, and carbon budget. New developments of LIGNUM model in this study were the incorporation of a biochemically-derived photosynthesis submodel; nested time steps for simulating physiological processes, structural development, and annual biomass production; incorporation of field-measured weather data for modeling the response of physiological processes to environmental variation; and application of a Monte-Carlo voxel space submodel for simulating the stochasticity of tree growth and improving computational efficiency. A specific parameter system was applied for modeling P. deltoides growth in the central Missouri, USA, environment. This adaptation of LIGNUM was applied on modeling growth of P. deltoides in a short-rotation agroforestry practice. The simulated height and biomass growth were close to field observations. Visualization of simulation results closely resembled the trees growing in an open site. The simulated response of tree growth to variations in photon flux input was reasonable. The LIGNUM model may be used as a complement to field studies on P. deltoides in short-rotation forestry and agroforestry.


Introduction 1.Eastern Cottonwood, a Versatile Short-Rotation Forestry Tree
Eastern cottonwood (Populus deltoides Bartr.ex Marsh.) is a fast-growing deciduous tree common in eastern North America from Quebec in the north to Florida in the south.Both pure P. deltoides and its hybrids with other Populus spp.have been intensively studied for short-rotation forestry in North America (Alig et al. 2000).It has been estimated that short-rotation poplar and cottonwood plantations could supply 40% of the hardwood pulp needs of the U.S. and reduce current pressure for exploitation of natural forests for pulpwood (Alig et al. 2000).There is an increasing interest in poplar wood for use in lumber and oriented strand board, which have long service life and, thus, carbon sequestration potential (Pallardy et al. 2003).Introduced P. deltoides is common in short-rotation forestry in India (Puri et al. 1994, Singh et al. 1997).It has also become a prime agroforestry species in Northern India (Singh 1998).Silvopastoral practices combining P. deltoides with a mixed understory of native grasses and legumes are used for erosion control in New Zealand (Guevara-Escobar et al. 2000).
The usage of P. deltoides in Northern-American agroforestry practices is also increasing.
Populus deltoides is well-suited to flood plain forests.It grows best on moist, well-drained, fine sandy loams or silt loams (Fig. 1a).Its widespreading root system and rapid growth make it useful for soil erosion control along river banks and flood plains.The most important feature of P. deltoides is its resistance to flood damage (Fig. 1b).During the 1993 floods of the Mississippi and Missouri rivers, the biggest in a century, mortality of P. deltoides was low in areas that were under floodwaters for less than a month, and at least 97% of the large (DBH > 20 cm) trees survived even 250 days of flooding (Larsen, unpublished data).Therefore, P. deltoides seems to be a promising species for flood plain forestry practices and it is one of the main species in the flood plain agroforestry projects of the University of Missouri Center for Agroforestry.In this contribution, we describe the new application of the functional-structural tree model LIGNUM (Perttunen et al. 1996, 1998, 2001, Perttunen and Sievänen 2005), which was developed for P. deltoides in this context.

Functional-Structural Tree Model LIGNUM
Physiological processes and morphology are two critical properties that must be accommodated to plant growth models (Kurth 1994, Perttunen et al. 1996).Both process-based and morphological modeling approaches have advantages and limitations (Sievänen et al. 2000).The processbased models do not have reference to the threedimensional crown structure, which may affect radiation interception and carbon allocation.Morphological models, on the other hand, are limited in their capacity for describing tree physiology pertinent to growth processes.These two basic approaches are to a certain degree complementary and can be combined into functional-structural tree models (Kurth 1994, Sievänen et al. 2000).
Several approaches for functional-structural plant modeling have been presented in special issues of Annals of Forest Science 57(5/6) in 2000, New Phytologist 166(3) in 2005, and Functional Plant Biology 35(9/10) in 2008.
The LIGNUM model (Perttunen et al. 1996, 1998, 2001, Lo et al. 2001, Perttunen and Sievänen 2005) is designed to simulate the growth relationships among organs within a tree, from the uppermost leaves to the root system.The model has been developed primarily to enhance understanding of mechanisms underlying tree growth and a tree's response to different growing environments.It is used to explore the relationship of physiological processes within a tree to its crown structure and growth.
The key for integrating morphological and process-based models in LIGNUM is to consider a tree in terms of repeated basic units and use them to model both metabolism and spatial structure (Perttunen et al. 1996).A tree is treated as a collection of a large number of basic units that correspond to organs in a real tree (Perttunen et al. 1998).Sievänen et al. (2000) called these structural modules idealized elementary units (IEU); they are a common feature of functional-structural plant models.Four structural modules are used to represent the tree crown in LIGNUM (Fig. 2): tree segment, shoot-forming bud, axis, and branching point.These units describe both the architectural development of the tree crown in space and the carbon allocation in terms of the metabolism taking place in these units.Physi-ological processes, including photosynthesis and respiration, can also be explicitly related to the tree structures in which they occur.
In the earlier versions of the model, annual radiation interception, which was assumed to be linearly proportional to the annual photosynthetic production, was collected based on detailed tree structure (Perttunen et al. 1996(Perttunen et al. , 1998)).The carbon budget in LIGNUM linked photosynthesis, respiration, and tree growth.Structural development was implemented by L-systems (Perttunen and Sievänen 2005).The root system in LIGNUM was modeled with one parameter denoting its mass (Perttunen et al. 2001).
The earlier LIGNUM versions combined a very exact description of tree structure and radiation attenuation within a tree canopy with a rough estimate of incident radiation and photosynthetic production.The annual integral of incident radiation was used by assuming standard overcast distribution throughout a year (Perttunen et al. 1998), and only annual whole-tree photosynthetic production was calculated (Perttunen et al. 1996).This approach produced good results when simulating early development of slow-growing trees like Scots pine (Pinus sylvestris L.; Perttunen et al. 1996Perttunen et al. , 1998)), jack pine (P.banksiana Lamb.; Lo et al. 2001), and sugar maple (Acer saccharum Marsh.; Perttunen et al. 2001).
In functional-structural modeling, it would be desirable to apply a process-based CO 2 exchange model for estimating C available for growth, and use variable weather data input in order to better simulate the growth responses of a tree to environmental variability.Unfortunately, the submodel used for estimating radiation interception in earlier LIGNUM versions, i.e. pairwise comparisons between all leaves of a tree for defining mutual shading (Perttunen et al. 1998), makes short time step simulations of weather and CO 2 exchange impossible because of extreme computational complexity even for small trees.Thus, the general framework of the LIGNUM model was retained when applied for growth simulation of P. deltoides, but submodel improvements were made to enable short time step simulation of physiological processes of this fast-growing species.

Objectives
The aim of this contribution is to describe adaptation of the functional-structural tree model LIGNUM to a new broad-leaved species with submodel improvements and parameter adaptation.This adapted model has been applied to simulate the growth and yield of Populus deltoides under a flood plain environment.The specific objectives of the project were to: 1) develop a nested time step system with different update periods for physiological processes and structural development for accommodating the fast-growing property of P. deltoides; 2) develop a submodel for using real, field-measured weather data as input for modeling the response of physiological processes to the environment; 3) apply a biochemically-derived photosynthesis model (von Caemmerer and Farquhar 1981, Long 1991, De Pury and Farquhar 1997) for estimating leaf CO 2 exchange; 4) improve the computational efficiency of photosynthetic photon flux interception in LIGNUM model for enabling short time step process simulations; 5) improve the estimation of C allocation to roots and include the first approximation of C rhizodeposition (Högberg and Read 2006) to LIGNUM; and 6) apply a parameter system specific to P. deltoides growing in a mid-Missouri, USA, flood plain environment.
These adaptations can be used for further extensions of LIGNUM development in simulating the growth and yield of P. deltoides in agroforestry systems and short-rotation plantations.In order to reach objective 4, the voxel space approach (Sievänen et al. 2008) was modified for estimating photon flux interception within a broad-leaved tree canopy.

Model Flow Diagram
The modeling time step is the time period for an individual cycle of the process being modeled.The modified LIGNUM model applies three nested time steps in tree growth simulation: an annual time step, a structural time step, and a physiological time step.The annual time step matches the natural annual cycle of tree growth.The intermediate structural time step is used to capture the structural developments of P. deltoides during a growing season.In the flow of the entire model, canopy structure, root biomass, and root respiration are updated every six weeks so that its changing effects on photosynthetic photon flux density within the tree canopy, photosynthetic production, and tree growth may be captured.The physiological time step functions at the scale of a fraction of a day.
Fig. 3 shows the three nested model loops.The inner physiological loop collects the intercepted photon flux, photosynthesis, and leaf respiration every 30 minutes and sums up total photosynthesis for use in subsequent structure updates.The medium structural loop integrates the net CO 2 assimilation from the inner loop to obtain the net biomass growth after subtracting for respiration.The net biomass growth is allocated into each part of the tree for updating the tree structure.The structure update is implemented four times each year.The outer loop encompasses all these procedures within the annual growth cycle.

Intermediate Time Step: Derivation of Tree Structure
In the overall model configuration, the canopy structure is updated periodically during a growing season so that it may affect photon flux density patterns within the tree crown structure, as well as photosynthetic production and tree growth.
The structural time step is a fraction of the annual cycle, selected so that canopy expansion during the time period is detectable.At the same time, the structural time step is long enough to simplify the model calculations.
In natural tree crown development, each bud can develop into a new branch in a growing season, which would include several segments, branching points, and a new apical bud.Each branching point carries one lateral bud, which creates a leaf and an axillary bud, which may lead to further branch development in the following growing season.New leaves emerge from branch tips one by one with their supporting segments.Naturally, the derivation of segments is part of a continuous process resulting in tree growth.This continuous process is computationally divided into several phases, the structural updates.Only the apical bud of each branch is set as active during a single growing season and serves as the basis for deriving new segments in each structure update step.The apical bud contains several leaf primordia that form a leaf and a dormant axillary bud at the leaf base.These dormant buds become active during the subsequent growing season (Figs. 3  and 8).All leaves are shed at the end of the growing season.Root dynamics are also modeled at the intermediate time step.The details of C dynamics in roots are given in section 2.5.3.research articles

Short Time Step: Interception of Photon Flux and Photosynthesis
The  so that there was a photosynthetically important change in unshaded Q because of changes in solar position.Net CO 2 assimilation is accumulated and allocated into tree structure in the structural update step.Computational efficiency and accuracy were considered in the determination of time steps for physiological processes and structural update in addition to biophysical criteria.

Incident Photosynthetic Photon Flux Density
Real weather data were used as the input for the model.The measured Q was divided into diffuse (Q d ) and direct (Q b ) components using the relationship between measured and potential global radiation (Weiss andNorman 1985, Nygren et al. 1996).The model was parameterized by Weiss and Norman (1985) in Nebraska under similar climatic conditions that prevail in mid-Missouri.An object on the ground receives diffuse flux from all directions of the visible hemisphere, but the relative importance of different sky regions depends on geographic and temporal atmospheric conditions.The distribution of photon flux in the whole sky was by means of the Firmament submodel (Perttunen et al. 1998(Perttunen et al. , 2001)), which deals with direct and diffuse fluxes separately.The details of the submodel are described in Appendix 2.

Voxel Space
Interception of direct and diffuse photosynthetic photon flux is determined by incident photon flux, mutual shading of leaves, and the optical properties of leaves.In earlier versions of LIGNUM, radiation transmission in a tree segment was calculated according to Oker-Blom and Smolander (1988), and mutual shading among all segments was estimated using analytical geometry (Perttunen et al. 1998(Perttunen et al. , 2001)).The voxel space model was introduced into LIGNUM to reduce computational complexity associated with the radiation interception model.Instead of calculating the mutual shading of each leaf in the tree canopy, voxel boxes are used as the calculation unit, which simplifies calculations substantially.
Voxel space is the space in which the tree grows.The 3D tree growing space is divided into small imaginary cubic boxes called voxels.Just as a pixel is the basic pictorial planar element in a two-dimensional system, a voxel is the basic volumetric element in a three-dimensional system.In an individual tree simulation, the size of voxel space is extended to accommodate the growing tree.The number of voxels increases with tree growth, which saves calculation redundancy in juvenile tree simulations.Interception of Q is calculated based on the voxels.The 3D tree structure generated by LIGNUM defines the position of each leaf within the voxel space.Every tree component occupies one or more voxels (Fig. 4).
The voxel acts as the unit for calculating photon flux interception.All leaves in one voxel are added up as the content of that voxel.The greater the leaf area is in a given voxel, the greater the likelihood of that voxel intercepting incoming photon flux.A beam might be intercepted by one of the voxels in its pathway or go completely through the tree growing space.The bigger the voxel, the less accurate the resolution of voxel space becomes.The voxel size is selected according to the branch elongation speed to simplify the calculation of photon flux interception while assuring that the modeling results reflect reality.Given a standard overcast sky (Ross 1981), the voxel space receives diffuse flux from each sector in the sky in a specific direction.The direct flux is emitted from the apparent direction of the sun.Intercepted Q b and Q d in each voxel are summed to yield the total incoming photon flux for the leaves in that voxel.

Big Leaf
There may be one or more leaves with different area and angle, and/or tree segments in a given voxel.Each leaf has an angle α between the leaf's normal direction, or the direction perpendicular to the plane of the leaf blade, and the incoming photon flux.The interception of photon flux by the leaf is the result of incoming flux multiplied by cosα.In order to simplify the cosine correction of photon flux interception, one big leaf is generated to represent all leaves in a given voxel.The normal of the big leaf determines the general direction of exposure for all leaf components in the voxel.The big leaf normal is calculated as the leaf-area-weighted sum of all individual leaf normals in the voxel: where N is the general big leaf normal; LA i is the area of leaf i, and N i is the normal of leaf i.The leaf area of the big leaf is the total area of the leaves in the box.Let α N denote the angle between the normal of the big leaf and the incoming photon flux direction; the photon flux received by the big leaf is calculated as the product of the original flux density multiplied by cosα N .When the area of the big leaf is bigger than the projected area of voxel A v , the big leaf is folded to several layers to fit in the voxel.

General Rules for Photon Flux Interception
The stochastic Monte-Carlo voxel space model is used in LIGNUM to calculate the mutual shading among voxels in the path of photon flux.The interception of photon flux is affected by whether the flux passed through other voxels before reaching the focal voxel.The voxels in front of the focal voxel in the photon flux path are called shading voxels.The calculation of the intercepted total photon flux, Q i , proceeds according to the following steps: I) Mutual shading is based on the relative geometric position of voxels filled with leaves.For each voxel, a path is traced towards each sky sector and a list of voxels on that path is generated.Each path is used to calculate a fraction of the Q d from each sky region.Solar position determines the direction of Q b .The path corresponding to this direction is used to calculate the nents may either be intercepted by shading voxels or they reach the focal voxel.III) The big leaf is used to describe optical properties of leaves (Monteith 1965, Ross 1981).IV) The Monte Carlo voxel approach is employed for determining whether or not a voxel intercepts the photon flux.In so doing, the projected big leaf and voxel surface areas perpendicular to flux direction from each sky region n -A l (n) and A v (n), respectively -are calculated (Fig. 5).The probability (p) of the directional flux interception within the voxel is the ratio A l (n) / A v (n).The Bernoulli process (Bertsekas and Tsitsiklis 2002) is used to determine whether or not a flux is intercepted by the voxel i with the probability of p i .

Bernoulli p probability p probability p ( ) -
The Bernoulli(p i ) function has result either 1 or 0.
The result 1 for the Bernoulli process means that the photon flux is intercepted by the voxel.V) For a particular voxel, the absorbed photon flux is a fraction of the incoming photon flux.When exposed to photon flux, a single layer of leaves is considered to absorb 80% of the incoming photosynthetic photon flux, with 10% being transmitted through the leaf, which may be intercepted and absorbed by another layer of leaves (Monteith 1965, Ross 1981).Reflection and scattering of the photosynthetic photon flux in the canopy are not considered in the current model version.In the case of non-interception, the whole flux passes through the voxel.When there is more than one layer of leaves in the voxel, 90% of the Q is considered to be intercepted by the voxel biomass.VI) The total photon flux intercepted by a focal voxel is calculated as the sum of all fluxes -direct and diffuse -from each sector.The intercepted photon flux is shared equally by all leaves in the voxel.
Four interception scenarios are possible in the pathway of the photon flux (Fig. 6): 1) the photon flux may penetrate through all the shading voxels and reach the focal voxel unaltered; 2) the photon flux may be transmitted through the foliage of one shading voxel, which has one layer of leaves and the focal voxel receives 10% of incident Q; 3) the photon flux may be absorbed by the foliage of one shading voxel, which has more than one layer of leaves and the focal voxel is considered completely shaded; or 4) the photon flux may be transmitted through the foliage of one shading voxel and be intercepted by another shading voxel, resulting in complete shading of the focal voxel.

Standard Diffuse Photon Flux Interception
Each calculation of leaf photosynthesis with 30-minute interval is implemented based on the current tree structure.However, since tree structure remains the same during a structural time step, each voxel gets a similar fraction of diffuse photon flux during this period.Thus, the standard diffuse photon flux interception is applied for simplifying the calculations.Within each tree structure update cycle, diffuse photon flux interception by each voxel is only calculated in the first photosynthetic process with new structure to create a standard distribution of intercepted diffuse photon flux.The standard distribution describes the diffuse photon flux intercepted in each voxel as a fraction of total diffuse photon flux incident on the canopy (Appendix 2).For the other photosynthetic processes within the same structural time step, these fractions are used to calculate the interception of diffuse photon flux in each voxel.

Photosynthesis Submodel
Leaf photosynthesis (P) was modeled as a function of intercepted total photon flux (Q i ) and ambient temperature (T a ) according to Farquhar's approach (von Caemmerer and Farquhar 1981, Long 1991, De Pury and Farquhar 1997, see Appendix 3).The photosynthetic production calculated is the net CO 2 assimilation of a leaf (A) after subtracting leaf respiration.Net CO 2 assimilation is integrated over the leaf area of the whole tree; at night time, A is negative because of foliage respiration.Net primary production (NPP) is estimated as the difference between foliar net CO 2 assimilation and respiration by the woody parts of the tree.

Lindenmayer System for Eastern Cottonwood
Lindenmayer systems (L-systems) are essentially parallel rewriting systems operating on strings, i.e. sequences of symbols.An L-system is defined by an alphabet of symbols and a set of rules called productions, each rule replacing a symbol 'a' by a string 's'.The rewriting begins from the initial start string called axiom.The formal definition of the classes of L-systems, string rewriting, and interpretation can be found in Prusinkiewicz and Lindenmayer (1990) and Prusinkiewicz et al. (1997).
In the plant modeling context, the symbols represent both the botanical entities (internodes, buds, flowers, leaves etc.) of the growing organism and the topological and the geometrical information how the entities are connected with each other.The development of the branching structure of a plant is achieved by applying the rules to the given axiom and then repeating the rewriting process in discrete time steps.L-systems and their extensions have been applied to tree architecture (Kurth 1999, Mutke et al. 2005), development of a shrub under different environmental conditions (Salemaa and Sievänen 2002), root growth (Mech and Prusinkiewicz 1996), and plant-insect interaction (Hanan et al. 2002), among other modeling questions.
The LIGNUM model has been interfaced with the language L (Perttunen and Sievänen 2005), which is an extension of L-systems.Based on the definition of the language L (Prusinkiewicz et al. 1999), Karwowski (2002) created the original parser for the language and further implemented the L+C language (Karwowski and Prusinkiewicz 2003), including features not present in L such as fast linear time information transfer in the simulated plant.
The key idea in L-systems is the rewriting of symbols.The language L follows the same idea.A module (corresponding to a symbol) in L has a name and it can take any number of arguments of any type in the C++ programming language (Stroustrup 1997).In addition, L can embed C++ for computations.A rule in L consists of a predecessor module possibly with its context (informally, the context of a module means what is on its left or right hand side) and a production defining the successor string.A special module Start corresponds to the axiom.For example, a simple L-system of two modules A() and B() demonstrating the language L can be defined as follows: The first four strings of modules after first four derivation steps produced by the system are, starting from the axiom: A(), B()A(), A()B()B()A() and B()A()A()B()A()B()B()A().Note the parallel rewriting.
To model higher plants like trees, L must have modules that represent the various plant parts, their geometrical aspects and their topology, especially the branching structure of the plant.
For the geometrical interpretation of an L-system defined in L, a turtle interpretation (Prusinkiewicz and Lindenmayer 1990) has been implemented.A turtle is simply a 6-tuple (x,y,z,  H,  L,  U) representing its state.The Cartesian coor- dinates (x,y,z) determine turtle's position in three dimensional space.The three perpendicular unit vectors (  H,  L,  U ) specify turtle's heading, direc- tion to the left and the direction up respectively.They define the current orientation of the turtle and satisfy the vector cross product equation 7).In the initial orientation at the origin, the turtle heads up.
Three modules in L can change turtle's orientation in space: To move the turtle to the direction pointed to by the heading by length s, the module F(s,…) has been defined.The ellipsis denotes that F can take additional arguments.To model branching, the two modules SB() and EB() delimit a branch denoting the beginning and the end of a branch respectively.
As the module F can be interpreted as the cylindrical tree segment (Fig. 2), another module B(…) for the bud is needed to create a tree-like structure in L that can be converted with the turtle interpretation to LIGNUM for metabolic processes.The algorithms for the two-way communication linking the LIGNUM model and the L-system for considering interactions between tree structure and functioning have been described in detail earlier (Perttunen and Sievänen 2005).The L-system for P. deltoides is listed in Appendix 5. Specifically, the module F(s,o) denotes a tree segment of length s with branching order o.The module B(PoplarBudData) denotes a bud with information passed in PoplarBudData data structure between LIGNUM and the L string (for example, if the bud is active, dormant or dead, its branching order and its position in space).
The initial tree consists of a 10-cm-long segment, an active apical bud and two active lateral buds (Appendix 5, the module Start).After the integration of photosynthesis and respiration (Fig. 3) two rewriting steps of the L-system string are needed for each structure update step: one for segment elongation and the second for formation of axillary buds.Turn, Pitch, Roll, and F in the L string creates a trajectory for the turtle in space that forms the tree architecture (Fig. 9).
The rule F(s,o) < B(d) (Appendix 5) defines the elongation.If there is a segment in the left context (denoted by <) of an active bud, the bud is rewritten to a new initially 30-cm-long segment followed by Split() marker symbol and to itself.Depending on the branching order of the active bud, the rewriting implements either the elongation of the main axis or side branches (Fig. 8).
After the elongation, the L string is converted to LIGNUM structure (Fig. 2) for carbon allocation (Fig. 3) to determine dimensions of the new tree segments, calculate the induced diameter growth, and determine if the buds are dead or alive according to available photosynthates (Eq.4).The resulting new tree structure is transmuted back to an L string.
After the carbon allocation the rule F(s,o) > Split() (Appendix 5) forms the initially dormant lateral buds.The segments with Split() marker in their right context (denoted by >) are rewritten to a series of segments of equal length (cf.Fig. 15) and new dormant lateral buds are generated with azimuthal distributions and bifurcation angles according to Figs. 13 and 14.Only the three first cases for the main axis and branches are shown in Appendix 5 because similar repetitive fifty rewrites for new segments up to 1 m for both main axis and branches would take too much space.Finally, the last rule Split() (Appendix 5) removes these marker modules from the L string during the formation of lateral buds due to the parallel rewriting.
The final L string is again converted to the LIGNUM structure for leaf insertion to model a new photon flux interception step (Fig. 3) and metabolism in the new crown architecture.There are four structure updates during a growing season, with a 6-week-interval between them.In the following growing season, all dormant lateral buds resume activity and have the potential to develop a new branch.The structural development steps of a 3-year-old P. deltoides produced by the L-system without metabolic processes are shown in Fig. 9.

Carbon Allocation
Along with structural development, tree biomass is allocated in each part of a tree.Biomass increment (G) is determined as the difference between the net C assimilation by leaves (P); respiration (R t ) in the stem, branches, and roots; and rhizodeposition of C (D): Carbon is used as the currency of energy flow in the program.Photon flux interception in each part of the tree is calculated with self-shading among tree compartments.The net CO 2 assimilation is the leaf photosynthesis minus leaf respiration.Whole-tree net C assimilation (A t ) is calculated by integrating net C assimilation by individual leaves (A l ) over the whole canopy during a calculation period.The accumulated net C gain after subtracting the total maintenance respiration of woody parts is the C available for tree growth.
The net C gain is then allocated to biomass all over the tree in the form of segment elongation, diameter increase, and foliage and root growth.Growth respiration, which is assumed to be proportional to biomass production, is subtracted to calculate net increase in tree size.Tree growth leads to a change in tree structure, which in turn affects photon flux interception in the next growth cycle.The relationship between NPP, the wholetree net C assimilation and the maintenance and growth respiration of other tree compartments is depicted as: Fig. 9.The structural updates from 1-1 (1 st update of the 1 st year) to 3-1 (1 st update of the 3 rd year) without metabolism for Populus deltoides as defined by the L-system in Appendix 5 corresponding to the schematic presentation in Fig. 8.The elongation phase (the rule F(s,d) < B(d) in Appendix 5) was changed to create new segment lengths between 2 and 4 cm and between 1 and 2 cm in the main axis and branches, respectively.S: The initial 10-cm-tall tree with the active apical and two active lateral buds.
where R mi is the maintenance respiration of the ith compartment of the tree, which is proportional to the existing biomass in each compartment, and k is the proportion of growth respiration (R gi ) out of the C available for growth.The factor k is assumed to be equal to all tree compartments and R gi is, consequently, proportional to biomass production: The growth respiration of foliage is embedded in Farquhar's photosynthesis model according to eq. (A3-1) in Appendix 3. The same function (eq.A3-12 in Appendix 3) was used for leaf night and daytime respiration rate.Growth occurs only when the whole-tree net C assimilation exceeds the maintenance respiration demand by woody organs and fine roots.If maintenance respiration exceeds the whole-tree net C assimilation during a six-week structure update period, the tree dies.Dormant season respiration is not computed because of the low temperatures prevailing from November through March in mid-Missouri where the field data were gathered (Table 1).The net primary production is allocated (Perttunen et al. 1996, 1998, Perttunen et al. 2001): where iW n (kg) is the C allocation to new tree segments and buds, iW o (kg) C allocation to the woody biomass increment by radial growth of existing tree segments, iW rc (kg) is the C allocation to roots including rhizodeposition, and λ is a unitless parameter to control the sizes of new tree segments.We assume that all NPP is allocated to tree growth during the structure update, i.e. there is no carbon storage pool for photosynthates.The total demand for photosynthates caused by a new tree segment at the time of its emergence is not known.It can be calculated only by traversing the tree from the tip of the branch to base of the tree and assess the induced diameter growth.Thus, given certain sizes of new segments the eq.( 6) does not necessarily hold.The balance of the eq.( 6) can be solved for λ iteratively with Bisection method (Press et al. 1992).During the iteration λ alternately elongates and shortens the segments and generates demand for diameter growth throughout the tree according to pipe model principle (Shinozaki et al. 1964) defined for P. deltoides (Appendix 4).Also, the C allocation to roots iW rc depends on new foliage and is therefore controlled by λ (Section 2.5.3).Above-ground structural development and the C allocation are affected by branching, by the relative sizes of tree segments and radiation environment.The effect of the relative size of a segment is described by its vigor index (Perttunen et al. 2001, Nikinmaa et al. 2003).The greater the vigor index value, the more NPP is allocated on the segment.Thus, the thicker the axis the more relative growth potential it has.The effect of branching is described by the Gravelius order (Gravelius 1914) of the tree segment.The higher the Gravelius order, the less NPP is allocated on the segment.A tree compartment that absorbs more radiation grows faster than one that is more shaded by other compartments.The radiation index acts as a factor in segment elongation along with vigor index and Gravelius order.The radiation index ranges from 0 to 1 following the photosynthetic production curve as a function of absorbed photon flux (Fig. 10).The radiation index is set to 1 for segments under saturating photon flux; to 0 for segments under a photon flux density below the radiation compensation point of photosynthetic rate; and from 0 to 1 according to the radiation index curve for all segments growing in intermediate radiation environments (Fig. 10).The effects of vigor index, Gravelius order, and radiation factor were multiplicatively applied to segments' growth in C allocation within a tree as described in Appendix 4.

Root Biomass
Root dynamics for P. deltoides include birth and senescence, respiration distinctly for fine and coarse roots, and rhizodeposition, i.e. the loss of C from roots to soil.It is the sum of several processes: 1) root cap and border cell loss; 2) death and lysis of root cells (cortex, root hairs etc.); 3) C flow to symbiotic associates in soil (e.g.mycorrhizae); 4) leakage of solutes from living cells (root exudation); and 5) insoluble polymer secretion from living cells (mucilage) (Jones et al. 2009).At an annual time step, part of the root mass, W r , dies as defined by the parameter s r .The living root mass after a growing season (W rnew ) is: At each structural update the new roots are created as follows.First, we assume that the share of new carbon allocated for the roots (iW rc ) in the carbon balance is only partially used to growth due to the rhizodeposition.We assume that 30% of C allocated to roots (defined by 1-a rd ) is used to create new root mass (∆W r ): where a rd is the proportion of rhizodeposition out of total C allocated to roots.To allocate C for roots, we assume that new foliage requires new roots defined by the parameter a r .Secondly, the roots that died at the end of the growing season must be replaced.Thus, as there are four structure updates per growing season, the C requirement for new roots is: where W fnew is the mass of the new foliage created during the structure update.Finally, for respiration, we assume that 50% of the root mass is fine roots, defined by a cr (Coyle and Coleman 2005).Thus, respiration of the roots (R mr ) in the intermediate time step is: where rc mr and rf mr are annual maintenance respiration rates for coarse and fine roots, respectively.
3 Field Data The soil nutrition, especially N content, strongly decreases with increasing depth (Lu 2006).Parameterization data were measured in an alley cropping experiment (Fig. 11).Trees were planted with 6 × 18 m spacing in April 2001 in association with white clover (Trifolium repens L.).The P. deltoides clone was a selected Midwestern industry clone.No within-row canopy closure had occurred before the measurements in summer 2002 and 2003.Thus, the alley cropping was considered to be an open-growth plantation.

Weather Data
Weather data were collected using an automated weather station (Campbell Scientific, Logan, UT, USA).The weather station recorded global radiation, photosynthetic photon flux density, net radiation, wind speed, soil temperature, air temperature, air humidity, and rainfall.A Li-Cor LI-190SB Quantum sensor (Li-Cor, Inc., Lincoln, NE, USA) was used for measuring photon flux density; Kipp & Zonen CM6B pyronometer with heating/ventilation unit (Kipp & Zonen BV, Delft, The Netherlands) was used for global radiation; and HMP45C-L integrated temperature/relative humidity probe (Vaisala Ltd, Helsinki, Finland) for air temperature.Photosynthetic photon flux density and temperature were used in the current model, and global radiation was applied for estimating direct and diffuse components of incident photon flux density according to the empirical model of Weiss and Norman (1985; see Appendix 2).Other weather data were not used in this study.All radiation data were cosine corrected.
Data were measured every 10 seconds and the average values were recorded at 10-minute intervals by a CR23X data logger (Campbell Scientific).The date for weather data was recorded as Julian dates.The logger time was converted from central time to real solar time for ease of calculations.

Branching Pattern
Two-and three-year-old trees were sampled for branch structure.Three trees with large, average, and small diameter at breast height (DBH) were selected.The measurements included branch base diameter, branch bifurcation angle, internode seg- ment diameter, and segment length.The data were organized by branching order that was set following developmental topology (Berntson 1996).The lowest order was the tree stem, which was set as order 0. The main branches, or those attached directly to the stem, were assigned to order 1 and the order number increased with increasing branch bifurcation.The tree heights were also measured for trees from two-to five-years old for model validation.Woody biomass (stem, branches, and roots) of 4-year-old trees was estimated from the stem basal radius applying the model developed by Nygren and Pallardy (2008).
Height of branch insertion point, diameter of main branch, diameter of stem below main branch, insertion angle of the main branch, and the angle of branch bending were recorded for each main branch (Fig. 12).Detailed structural measurements were implemented for twelve main branches selected from each tree.Sampling was stratified by selecting four branches from the lowest third of crown length, four from the middle third, and four from the uppermost third.The structural measurement in branch included branch basal diameter, diameters in each bifurcation point, length of internode segment, and bifurcation angle of sub-branches (Fig. 12).The bifurcation was tracked in each selected main branch until the last order of the branch.Onesided surface area of selected matured leaves was measured by scanning.

Parameters of the Photosynthesis Submodel
The photosynthesis submodel was parameterized based on unpublished data of Pallardy.The model parameters were determined using eastern cottonwood saplings derived from cuttings from a dense plantation adjacent to the alley cropping stand used for field measurements of this study.The cuttings were grown in a greenhouse and transferred to a growth chamber for gas exchange measurements.The growth conditions of the cuttings and the measuring procedure were similar to those applied by Dowell et al. (2009).Gas exchange of mature leaves was measured at 25°C ambient temperature and various CO 2 concentrations in the sample cuvette with a Portable Photosynthesis System (LI-6400, Li-Cor, Inc., Lincoln, NE, USA).Generation and analysis of A l -C i curves was conducted following the procedures recommended by Long and Bernacchi (2003).
The parameters determined from the data included maximum rate of carboxylation (V cmax25 ), maximum rate of electron transport (J max25 ), photosynthetic rate under ambient CO 2 concentration and saturating radiation (P max25 ), and dark respiration rate (R d25 ).The cuttings represented three clones of eastern cottonwood.Because no significant differences in gas exchange characteristics between clones were detected (Lu 2006), we used the parameters derived from the clone with V cmax25 and J max25 close to the average (clone 1112 of the Missouri Department of Conservation).

Weather Conditions
The growth environment of the studied P. deltoides is in the temperate zone of North America.

Branching Angles
Branching angle included both horizontal and vertical angle.In the case of the angle between the stem and a main branch, the horizontal angle was recorded as azimuth, i.e. clockwise from North, and vertical angle was the inclination referenced to horizon.The angles between a main branch and subbranches were measured with respect of the main branch.The branch angles affect branch development and partly determine structural pattern of the tree crown.Although the main branch azimuth appeared to have a small southern bias (Fig. 13), the branch azimuth distribution was considered to be initially random at the time of bud emergence.Main branch inclination angles ranged mostly from 30° to 70° with an average value of 44° and a standard deviation of 14° (Fig. 14).The above features of branch angles provided the branch derivation pattern for P. deltoides in LIGNUM.

Growth Dynamics
Branch growth dynamics include branch elongation, leaf emergence, and maturation for new branches.New branches exhibited active growth from April to October, with final length averaging about 0.7 m (Fig. 15b).New leaves emerged along with branch elongation.The average leaf numbers increased from May and reached a maximum number of ca.35 before leaf fall began in August (Fig. 15a).
The time interval between emergence of two     concatenated leaves was about 1 to 2 days.It took about 10 days from emergence to full leaf maturity.The leaf mass may be expressed in terms of leaf area and specific leaf area.The specific leaf area varied during the leaf development; the measured value for mature leaves was 15 m 2 kg -1 .

Model Parameterization
LIGNUM was parameterized for P. deltoides based on field data and a literature review for published data (Table 2).Sapwood maintenance respiration ratio was set to be 0.015 kg [C] kg -1 [C] y -1 and growth respiration was set to be 0.125 kg [C] kg -1 [C] of new growth (Rauscher et al. 1990).Fig. 16 shows the net CO 2 assimilation rate estimated by Farquhar's approach (Appendix 3) at temperatures 15°C, 25°C, and 35°C using the CO 2 exchange parameters from Table 2.The temperatures shown cover the typical range of growing season temperatures in mid-Missouri.The radiation index curve (Fig. 10) was converted from net CO 2 assimilation curve at 25°C (Fig. 16) and was used in segment C allocation submodel.
The structural modeling time step and voxel size were determined by branch elongation (Fig. 15).The average length of segment internodes was 0.055 meter.The selection of voxel size in LIGNUM model was related to segment length so that voxel size was bigger than a single segment in order to contain one or more segments completely; otherwise no gain in computational complexity would have been achieved.At the same time, the resolution of voxel space was set small enough for accurately describing the tree structure.Branch growth rate during a growing season was stable and almost linear (Fig. 15).It took about six weeks for a branch to grow 0.2 m.Thus, the structural time step and voxel size were harmonized by setting the structural time step to 6 weeks and voxel side to 0.2 m so that, on average, a new layer of voxels was created in the surface of the growing canopy during each structural update step.This resulted in a reasonable computational load.A 5-m-high tree filled about 832 voxels, 130 of them with leaves inside, and it took about eight minutes to simulate tree growth to this size.
Carbon allocation to the root system and rhizodeposition were parameterized using literature data.First, it was assumed that fine roots comprise half of the root mass of young P. deltoides (Coyle and Coleman 2005), and root senescence corresponds to fine roots only.These fine roots must be replaced together with the need of fine root increment for fulfilling the water and nutrient requirements of the growing foliage.Coleman et al. (2000) derived the following equation for estimating soil respiration rate (r s ) in a P. deltoides plantation as a function of soil temperature (T s ): where T s is in °C and r s in µmol m -2 s -1 .Applying the average growing season soil temperature in our study site, 19.24 °C, to Eq. ( 11) gives soil respiration rate 2.537 µmol m -2 s -1 ; thus, ca.46 mol m -2 or 552 g m -2 in a growing season.Further, it has been estimated that fine root respiration in a P. deltoides plantation accounts for 20% of soil respiration (Horwath et al. 1994) and fine root density of P. deltoides is 356 g m -2 (Coleman et al. 2000).Combination of these factors with the estimate of soil respiration derived from Eq. ( 11) gave the average fine root respiration rate of 0.3 g [C] g -1 [BM], which was used in the LIGNUM model (Table 2).We recognize that this is a rough first approximation of fine root respiration in the LIGNUM model and its implications will be discussed in Chapter 5.2.Rhizodeposition was determined by closing the C balance, and it was set to be 70% of C allocated to the roots (cf.Högberg et al. 2002).

Simulation Results and Model Validation
The simulation results included tree height, and total stem, branch and root biomass.LIGNUM may also produce additional information such as branch length, and number of leaves and segments on request.In addition to numerical data, simulation results are presented as a tree image.The  2.
tree visualization provides heuristic expression of the simulated results, which can be compared and validated with the visual appearance of real tree growth.The results of LIGNUM may be evaluated by examining the appearance and the general shape of the tree simulated with computer graphics and by comparing the numerical data generated by the model with the field data.
Comparison between selected tree characteristics simulated with LIGNUM and measured in the field site is presented in Table 3.Because LIGNUM includes stochastic estimation of the attenuation of photon flux density within the canopy, each simulation with a given parameter set gives slightly different results.Thus, all simulations were run five times and we reported the mean and standard deviation.The field data were from trees of different size growing in the site where parametrization data were collected.There was a fairly good agreement for height between the simulated results and field observations, although there was a slight underestimation in year 2 and slight overestimation in years 4 and 5 (Table 3).The simulated woody biomass for 4-year-old P. deltoides was higher than the mean in the field, 5.21 kg per tree, but lower than the maximum observed, 8.36 kg per tree.The simulated biomass increment pattern and development of root:shoot ratio are shown in Fig. 17.
We also simulated the development of a partially shaded young P. deltoides receiving 75% of the incident photon flux density (Table 4, Fig. 18).Height growth was only slightly reduced in comparison to the completely unshaded tree (Table 3).The estimated woody biomass at 4-year-old trees, 5.02 kg per tree, was very close to the field average.In a nearby dense P. deltoides plantation, the average height of three P. deltoides clones varied from 8.01 to 9.12 m and woody biomass from 6.44 to 6.95 kg per tree in 5 years after planting (Dowell et al. 2009).These height measurements were in close agreement with the simulation under partial shading but biomass production was lower than simulated (Table 4).
The simulated height growth followed the same pattern as the simulated total leaf area of a tree (Fig. 18).No reduction in either height or leaf area growth rate was visible during the 8-year- simulation.This development is reflected in the visual representation of an 8-year-old P. deltoides in an open growth environment (Fig. 20b).Fig. 19 visualizes the LIGNUM simulation results for two growing seasons, starting from a 10-cm-long cutting through the fourth structure update of the second growing season.Leaf primordia in apical buds generate leaves but dormant axillary buds in leaf base do not form new branches before the following growing season.Inclusion of metabolic processes, like photosynthesis and respiration, and C allocation with segment split makes these LIGNUM simulated trees bigger and more variable in shape than the trees based on deterministic L-system rewriting rules only in Fig. 9.This emphasizes the importance of combining tree structure and functions in growth modeling.
Fig. 20 shows simulated 4 and 8 years old P. deltoides trees growing in an open site under climatic conditions of mid-Missouri.The images of the simulated tree showed strong apical dominance.Each new branch attached several new leaves, which split the entire branch into several segments.The axillary bud at the base of each leaf derived a new branch in the next growing season.Branches in lower position had fewer leaves because of the self-shading by upper crown.Stem diameter tapered from base to tip of the tree.The modeled structure of P. deltoides appeared to satisfactorily reflect the growth pattern in the field (Fig. 11).The crown shape developed in LIGNUM simulation was close to the actual canopy form.
The LIGNUM model was also run with reduced radiation level for a 4-year-old tree growth by computationally shading sky sectors (Perttunen et al. 2001).Under Q reduced to 50% of measured photon flux density, the simulated tree had much less leaf area (Fig. 21) compared to such a tree without shading (Fig. 20).Simulated tree height with 75% Q was 6.26 m, which is ca.10% less than the height (6.95 m) at full sun at the same age.Under Q reduced to 50%, the simulated tree attained a height of 5.37 m in four years.The total leaf area under 50% shade was 29.9 m 2 as opposed to 70.9 m 2 at full sun (Fig. 21).Total woody biomass at 4-year-age was 8.32 kg under full sun and 2.79 kg under 50% shading.

Simulation Running Time
Tree growth was simulated during a growing season from May to October.The length of the simulated growing season was 180 days.Because photon flux interception and leaf photosynthesis were computed every 30 min, there were 8 832 calculations (184 × 24 × 2) for direct photon flux interception and leaf photosynthesis for each voxel during a growing season.In addition, one diffuse photon flux interception was calculated in    each structural time step, or four times in a growing season, for all photon fluxes coming from all sectors of sky.Conducting these calculations by pair-wise comparisons between all leaves of a tree as in earlier LIGNUM versions (Perttunen et al. 1996(Perttunen et al. , 1998(Perttunen et al. , 2001) ) would have been impossible.Thus, the voxel space approach was used.
The number of voxels was significantly lower than the number of leaves (Fig. 22).Calculation of diffuse photon flux involved all sky sectors (801 in simulations presented in this contribution) but it was conducted only once for each structure update cycle.Direct photon flux that originated from one sky sector only was computed every 30 min.Thus, the computational time depended on the product of the number of voxels containing tree foliage and mean number of shading voxels on the direct beam's route.Because multiplying constant numbers does not increase the complexity of simulation, the computational time of the voxel space model was directly proportional to the total number of the voxels comprising the tree growing space, as opposed to the computational time in pairwise comparison, the latter of which was equal to the second power of the number of leaves.For example, a 5-year-old tree had 8083 leaves and 218 leaf-filled voxels (Fig. 22).The pairwise comparison of leaves (Perttunen et al. 2001) would have required 8,083 2 = 65,334,889 computational operations for each estimation of photosynthesis.A reasonable assumption is that a voxel is shaded on average by 7 voxels (Lu 2006) giving 7 × 218 = 1,523 computational operations for the same task in the voxel space.Because the other operations take about the same time in the pair-wise comparison and voxel model, the computational complexity of the former is ca.42,814 times higher in this example.

Sensitivity Analyses
Nine parameters were subjected to sensitivity analysis: maximum rate of electron transport at 25°C (J max25 ); maximum rate of carboxylation at 25°C (V cmax25 ); day-time leaf respiration rate at 25°C (R d25 ); shoot maintenance respiration rate (r ms ); wood basic density (ρ w ); specific leaf area (SLA); root to foliage C allocation ratio (a r ); ratio of rhizodeposition to total C allocated to roots (a rd ); and coarse root to total root biomass ratio (a cr ).The leaf gas exchange parameters were measured from seedlings grown in a greenhouse and it is, thus, important to know the sensitivity of the model to them.Morphological parameters and respiration rate might also affect tree growth in terms of resultant biomass and volume.The sensitivity of LIGNUM to changes in these parameters was evaluated based on their effects on simulated tree height, foliage mass, root mass, and total woody biomass.
The parameters were varied by ±25% of their original value.In most cases, the output was affected similarly by variation in a given parameter (Table 5): all either increased or decreased.LIGNUM was least sensitive to sapwood maintenance respiration, coarse root to total root biomass ratio, and wood basic density.The increase of the values of parameters J max25 and V cmax25 strongly increased the output values; the increases of R d25 , SLA, and a r reduced the output values.
Combined sensitivity analyses were accomplished by applying two parameters together (Tables 6 and 7).Four parameters with strong Table 5. Sensitivity analysis of LIGNUM to variations in values of physiological and structural parameters of Populus deltoides.The numbers are percentage change in an output variable of four-year-old P. deltoides in comparison to results using parameter values in Table 2. Symbols: J max25 = maximum rate of electron transport at 25°C; V cmax25 = maximum rate of carboxylation at 25°C; R d25 = day-time leaf respiration rate at 25°C; r ms = shoot maintenance respiration rate; ρ w = wood basic density; SLA = specific leaf area; a r = root to foliage C allocation ratio; a rd = ratio of rhideposition to total C allocated to roots; a cr = coarse root to total root biomass ratio.Table 7. Sensitivity of LIGNUM to variation in parameter combinations using the total wood mass (stem, branches, roots) of four-year-old Populus deltoides as an index.Each two parameters are combined as a pair for sensitivity analysis.Numbers in the matrix indicate the the total wood mass change with 25% increases in column head parameters and 25% reductions in row head parameters.The value in the cell is the increase or reduction rate in percentage compared to the total wood mass simulated using original parameter values from Table 2. Symbols: J max25 = maximum rate of electron transport at 25°C; V cmax25 = maximum rate of carboxylation at 25°C; R d25 = day-time leaf respiration rate at 25°C; and SLA = specific leaf area.-45.0%Table 6.Sensitivity of LIGNUM to variation in parameter combinations using the total wood mass (stem, branches, and roots) of four-year-old Populus deltoides as an index.Each two parameters are combined as a pair for sensitivity analysis.Numbers in bold in the matrix indicate the total wood mass change with 25% increases in both tested parameters; the numbers in italics in the matrix indicate the total wood mass change with 25% reductions in both tested parameters.The value in the cell is the increase or reduction rate in percentage compared to the total wood mass simulated using original parameter values from Table 2. Symbols: J max25 = maximum rate of electron transport at 25°C; V cmax25 = maximum rate of carboxylation at 25°C; R d25 = day-time leaf respiration rate at 25°C; and SLA = specific leaf area.
25% + + 25% -- effect on output (J max25 , V cmax25 , R d25 , and SLA) and a parameter with weak effect (ρ w ) were tested.Each of the five parameters included was tested in conjunction with another parameter, yielding four possible scenarios of combined variation (Tables 6 and 7).The single parameter sensitivity analyses revealed that changes in model parameters resulted in parallel changes in model output such as height and biomass.Thus, the combined sensitivity analyses were performed using only the total woody biomass as a representative output variable to test the effects of combined parameter variation.
The combined test for sensitive parameters indicated that changes in two parameters in the same direction could strongly influence tree growth; e.g.increase in both J max25 and V cmax25 enhanced biomass production by more than the additive effect of single parameter changes (Tables 5 and 6).The effect of combined reduction in these parameter values had a smaller effect on biomass than the additive effect of single parameter changes, though more than the effect of changes in either parameter alone.A change in R d25 together with J max25 significantly reduced the positive effect of the latter on growth (Table 6).The most significant parameter pair was the combination of increases in both J max25 and SLA, which led to a 151% increase in tree biomass (Table 6); also the reduction of these parameters together resulted in the largest biomass reduction.
Changes of two parameters in the opposite direction might cancel each other's effect on tree growth; this was most obvious between J max25 and SLA (Table 7).However, in most cases one of the parameters in the pair was dominant and opposite change in the other parameter only reduced the effect of the dominant one; an example of such a pair is J max25 and V cmax25 (Table 7).

New Features in the LIGNUM Model
Since the first publication of the functionalstructural tree model LIGNUM for Scots pine (Perttunen et al. 1996), it has been under constant development (Perttunen et al. 1998, Perttunen andSievänen 2005) and it has been adopted to new tree species (Lo et al. 2001, Perttunen et al. 2001).These developments have shown that the basic modular structure of LIGNUM is flexible and it provides a solid basis for modeling functional-structural development of both conifers (Perttunen et al. 1996, 1998, Lo et al. 2001) and deciduous trees (Perttunen et al. 2001).However, only the recent inclusion of a new submodel for structural development based on L-systems (Perttunen and Sievänen 2005) and introduction of the voxel space radiation interception submodel (Sievänen et al. 2008) have added fundamentally new features to LIGNUM beyond reparameterization for a new species (Lo et al. 2001, Perttunen et al. 2001).In this contribution, we have described a new version of the LIGNUM model that includes several basic developments in addition to adaptation to P. deltoides.These developments include application of a new leaf CO 2 exchange model, use of real weather data, the system of nested time steps, inclusion of the first approximation of fine root dynamics and rhizodeposition, and formulation of the voxel space submodel for computing the interception of photon flux in broad-leaved trees.
In previous LIGNUM studies, a linear relationship with intercepted radiation was used to predict photosynthetic production of a tree in order to simplify the modeling process (Perttunen et al. 1996).There are also tree growth simulation models using photosynthesis-Q response curves to calculate the tree photosynthesis more precisely (Rauscher et al. 1990).A biochemically-derived photosynthesis model (Farquhar et al. 1980, von Caemmerer and Farquhar 1981, Long 1991, De Pury and Farquhar 1997, Medlyn et al. 2002), now commonly called Farquhar's model, was applied in this study to reflect the photosynthetic response to varying photon flux density affected by internal crown shading.Farquhar's model yielded a non-linear relationship between photon flux density and net CO 2 assimilation rate resulting in a stronger physiological basis for modeling primary production of a tree.
Photosynthetic production was calculated every 30 minutes with concurrent absorbed photon flux density, temperature, and solar position.The 30-minute-photosynthesis calculation with real weather data helped to connect the model and will vary from case to case even with the same set of parameters.However, the average result will reflect the actual growth situation (Figs. 17 and 18, Table 3).
Tree roots were considered only in terms of root biomass in the earlier LIGNUM versions (Perttunen et al. 1996, 1998, 2001, Perttunen and Sievänen 2005).However, coarse and fine roots are functionally different: the former are part of the transport system of a tree while the latter are metabolically active and forage the soil for water and nutrients.It has been shown that fine roots are highly dynamic in nature (Block et al. 2006) and that rhizodeposition of C may form an important proportion of C allocation within a tree (Högberg and Read 2006).We did not have any direct measurements on these factors but the fine root dynamics and rhizodeposition were modeled according to a literature review (Horwath et al. 1994, Coleman et al. 2000, Högberg et al. 2002, Coyle and Coleman 2005, Block et al. 2006, Högberg and Read 2006).While this approximative module helped to achieve a realistic C balance of young P. deltoides, it will require further study for proper parametrization, especially because the model seems to be sensitive to the root to foliage C allocation ratio (Table 5).
These new developments were integrated to the original LIGNUM model.The voxel model was integrated to Firmament submodel on incident radiation (Appendix 2; Perttunen et al. 1998), vigor index (Nikinmaa et al. 2003) was applied for modeling apical dominance, tree structure was derived using the L-systems-LIGNUM combination (Perttunen and Sievänen 2005), and tree structure was based on the elementary units developed for sugar maple (Perttunen et al. 2001) with slight modification.These existing LIGNUM submodels were re-parametrized for P. deltoides based on field data and literature review.

Comparison of Simulation Results with Field Data
Modeling results for tree height were close to field measurements (Fig. 18).The one-year result for poplar growth from the ECOPHYS model is 1.10 m (Rauscher et al. 1990), which is close to the LIGNUM simulation results (1.06 m; Fig. 18) in local environment so as to better reflect the effect of environment on tree growth.Photosynthetic light reactions respond within seconds to changes in photon flux density, which also varies rapidly on partially cloudy days.Thus, even a 30-min time step may be considered long but it provides a reasonable compromise between the number of calculations and precision.During a 6 week structure update period, the 30 min environment simulation is expected to give a realistic variation in cloudiness.The 30-min-step is short enough to fully cover the changes in radiation environment because of the changes in solar position during a day.In earlier LIGNUM versions, the environment was described only by means of the annual integral of global radiation assuming standard overcast distribution throughout a year (Perttunen et al. 1996(Perttunen et al. , 1998) ) and by masking sectors of the sky for simulating the effect of forest gap (Perttunen et al. 2001).Incorporation of real weather data is expected to increase realism in the description of the growing environment of a tree in LIGNUM.The use of an annual time step was a reasonable and obvious choice for earlier LIGNUM versions dealing with relatively slowly-growing tree species, given that tree growth occurs via a process cycle repeated every year (Perttunen et al. 1996(Perttunen et al. , 1998(Perttunen et al. , 2001)).A short time step simulation has been applied in juvenile poplar growth in the ECOPHYS model (Rauscher et al. 1990) to provide a more detailed update of tree structure.Short time step simulation was also considered necessary in simulating the fast-growing P. deltoides.The application of two nested short time steps enabled the model to update the tree canopy structure more frequently and adjust accordingly photon flux interception in the model.The length of the time step for structure update was intermediate to that of the photosynthesis and the seasonal time steps.
The Monte-Carlo voxel space model, which is applied in this LIGNUM version, is stochastic in nature.The chance of radiation interception is related to the leaf area density in the voxel.The stochastic Monte-Carlo model accommodates the inherent variety in the growing tree in nature.Every tree is unique, and its growth will be affected by both intrinsic and environmental factors.The results of tree growth simulations this study for one-year-old P. deltoides grown in an open site.Compared to the detailed calculations of geometrical shape and size of each leaf required in ECOPHYS, the use of voxel space in LIGNUM simulation greatly simplified the modeling process.Furthermore, such modeling simplifications enable the LIGNUM model to simulate and predict growth for older trees than more detailed models.
Populus deltoides is known to be a fast-growing and shade-intolerant species (Cooper 2002).The system of nested time steps for accommodating fast growth and modeling efficiency appeared to produce realistic results both in quantitative terms in height and biomass growth (Figs.18 and 19, Table 3) and visually by generating trees (Fig. 20) that closely resembled the trees in the open study site (Fig. 11).Visualization results also showed that application of the vigor index (Nikinmaa et al. 2003) together with the Gravelius order of an axis created realistic apical dominance.Branching pattern was also realistically replicated in visualizations (Fig. 20).
Visualization also showed that shade intolerance of P. deltoides was adequately described in the LIGNUM model (Fig. 21).The sensitivity analysis revealed that LIGNUM was sensitive to the maximum electron transport rate, J max25 , used in the photosynthesis submodel (Table 5).This feature seems realistic given the shade intolerance of P. deltoides: sensitivity to J max25 indicates that even slight shading may substantially reduce photosynthetic production.
The parameter set for fine root dynamics and rhizodeposition was partly estimated from published work and partly fit to yield an output comparable with field observations.The model seems to be sensitive for two of the parameters of root dynamics, the root to foliage C allocation ratio and the ratio of rhizodeposition to total C allocated to roots (Table 5).Thus, these parameters will require further study.Although the parameters of the CO 2 exchange model were determined in a growth chamber using greenhouse-grown seedlings, they may be considered to be quite accurate.Further, the attenuation of solar radiation within the tree canopy may be assumed to be computed fairly accurately as the model is based on solid theory of radiation penetration into a plant canopy.Thus, C allocation to roots may be at the moment the weakest point of the model, although better justified than the earlier constant proportion of tree biomass.Further, the model does not consider the effects of drought on leaf CO 2 exchange.Dry periods are common in the field site in mid-Missouri in late summer and, thus, the exclusion of drought effects may explain the high estimate of tree biomass without shading (Table 3).
Even with a large amount of field measurements, the tree structure and physiological growth rules were simplified for the sake of modeling efficiency.The branches were assumed to be evenly distributed in all azimuth directions.Carbon allocation for the shoot applied empirical rules in segment elongation and thickening such as vigor index and radiation index.Aboveground environment of a tree was described as a voxel space.In spite of these simplifications, exclusion of the potential drought effects, and the provisional literature-based C allocation submodel to roots and rhizodeposition, LIGNUM model adapted to P. deltoides produced satisfactory convergence between simulation results and field data as indicated above.

Applications
In this study, the improved LIGNUM model was applied for simulating growth of P. deltoides in short-rotation applications in flood plain areas of mid-Missouri.The results on individual tree growth in an open site are especially useful for evaluating the suitability of P. deltoides for agroforestry designs with widely-spaced trees like the one used for measuring model parametrization and validation data in this study (Fig. 11).This kind of relatively open tree spacing is common in North American agroforestry (Garrett and Harper 1999, Gillespie et al. 2000, Miller and Pallardy 2001, Udawatta et al. 2005).Simulation models provide a research tool for complementing -but not replacing -traditional field experiments in forestry and agroforestry, which are more time consuming, labor intensive, and susceptible to adverse environmental events.
One of the purposes of tree growth modeling in forest research and management is to predict future tree growth and to answer questions regard-ing the potential effects of certain management practices on a forest over time (Pretzsch et al. 2008).Tree growth is a long-term process, results of which cannot be observed instantly.Simulation models may be used to predict the effects of forest management on future growth.Functionalstructural tree models like LIGNUM may also be used as a researchers' "workbench" for examining different physiological and ecological hypotheses on tree growth and yield when designing field and controlled-environment experimentation (Sievänen et al. 2000).
Weather factors were incorporated into the LIGNUM model in this study.Thus, LIGNUM may be used for predicting possible effects of global climate change on trees, although partial re-parameterization may be needed even for P. deltoides for simulating the physiological effects of elevated atmospheric CO 2 concentration and temperature acclimation (Turnbull et al. 2002, Bernacchi et al. 2003, Davey et al. 2006).
Because the current version of LIGNUM model was adapted for simulating the behavior of P. deltoides in a short-rotation production system, it was assumed that the heartwood formation in the juvenile P. deltoides was negligible (Bruce Cutter, University of Missouri-Columbia, pers.comm.).It is necessary to include heartwood formation in modeling the growth of older P. deltoides trees.Further, LIGNUM does not include reproductive allocation and other physiological processes that affect development of mature trees.Thus, the LIGNUM model adapted to P. deltoides must be used only for its original purpose, i.e. studying short-rotation systems.
Although the simulation is intended to facilitate tree growth research by easing the burden of field experiments, extensive field work will always be necessary to understand tree growth patterns.The morphological analysis requires data on crown structure including branch bifurcation, branch direction, branch angle, segment length, and leaf size.Information about branch diameters and biomass distribution throughout the tree is needed for modeling C allocation.

Concluding Remarks
Integrated with the original LIGNUM model framework, the version adjusted for P. deltoides included new features in the modeling environment.These were the biochemically-derived photosynthesis submodel; nested time steps for allowing reasonable resolution for simulating physiological processes, structural development, and annual biomass production; incorporation of field-measured weather data for modeling the response of physiological processes to the environmental variation and calibrating simulation results by field observation; and application of the Monte-Carlo voxel space submodel for both simulating the stochasticity of tree growth and improving computational efficiency.The new simulation applied the specific parameter system for P. deltoides growing in mid-Missouri flood plain environment.Besides the empirical branching parameters, CO 2 exchange parameters, and shoot growth rules, the model included speciesspecific growth patterns, such as radiation interception assumption, C allocation principles, and literature-obtained parameters.
The LIGNUM model has proven to be a general one for tree growth simulation for both conifer and deciduous species, whether slow-growing or fastgrowing in nature.The simulated height and biomass growth of P. deltoides adequately matched with field observations although the simulations tended to produce trees that resembled the biggest trees in the field site.Visual representation of the simulated trees facilitated interpretation of the simulation results as well as encouraged user confidence in the modeling approach.The simulated response of tree growth to variations in photon flux input was logical.The LIGNUM model version presented may be used as an aid for predicting tree growth and yield, economic analyses, and management decision making concerning P. deltoides in short-rotation systems.However, it cannot be applied for analyzing growth and yield of mature P. deltoides because heartwood formation and physiological processes affecting mature tree growth are not included in the current version.Improving the precision of root functions and including root architecture together with simulation of drought effects are obvious needs for further model development.

Appendix 2. Computation of incident radiation.
The distribution of photon flux in the whole sky is set in LIGNUM by means of the Firmament submodel (Perttunen et al. 1998(Perttunen et al. , 2001)).In the LIGNUM version described in this contribution, Firmament was modified to handle direct and diffuse fluxes separately.The upper hemisphere is divided into several sectors in both horizontal and vertical directions.The division is implemented so as to have as equal solid angle as possible.The width of sectors in the same azimuth is the same, and there are different numbers of sectors at each inclination level.Both the number of inclinations and the number of azimuth can be set in the model.Let N φ denote the number of inclination divisions (vertical) and N Θ the number of azimuth divisions (horizontal).The number of the azimuths defined in the model is the mean number of sectors contained within the inclination zones.There is also a zenith sector at the top of the hemisphere pointing directly upward.The total number of sectors in the hemisphere is equal to N φ × N Θ + 1.The total area of the upper hemisphere is 2π; thus the area of the zenith sector is equal to 2π / (N φ × N Θ + 1).The width of the inclination zone is equal to (π / 2 -angle of zenith sector) / N φ .
The diffuse photon flux density from each sky sector is distributed according to the zonal brightness of the standard overcast (Ross 1981): where φ is the elevation of the sector positioned in the sky hemisphere, and d φ is the fraction of diffuse flux from sector φ out of the total diffuse flux received by the sphere.The photon flux density received by each sector increases with the increase of its altitude.The equation of d φ is scaled so that the flux on the horizontal plane is equal to 1: Both direct and diffuse components of Q reach the tree through the sky sectors.Diffuse flux occurs in every sector, and the direct flux in one sector.The photon flux available for photosynthesis at a particular point is the sum of the photons intercepted through all sky sectors.
The value of Q measured in the field is the sum of direct and diffuse components.The ratio (Φ b ) of the direct component of photon flux to the measured total photosynthetic photon flux density was estimated (Weiss andNorman 1985, Nygren et al. 1996): where R d is daytime respiration rate of leaves (µmol m -2 s -1 ), i.e. respiration related to metabolic processes other than photorespiration; Γ * /C is the ratio of the rates of photorespiration and carboxylation; C i is the intercellular CO 2 concentration in a leaf (µmol mol -1 ), corrected for solubility from ambient CO 2 concentration, C a , at temperature T relative to 25°C (Long 1991): V c is the carboxylation rate (µmol m -2 s -1 ); Γ * is the CO 2 compensation point of photosynthesis in the absence of daytime respiration (µmol mol -1 ).The CO 2 compensation point is estimated: where K c is the Michaelis constant for CO 2 (µmol mol -1 ); V omax is the maximum ribulose bisphosphate (RuP 2 ) saturated rate of oxygenation (µmol m -2 s -1 ); K o is the Michaelis constant for O 2 (mmol mol -1 ); O i is the intercellular concentration of O 2 in the leaf (mmol mol -1 ), corrected for solubility from ambient O 2 concentration, O a , at temperature T relative to 25°C (Long 1991): V cmax is the maximum RuP 2 -saturated rate of carboxylation (µmol m -2 s -1 ) at temperature T corrected relative to the rate at 25°C (Yin et al. 2004): where V cmax25 is the V cmax value at temperature 25°C and R is gas constant (8.314J mol -1 K -1 ).The relationship between V omax and V cmax is (Long 1991): Farquhar et al. (1980) predicted an abrupt change in the dependence of photosynthetic rate on intercellular pressure of CO 2 as the limitation on photosynthesis changed from RuP 2 -saturated carboxylation rate to RuP 2 regeneration limited rate.The RuP 2 -saturated rate of carboxylation, W c (µmol m -2 s -1 ) is: The RuP 2 -limited rate of carboxylation, W j (µmol m -2 s -1 ) is: where J is the potential rate of electron transport (µmol m -2 s -1 ).It depends on absorbed photon flux (Q; µmol m -2 s -1 ): where J max is the light-saturated potential rate of electron transport (µmol m -2 s -1 ).An optimum response of J max to temperature is described (De Pury and Farquhar 1997, Medlyn et al. 2002, Yin et al. 2004): where S J is the entropy term (650 J K -1 mol -1 ) and D J is the energy of deactivation (200,000 J mol -1 ).Thus, carboxylation rate V c in eq.(A3-1) is: where R d25 is the respiration rate at reference temperature 25°C; R d is the respiration rate adjusted to measured ambient temperature T; Q 10 is the temperature coefficient (Q 10 = 1.78;Turnbull et al. 2002), which is defined as the change of respiration rate with temperature increasing by 10°C.
For the diameter growth, let's assume that unit area of sapwood can support foliage mass defined by parameter Y c .Thus, given the leaf area LA s in the new segment, the sapwood requirement SA ws of the new segment becomes: where SLA ws is the specific leaf area.
For the older segments, we assume that a target segment below a branching point must match the sapwood areas of the attached segments immediately above, connected to the same branching point.More precisely, denote the cross sectional area of the target segment with SA wt and the sum of the sapwood areas of the attached segments SA wu .Then the new cross sectional area SA wtnew of the target segment is: We assume there is no sapwood senescence in young P. deltoides trees up to 8-year-old.The equations for SA ws and SA wtnew define the pipe model for P. deltoides and prevents segments from shrinking if the sapwood required, SA wu , is already supported.

Fig. 1 .
Fig. 1.Top: Populus deltoides growing in sand deposits left by the 1993 centennial floods of the Missouri River in Plow Boy Bend, Moniteau County, Missouri, USA, in May 2001.Bottom: Populus deltoides growing in a Missouri River flood plain area frequently flooded in the springs; Franklin Island, Howard County, Missouri, May 2001.
physiological time step takes place within the structural development cycle.Physiological processes, including interception of photon flux, leaf photosynthesis, and leaf respiration are affected not only by tree structure but also by the radiation environment and air temperature.Real weather data are utilized in LIGNUM simulation.Photon flux interception and leaf photosynthesis are calcu-lated in the model with 30-minute-interval under the variable photosynthetic photon flux density (Q) affected by solar position and cloudiness.The interception of photon flux is calculated applying a voxel space model.Photosynthesis and leaf respiration are calculated in terms of intercepted Q, air temperature, and leaf area (Appendix 3).Short-time calculation of photon flux interception and photosynthesis enables simulation of the nonlinearity of the photosynthetic response to Q.The length of the physiological time step was selected

Fig. 3 .
Fig. 3. General model structure for individual tree growth in LIGNUM for P. deltoides.The inner loop (dotted line) is the physiological time step for interception of photon flux and photosynthesis repeated every 30 minutes.The intermediate loop (dashed line) is the structure update cycle with a longer interval, called structural update time step or medium time step.Structural update process includes an expansion step and a split step.Axillary buds from segment split are dormant and would be reactivated at the beginning of next growing season.The outermost loop (solid line) is the growing season cycle.

Fig. 4 .
Fig. 4. The voxel space divides the growing space of the tree into cubic voxels.Tree compartments including leaves and segments are added into each voxel.The more tree compartments a voxel contains, the greater the possibility that the photon flux will be intercepted in the voxel.(Figure courtesy of Mika Lehtonen).

Fig. 5 .
Fig. 5. Interception of photon flux in voxel space.Left: A general voxel with a tree segment and leaves.A photon flux entering the voxel from a sky region may pass through the voxel (thick line) or be intercepted by tree compartments (thick line turning thin).Right: Projection of the voxel (polygon) and big leaf area (shaded area in polygon) to the plane perpendicular to the photon flux from the sky region.

Fig. 6 .
Fig. 6.Four scenarios for interception of an incoming photon flux: (1) penetration to the focal voxel without contact with foliage in shading voxels; (2) transmission through the foliage of a shading voxel, which has one layer of leaves; (3) absorption by the foliage of one shading voxel, which has more than one layer of leaves; or (4) transmission through the foliage of one shading voxel and interception by another shading voxel.
Turn(α): Turn by rotating around  U by angle α.Pitch(α): Pitch by rotating around  L by angle α.Roll(α): Roll by rotating around  H by angle α.

Fig. 7 .
Fig. 7. Controlling the orientation of the turtle with three unit vectors  H ,  L and  U.The sequence of modules

Fig. 8 .
Fig. 8. Schematic presentation in 2D of the two first structural updates defined by the L-system in Appendix 5.For simplicity, assume that after the elongation phase the allocation of photosynthates results segment lengths between 2 and 4 cm in the main axis, and between 1 and 2 cm in the branches.S: The initial 10-cm-tall tree with an active apical and two active lateral buds (the rule Start in Appendix 5).1-1E: The 1 st elongation phase, the generation of 3 new 30-cm-segments (the rule F(s,d) < B(d)).1-1: The 1 st structure update of the 1 st year complete.The main axis segment is split into 3 and the two segments in the branches are split into 2 segments each (The 2 nd and the 3 rd rewrite in the rule F(s,o) > Split()).The new lateral buds (marked grey) remain dormant during the year they are created.1-2E: The second elongation phase.1-2: The 2 nd structure update of the 1 st year complete.

Fig. 10 .
Fig. 10.Radiation index curve for Populus deltoides as a function of relative photosynthetic photon flux density (Q rel ).The Q rel value determines the radiation index in segment elongation in the carbon allocation submodel.The Q rel threshold point is 0.036, under which radiation index is equal to 0 (radiation compensation point of photosynthesis).The curve is borrowed from net CO 2 assimilation curve at 25°C in Fig.16by setting radiation index value 1 to the approximate saturation point of photosynthesis (460 µmol m -2 s -1 ).

Fig. 11 .
Fig. 11.The field site used for measuring parameterization and validation data for LIGNUM model adapted for Populus deltoides in August 2003.The experiment was situated in the Horticulture and Agroforestry Research Center of the University of Missouri-Columbia in Howard County, Missouri, USA.

Fig. 12 .
Fig. 12. Measurements of branch structure.H: height of insertion point; L: segment (internode) length; D0: stem diameter below a branching point; D1: stem diameter above branching point; D2: diameter of main branch; d0: diameter of main branch below a branching point; d1: diameter of main branch above a branching point; d2: diameter of subbranch; α: inclination of main branch; and β: branch bending angle.

Fig. 13 .
Fig. 13.Azimuth distribution of main branches in Populus deltoides growing in an open site.The strength in each direction is the frequency of branch azimuths.

Fig. 14 .
Fig. 14.Inclination distribution of main branches in Populus deltoides growing in an open site.Inclination was measured relative to horizon.

Fig. 15 .
Fig. 15.Branch growth dynamics of Populus deltoides in an open site: a) evolution of the average leaf number per branch during a growing season; b) the average elongation of new branches during a growing season.Error bars denote standard deviation.

Fig. 16 .
Fig. 16.Net CO 2 assimilation rate in mature leaves of Populus deltoides as a function of photosynthetic photon flux density.The net CO 2 assimilation rate was computed using the model presented in Appendix 3 with parameters given in Table2.

Fig. 18 .
Fig. 18.Simulated changes in height and leaf area of Populus deltoides during eight growing seasons simulated using LIGNUM model.Each point is the mean of five runs of the LIGNUM model.

Fig. 19 .
Fig. 19.Simulation of Populus deltoides during two years with LIGNUM using L-system defined in Appendix 5 with metabolism included.Parameters have nominal values.S: The initial tree is as in Fig. 9. 1-1: first structure update of the first year.2-4: fourth structure update of the second year.The architectural development of the tree is now affected by the local radiation regime in the tree crown and the resource allocation in different parts of the tree (cf.Fig. 9).

Fig. 20 .Fig. 21 .
Fig. 20.Visualization of a simulated 4-year-old and 8-year-old Populus deltoides growing in a central Missouri environment without shading by other trees.Fig. 21.Visualization of growth response of Populus deltoides to shading: a simulated 4-year-old tree under 50% and 75% of measured photon flux density.

Fig. 22 .
Fig. 22. Increase of the number of leaves (N l ), voxels (N v ), and voxels with leaves (N lv ) as a function of tree age (a) in Populus deltoides simulations with LIGNUM.The leaf number trendline is N l = 40.5×a 3.1882 (r 2 = 0.98); voxel number trendline is N v = 42.6×a 2.2302 (r 2 = 0.97); the trendline for the number of voxels with leaves is N lv = 2.00×a 2.8693 (r 2 = 0.99).
x) is the area of the horizontal plane surface of the unit sphere, and d(x)sin(x) is the photon flux density on horizontal plane surface.When diffuse photon flux is set with Q d in the Firmament submodel, the diffuse brightness for each sector is Q d × d φ according to the sector position.
r b (PAR) and r g (PAR) are potential direct and global radiation of photosynthetically active wavebands; r obs is measured global short-wave radiation, and r g (SWR) is potential global short-wave radiation.The potential direct and global photosynthetically active radiations were calculated(Nygren et al. 1996): is subtracted from photosynthetic rate to yield net C uptake.The leaf respiration rate is adjusted by temperature with a temperature coefficient: Thermic growing season ranges from April to October but trees are leafless most of April and in this work we will consider the period from May through October as the growing season.Most precipitation falls in summer from May to August.Summers are hot and sunny.Water deficits typically develop in late summer because of high net radiation combined with high temperatures and relatively low rainfall.Midwinter temperatures vary around 0°C and freezing temperatures may occur from October to April.Selected monthly weather data for years 2002-2005 are shown in Table 1.

Table 1 .
Selected weather data for the field site used for measuring model parameterization and validation data; University of Missouri's Horticulture and Agroforestry Research Center in New Franklin, MO, USA (92°46'W; 39°01'N; 197 m altitude).

Table 3 .
Comparison between simulations with LIGNUM model and field data for height growth and wood biomass (above and below ground) of Populus deltoides with no shading.Parameters have the basic values given in Table2.Measured data are mean ± standard deviation for five trees, simulated data are mean ± standard deviation for five model runs.
Fig. 17.Biomass development of Populus deltoides during eight growing seasons simulated using LIGNUM model.Each point is the mean of five runs of the LIGNUM model.

Table 4 .
Comparison between simulations with LIGNUM model and field data for height growth and wood biomass (above and below ground) of Populus deltoides with 25% shading.Parameters have the basic values given in Table2.Measured data are mean ± standard deviation for five trees; simulated data are for five model runs.