Appendix B. The Markov chain Synthetic Turbulence Generation (STG) model. A pdf file of this appendix is also available.
The Markov chain STG model is a version of the forest seed dispersal model of Nathan et al. (2002b), adapted for grassland ecosystems. The general structure is presented in the Methods section and by Katul and Albertson (1998) and Katul and Chang (1999). In this appendix the estimation of the Eulerian flow statistics needed to drive the Lagrangian dispersion is described, followed by descriptions of the computational grid, numerical scheme, boundary conditions and closure constants.
Upon temporal and spatial averaging the conservation of momentum equations, and following the closure approximations of canopy flows, the Wilson and Shaw (1977) closure model simplifies to the following set of ordinary differential equations (ODEs):
Mean Momentum:  _{} 
Tangential Stress Budget:  _{} 
Longitudinal Velocity Variance:  _{} 
Lateral Velocity Variance:  _{} 
Vertical Velocity Variance:  _{} 
Where _{} is a characteristic velocity scale, C_{d} is a drag coefficient, a(z) is a leaf area density, _{j} = a_{j} × L_{ws} , with L_{ws} a characteristic length scale specified using the formulation given by Katul and Albertson (1998) and not permitted to increase at a rate larger than k, and a_{1}, a_{2}, a_{3} and C are determined so that the flow conditions well above the canopy reproduce wellestablished surface layer similarity relations. With estimates of the five constants (a_{1}, a_{2}, a_{3}, C, and ), the above five ODEs can be solved iteratively for the five flow variables _{}, _{}, _{}, _{}, and _{}, which are used to drive the Lagrangian model.
The computational grid
The computational flow domain was set from zero to 20 × h. The grid node spacing is z = 0.005 m. This grid density was necessary due to rapid variability in leaf area density close to the canopy top. Parameter values at the exact location of the seed are calculated by interpolation between the grid nodes, or extrapolation in the case of x_{3} > 20 × h. To ensure seeds do not exit the atmospheric boundary layer during the computation of dispersal trajectories, the vertical position of a seed is not allowed to increase above 800 m.
The numerical scheme
The five ODEs for the Wilson and Shaw (1977) model were first discretized by central differencing all derivatives. An implicit numerical scheme was constructed for each ODE with boundary conditions to be discussed in the following section. The tridiagonal system resulting from the implicit forms of each discretized equation was solved using the Tridag routine from Press et al. (1992; pp.42–43) to produce the turbulent statistic profile. Profiles for all variables were initially assumed, and a variant of the relaxation scheme described by Wilson (1988) was used for all computed variables. Relaxation factors as small as 5% were necessary in the iterative scheme because of the irregularity in the leaf area density profile. The measured leaf area density was interpolated at the computational grid nodes by a cubicspline discussed in Press et al. (1992; pp.107–111) to insure finite second derivatives of a(z). Convergence is achieved when the maximum difference between two successive iterations in _{} did not exceed 0.0001%. We checked that all solutions were independent of z (as described in Katul and Albertson 1998). Calculation of dt is described in Appendix A.
Boundary conditions and closure constants
Typically, the wellestablished flow statistics in the atmospheric surface layer provide convenient upperboundary conditions for closure models. The boundary conditions used are:

Where _{} is the standard deviation of any flow variable _{} (= _{}), A_{u} = 2.2, A_{v } = 2.0, and A_{w} = 1.4 (Panofsky and Dutton 1984).
The closure constants are dependent on the choice of the boundary conditions and are determined by assuming that in the atmospheric surface layer (z > 2h), the fluxtransport term is negligible and that _{}, _{}, and _{} become independent of z for nearneutral conditions. These simplifications result, after some algebraic manipulations (e.g., see Katul and Albertson 1998; Katul and Chang 1999), in the following relationships between A_{u}, A_{v}, and A_{w} and a_{2}, a_{3}, and C:

Where A_{q} = (A_{u}^{2} + A_{v}^{2} + A_{w}^{2})^{1/2}. The closure constant a_{1} is determined by noting that the eddydiffusivity is k × (z – d) × u_{*} in the surface layer. Hence, q_{1} becomes identical to k × (z – d) × u_{*} leading to a_{1} = 1/A_{q}. The above equations are the first analytic expressions relating closure constants to ASL boundary conditions for the Wilson and Shaw (1977) model as described by Katul and Albertson (1998) and Katul and Chang (1999). Table B1 summarizes the closure constants used resulting from our choice of A_{u}, A_{v}, and A_{w}.
Table B1. Closure constants used in the Markov chain STG model for A_{u } = 2.2, A_{v} = 2.0, and A_{w} = 1.4.
Closure constant
(Wilson and Shaw 1977)Value
A_{1}
0.30
A_{2}
1.58
A_{3}
20.8
0.07
C
0.12
C_{d}
0.20
LITERATURE CITED
Katul, G. G., and J. D. Albertson. 1998. An investigation of higher order closure models for a forested canopy. BoundaryLayer Meteorology 89:47–74.
Katul, G.G., and W.H. Chang. 1999. Principal length scales in secondorder closure models for canopy turbulence. Journal of Applied Meteorology 38:1631–1643.
Nathan, R., G.G. Katul, H.S. Horn, S.M. Thomas, R. Oren, R. Avissar, S.W. Pacala, and S.A. Levin. 2002b. Mechanisms of longdistance dispersal of seeds by wind. Nature 418:409–413.
Panofsky, H. and J. Dutton. 1984. Atmospheric Turbulence: Models and Methods for Engineering Applications. John Wiley, New York.
Press, W.H., S.A. Teukolsky, W.T. Vetterling, and B.P. Flannery. 1992. Numerical Recipes in Fortran. Cambridge University Press.
Wilson, J.D. 1988. A second order closure model for flow through vegetation. BoundaryLayer Meteorology 42:371–392.
Wilson, N.R., and R.H. Shaw. 1977. A higher order closure model for canopy flow. Journal of Applied Meteorology 16:1198–1205.