## Abstract

Tall vegetation and complex terrain create difficult conditions for measuring and modeling net ecosystem–atmosphere exchanges of carbon, water vapor, and pollutants. The instability of canopy flow regimes over complex terrain is critical for understanding what factors are essential to control exchanges between different canopy flow regimes. In this paper, an analytical criterion of instability of the terrain-induced canopy flows is derived from the simplified thermal-hydromechanical equations by nonlinear instability analysis. The stability of the terrain-induced canopy flows and an oscillation solution are predicted based on the instability criterion. It is found that the critical values of control parameters are determined by the terrain slope, drag coefficient, and leaf area density of vegetation.

## 1. Introduction

Terrain-induced (TI) flows often result in complicated situations in measurements of trace gas exchanges between land surface and atmosphere. For example, the eddy covariance technique (micrometeorological approach) has become the most important method for measuring the land–atmosphere exchanges during the past decades (Baldocchi et al. 2001; Baldocchi 2003). The eddy covariance methodology was derived under assumptions of a perfectly flat topography and homogeneous land cover (Goulden et al. 1996; Lee 1998; Massman and Lee 2002). However, most eddy flux tower sites are located in complex terrain with heterogeneous vegetation and lie outside of the idealized condition (Baldocchi 2008; Finnigan 2008). Many research groups have reported mismatches between nocturnal eddy flux measurements and soil/stem chamber measurements at their sites (Norman et al. 1997; Law et al. 1999; Janssens et al. 2001; Barford et al. 2001; Curtis et al. 2002; Cook et al. 2004; Bolstad et al. 2004; Miller et al. 2004; Launiainen et al. 2005; Van Gorsel et al. 2007; Myklebust et al. 2008). It is commonly agreed that the mismatch between eddy flux measurements and chamber measurements is caused mainly by thermal drainage flows on slopes (Massman 1987; Yi et al. 2000; Massman and Lee 2002; Aubinet et al. 2003; 2005; Staebler and Fitzjarrald 2004; Feigenwinter et al. 2004, 2008; Miller et al. 2004; Sun et al. 2007; Aubinet 2008; Finnigan 2008). During calm nights, canopy flows are separated by a superstable layer that is defined by an infinite value of the Richardson number due to the minimum wind speed and the significant potential temperature gradient that are near the level of maximum leaf density (Yi et al. 2005; Yi 2008). Therefore, the CO_{2} respired from soils and advected by cold drainage flows below the superstable layer cannot be properly accounted for by the sensors on the tower above the canopy. Sun et al. (2007) reported that the warm valley breeze also caused significant horizontal advection errors in eddy flux measurements, particularly in transition time periods between day and night. These observations have demonstrated that the land–atmosphere exchanges depend on the regimes of TI canopy flows. In this paper, we do not intend to solve the advection error problem for the eddy flux measurements but rather study only an important aspect of the advection issues: the instability of TI canopy flows.

TI flows are complicated because they result from multiscale interactions between irregular terrain and the overlaying atmospheric boundary layer (Mahrt and Larsen 1982; Mahrt et al. 2001; Finnigan and Belcher 2004; Ross and Vosper 2005; Katul et al. 2006; Poggi and Katul 2007a,b,c; Ross 2008). The small scale (∼10 km; several hills) of TI flows close to the land surface plays a more important role in the land–atmosphere exchanges than does the large scale (∼100 km; mountain wave) of TI flows aloft at higher levels. Particularly, the very small scales (∼1 km; a single hill) of TI flows within and immediately above canopies directly affect the land–atmosphere exchanges (Massman 1997; Froelich et al. 2005; Pypker et al. 2007). These TI canopy flows depend strongly on the local slope, canopy structure, and thermally driven forcing (Yi et al. 2005; Froelich and Schmid 2007; Sun et al. 2007; Heinesch et al. 2007, 2008; Prabha et al. 2008; Mao et al. 2008; Schaeffer et al. 2008; Myklebust et al. 2008). They are also modulated by the relatively large scale of TI flows (Turnipseed et al. 2003, 2004). The basic characteristics of these TI flows can be classified by steady states (e.g., katabatic and anabatic flows) and persistent oscillations (Doran and Horst 1981). Different steady states and persistent oscillations of these TI flows affect the land–atmosphere exchanges differently. The stability conditions of these steady states and oscillations are of considerable importance to understanding the switch mechanisms of land–atmosphere exchanges from one type of TI canopy flows to another.

The purpose of this paper is to derive the stability conditions of TI canopy flows. This is a complex topic. Although TI flows have been extensively studied in the past several decades (Fleagle 1950; Prandtl 1952; Manins and Sawford 1979; McNider and Pielke 1981; Mahrt and Larsen 1982; Chimonas 1993; Atkinson and Shahub 1994; Whiteman et al. 2004), concerns of air pollution and the impact of TI canopy flows on land–atmosphere exchanges have called for attention in recent years (Sun et al. 2007; Pypker et al. 2007; Hammerle et al. 2007; Froelich and Schmid 2007; Wang et al. 2005; Turnipseed et al. 2002, 2003; Mao et al. 2007). The instability analysis of TI canopy flows has not been reported except for some canopy wave analyses (Lee 1997; Pulido and Chimonas 2001). In this paper, we attempt to derive the analytical criterion of instability conditions of TI canopy flows based on a simplified model developed by Fleagle (1950) and further modified by McNider (1982). The governing equations of velocity and potential temperature in the simplified Fleage–McNider model were used and averaged for the along slope layer that is without vegetation. The present study takes into account the interactions between vegetation and canopy flows in the layer-averaged governing equations. The simplification of the layer-averaged approach sacrifices detailed information in a layer but retains nonlinear characteristics overall that create instabilities.

## 2. Governing equations

For the purpose of this paper, it is necessary to formulate an idealized canopy model that assumes a canopy that 1) has a uniform vertical leaf area distribution, 2) has a uniform vertical distribution of the drag coefficient of the canopy elements, and 3) is located along a uniform slope. For the sake of analytical convenience, the simplest form of governing equations is required. Fundamentally, we pose a hypothesis that *the regimes of canopy thermal-slope flows are governed largely by buoyancy and drag forces*. In methodology, we use the layer average concept that was developed by Fleagle and McNider. Fleagle (1950) conceived of a simple layer model that treated a drainage layer equivalent to a parcel being cooled on a slope. McNider (1982), in a more formal development, then incorporated ambient stratification. In his layer model, McNider simply assumed that the model equations are valid for a drainage layer with idealized conditions; that is, airflow is uniform along the slope except for a constant gradient of potential temperature consistent with the ambient stratification. In this paper, the assumptions of McNider’s layer model are invoked for a drainage layer within canopy. This is rational because the data measured at the Niwot Ridge AmeriFlux site (Fig. 3; also see Fig. 8 in Yi et al. 2005) indicate that wind speed and potential temperature were relatively uniform in a drainage layer within the canopy.

Thus, the model equations of canopy flows in a slope-aligned coordinate system are given by

where *u* and *θ* are the layer-average, along-slope velocity and the potential temperature, respectively; *g* is the gravitational acceleration; *α*, the slope angle; *c _{D}*, the drag coefficient;

*a*, the leaf area density; and

*γ*, the ambient lapse rate;

*θ*

_{0}(

*z*) =

*θ*

_{00}+

*γz*is the ambient potential temperature;

*θ*

_{00}is the potential temperature at

*z*= 0;

*L*

_{c}is the average heat flux for the canopy layer; and

*f*

_{u}(

*u*,

*θ*) and

*f*

_{θ}(

*u*,

*θ*) are two functions of

*u*and

*θ*, defined in Eqs. (1) and (2), respectively. The difference

*θ*−

*θ*

_{0}in Eq. (1) represents deviation of the flow-dependent and ambient potential temperature, respectively, from the constant along-slope

*θ*gradient. In this paper,

*θ*

_{0}will be approximated by

*θ*

_{00}in a drainage layer because change in ambient potential temperature within a drainage layer is ignored (<

*γh*/

*θ*

_{00}∼ 10

^{−4}, if

*γ*= 6.5°C km

^{−1},

*h*= 15 m, and

*θ*

_{00}= 290 K, where

*h*is canopy height). The essential difference between the present model and the Fleagle–McNider model is the treatment of the frictional term in Eq. (1). The friction force in the Fleagle–McNider model is assumed to be proportional to speed; that is,

*F*= −

*k*

*u*for case of bare slope surface, where

*k*is a proportional coefficient. The frictional term in this paper is taken to be

*F*= −

*c*|

_{D}a*u*|

*u*for the case of the vegetated slope surface (Mahrt et al. 2000; Yi et al. 2005; Yi 2008).

## 3. Instability analysis

In this section, we employ the self-organization theory, developed during the early 1970s (Nicolis and Prigogine 1977; Haken 1977; Yi 1995), to identify which parameters are control parameters for transition between different regimes of TI canopy flows.

### a. Steady states

Steady states are defined as stationary states, determined by

where the subscript *s* indicates steady state. The equations of steady states can be solved if *u*_{s} > 0 (downslope flow) or *u*_{s} < 0 (upslope flow).

In the downslope case, Eq. (3) becomes

and the steady state solutions are

Real solutions require that canopy flows are colder than environment: Δ*θ* = *θ*_{00} − *θ*_{s} > 0.

In the case of upslope, Eq. (3) becomes

and the solutions of steady states are obtained:

Physical solutions require canopy flow warmer than environment: Δ*θ* = *θ*_{00} − *θ*_{s} < 0.

### b. Linear stability analysis

#### 1) Linearized equations

To determine the stabilities of the steady states, one assumes that *f*_{u} and *f*_{θ} can be Taylor expanded in a formal power series of *u*′ = *u* − *u*_{s}, *θ*′ = *θ* − *θ*_{s} around steady states (*u*_{s}, *θ*_{s}). Omitting the nonlinear terms from the Taylor expansion, the linearized equations are obtained:

where

The principle of linearized stability (Nicolis 1995) indicates that the stabilities of steady states (*u*_{s}, *θ*_{s}) are equivalent to that of trivial solutions (*u*′ = 0, *θ*′ = 0). Therefore, the stabilities and properties of the steady states can be determined by examining the evolution of the linearized equations.

#### 2) Characteristic equations

Equation (11) admits solutions that depend on time exponentially:

where (*û*, *θ̂*) are initial perturbations and *ω* is a characteristic exponent. Substituting (13) into (11), one obtains

The time-dependent solutions of (13) must satisfy the characteristic equation

where *T* and Δ are, respectively, the trace and the determinant of the Jacobian matrix of (11):

#### 3) Cold canopy flows (Δ > 0)

where *σ* = *gγ*sin^{2}*α*/*θ*_{00}, and the eigenvalues of steady state (*u*_{s}, *θ*_{s})

where

is a critical value of the potential temperature difference between the environment and canopy flows. Real eigenvalues require

Under condition (22), the steady state (*u*_{s}, *θ*_{s}) in Eqs. (6) and (7) is a stable node in the phase space because *ω* ± < 0. It is noticed that steady states

are also mathematical solutions of Eqs. (4) and (5) under the condition of cold canopy flow (Δ*θ* > 0). It can be easily verified that the steady state (23) is an unstable node in the phase space because *ω* ± > 0.

The eigenvalues of the steady states in Eq. (20) becomes complex conjugates

under the condition

where *i* represents the unit imaginary number in Eq. (24). Under condition (25), the steady state in Eqs. (6) and (7) is a stable focus, whereas the steady state in Eq. (23) is an unstable focus in phase space. Because of the multiplicity of eigenvalue solutions, (24) are odd, and

where Δ*θ* = Δ*θ _{c}* is a bifurcation point (Hopf bifurcation; Nicolis 1995). The bifurcating solutions are time periodic. The period of oscillation of the canopy flows is determined by

#### 4) Warm canopy flows (Δ < 0)

Real *ω*^{±} must satisfy the condition

Similar to the case of cold canopy flow, another mathematical solution for the steady state (*u*_{s}, *θ*_{s}) under warm condition (30) also exists:

Under condition (30), the steady state (*u*_{s}, *θ*_{s}) in Eq. (31) related to downslope flow becomes an unstable node, whereas the steady state (*u*_{s}, *θ*_{s}) in Eqs. (9) and (10) related to upslope flows is a stable node in phase space.

It can be easily verified that −Δ*θ _{c}* is a bifurcation point for warm canopy flows. Beyond the bifurcation point, that is,

the steady state (*u*_{s}, *θ*_{s}) in Eqs. (9) and (10) related to upslope flows is a stable focus in phase space; the related eigenvalues are

The steady state (*u*_{s}, *θ*_{s}) in Eq. (31) related to the downslope flow is an unstable focus in phase space; the eigenvalues are

The time period of oscillation of canopy flows under condition (32) is given by

The stabilities and properties of canopy flows varying with the bifurcation parameter are summarized in Fig. 1. The upslope flows are stable under conditions of warm canopy flows (Δ*θ* ≤ −Δ*θ _{c}*), whereas the downslope flows become stable under conditions of cold canopy flows (Δ

*θ*≥ Δ

*θ*). Because canopy flow temperature is close to the ambient temperature (Δ

_{c}*θ*> |Δ

_{c}*θ*|), the canopy fluid system performs an oscillatory motion. The period of time oscillation is determined by Eq. (27) with the condition 0 ≤ Δ

*θ*< Δ

*θ*or by Eq. (35) with the condition −Δ

_{c}*θ*< Δ

_{c}*θ*≤ 0. Equations (27) and (35) can be combined in the form

here *τ _{c}* = 2

*π*/

*σ*is a time period at Δ

*θ*= 0, and |Δ

*θ*| < Δ

*θ*.

_{c}## 4. Concluding remarks

The switch conditions of TI canopy flows from the patterns of steady states to oscillation are derived by a nonlinear dynamics approach. The potential temperature difference, Δ*θ* = *θ*_{00} − *θ*, between the environment and canopy flow serves as an order parameter, also called a bifurcation parameter in dynamic system theory. The critical value of the order parameter Δ*θ _{c}* is found to be dependent on the slope

*α*, ambient lapse rate

*γ*, canopy drag coefficient

*c*, and leaf area density

_{D}*a*; that is, Δ

*θ*=

_{c}*λ*sin

*α*/(

*c*) in Eq. (21).

_{D}aThe bifurcation diagram is shown in Fig. 1. Beyond the bifurcation point (or under cold canopy flow conditions, Δ*θ* > Δ*θ _{c}*) downslope canopy flows are asymptotically stable, whereas upslope canopy flows are unstable. Below another bifurcation point (or under warm canopy flow conditions, Δ

*θ*< −Δ

*θ*), the stabilities of downslope with upslope canopy flows exchange with each other; that is, upslope canopy flows become stable, whereas downslope canopy flows become unstable. Because canopy flow temperature is close to ambient temperature (|Δ

_{c}*θ*| < Δ

*θ*), the eigenvalues of the characteristic equation become complex conjugates, and the oscillation behavior of canopy flows is bifurcated. Downslope canopy flows appear as a stable focus in phase space under the condition 0 ≤ Δ

_{c}*θ*< Δ

*θ*, and they become an unstable focus under a slight warm condition (−Δ

_{c}*θ*< Δ

_{c}*θ*≤ 0). However, upslope canopy flows exchange their stability with downslope flows from a stable focus under the condition −Δ

*θ*< Δ

_{c}*θ*≤ 0 to an unstable focus under the condition 0 ≤ Δ

*θ*< Δ

*θ*. At both bifurcation points (± Δ

_{c}*θ*), bifurcations satisfy the conditions of the Hopf bifurcation; hence, they are called Hopf bifurcations.

_{c}The time-oscillating period of canopy flows is determined by Eq. (36), depending on the gravity (buoyancy) factor *gγ*/*θ*_{00}, slope *α*, and the order parameter Δ*θ*. The minimum time period *τ _{c}* = 2

*π*/

*σ*, occurs under the condition that canopy flow temperature is equal to ambient temperature. The estimated value of

*τ*is about 80 min if

_{c}*g*= 9.8 m s

^{−2},

*γ*= 0.0065 K m

^{−1}, and

*α*= 5°, which is close to the value of 1.5 h obtained in the numerical case studied by Doran and Horst (1981). The interesting feature of the period given in Eq. (36) is that the time period increases with the increase of absolute difference of temperature between canopy flows and the environment, as shown in Fig. 2. The time period

*τ*will become infinitely large as Δ

*θ*→ Δ

*θ*. The increase of the time period with the absolute value of the order parameter (|Δ

_{c}*θ*|) is important in understanding why there is a change in the motion patterns in the canopy from oscillatory to steady state at the bifurcation points ± Δ

*θ*, where the time period becomes infinitely long.

_{c}The stability analysis by nonlinear dynamics approach is quite abstract. Figure 3 provides better visualization of how the order parameter controls the regime switch of canopy flows. Under cold canopy flow conditions Δ*θ* ≥ Δ*θ _{c}* (Fig. 3a), the magnitude of buoyancy forcing in Eq. (1) is constant with a fixed value of the order parameter (Δ

*θ*), while the direction of the buoyancy forcing is persistent downslope. The magnitude of the drag force is determined by the velocity-squared law as described in Eq. (1). The direction of the drag force is always opposite to that of canopy flows. Initially, the drag force is much smaller than buoyancy force if inflow speed is slow and hence the canopy flows are accelerated. Eventually, the drag is increased to equal the buoyancy force in magnitude and be opposite in direction, and then canopy flows move downslope at a constant speed

*u*

_{s}(i.e., steady state).

Under the warm canopy flow condition, Δ*θ* < −Δ*θ _{c}* (Fig. 3b), the buoyancy force is persistently upslope. Given downslope inflows in Fig. 3b, both drag and buoyancy forces are initially upslope and result in deceleration of canopy flows. Eventually, canopy flows will stop moving downslope at some point and begin moving upslope because of the net force. This net force acts on canopy flows in the upslope direction, which is why the downslope inflow is able to result in upslope flows over the lower part of the slope under the warm canopy flow condition. It is observed that the drag direction is downward as air moves upslope over the lower part of the slope (Fig. 3b).

Because canopy flow temperature is close to ambient temperature, the buoyancy force is weak and is neither persistent downslope nor persistent upslope. The well-organized behavior of the canopy flows results from competition between drag and buoyancy forces shown in Fig. 3c and can be understood by the self-organization theory. Imagine that a parcel of air within the canopy becomes colder than the surroundings by the random action of disturbances, and the parcel air moves downslope because of the buoyancy force. A weak drag force, low because of parcel speed, acts on the downslope-moving parcel in the opposite direction to that of the air parcel. The parcel air will slow down because of the drag force and eventually move back in the upslope direction as it becomes warmer than its surroundings by downslope adiabatic compression. As the parcel air, driven by buoyancy force, moves upslope, it will be decelerated by the drag force. Meanwhile, the buoyancy force will go from upslope to downslope as a result of the upslope adiabatic expansion cooling. Therefore, the oscillations predicted by |Δ*θ*| < Δ*θ _{c}* result from the cooperation between drag and buoyancy forces, in which buoyancy force is dominant and drag force is “slaved” by the buoyancy force. The ordered spatial pattern indicates that the individual oscillatory parcels are well self-organized (highly correlated to each other), and competitive fluctuations are dominant in the self-organization process to form the macroscopic pattern.

The instability conditions of three different TI flow regimes derived in this paper should prove useful in understanding problems of measuring and modeling ecosystem–atmosphere exchanges, especially the advection issues in the eddy flux measurements (Finnigan 2008). For example, the horizontal advection errors measured at the Niwot Ridge AmeriFlux site were opposite between two different TI flow regimes: daytime warm upslope flow and nighttime cold downslope flow (Yi et al. 2008; Sun et al. 2007). The derived instability conditions elucidate fundamental reasons as to why these different regimes change, and how the strengths of flow regimes depend on the slope, leaf area density, drag coefficient, and potential temperature deficit. To the author’s knowledge, investigation of the oscillation regime of canopy flows has been limited to theoretical and numerical studies and has not been reported by rigorous observations. We believe that the oscillation regime of TI flows is more likely to occur during the transition time between day and night. The importance of transition time between day and night was emphasized for the CO_{2} transport over complex terrain by Sun et al. (2007).

## Acknowledgments

This work was financially supported in part by the PSC-CUNY Faculty Research Award (Grant 60020-3819). The author would like to thank three anonymous reviewers for numerous comments that improved the manuscript.

## REFERENCES

**,**

_{2}flux measurements in nocturnal conditions: An analysis of the problem.

**,**

_{2}advection in a sloping forest.

**,**

_{2}storage and advection conditions at night at different CARBOEUROFLUX sites.

**,**

**,**

**,**

**,**

_{2}in a mid-latitude forest.

**,**

**,**

**,**

**,**

**,**

**,**

_{2}budget in and above a forest canopy.

**,**

_{2}fluxes at three forest sites.

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

_{2}advection patterns on wind direction on a gentle forested slope.

**,**

_{2}flux estimates by eddy covariance and chamber-based model.

**,**

**,**

_{2}and sensible and latent heat fluxes during a full year in a boreal pine forest trunk-space.

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

_{2}efflux in an annual semi-arid grass, Bromus tectorum.

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

^{13}C of respired CO

_{2}in a subalpine coniferous forest.

**,**

_{2}advection.

**,**

_{2}transport over complex terrain.

**,**

**,**

**,**

**,**

**,**

_{2}budget and advective contributions to measurements of net ecosystem–atmosphere exchange of CO

_{2}.

**,**

**,**

**,**

_{2}from a very tall tower.

**,**

**,**

_{2}exchange in a high-elevation, subalpine forest ecosystem.

**,**

## Footnotes

*Corresponding author address:* Chuixiang Yi, School of Earth and Environmental Sciences, Queens College, City University of New York, 65-30 Kissena Blvd., Flushing, NY 11367. Email: cyi@qc.cuny.edu