The linkages between soil moisture dynamics and convection triggers, defined here as the first crossing between the boundary layer height (hBL) and lifting condensation level (hLCL), are complicated by a large number of interacting processes occurring over a wide range of space and time scales. To progress on this problem, a soil–plant hydrodynamics model was coupled to a simplified ABL budget to explore the feedback of soil moisture on convection triggers. The soil–plant hydraulics formulation accounted mechanistically for features such as root water uptake, root water redistribution, and midday stomatal closure, all known to affect diurnal cycles of surface fluxes and, consequently, ABL growth. The ABL model considered the convective boundary layer as a slab with a discontinuity at the inversion layer. The model was parameterized using the wealth of data already collected for a maturing Loblolly pine plantation situated in the southeastern United States. A 30-day dry-down simulation was used to investigate the possible feedback mechanisms between soil moisture and convective rainfall triggers. Previous studies, which made use of surface flux measurements to drive an ABL model, have postulated that a negative feedback was possible, which could award the ecosystem with some degree of self-regulation of its water status. According to model simulation results here, this negative feedback is unlikely. However, drastic changes in external water sources to the ABL are needed for triggering convection when soil moisture is depleted. The apparent negative feedback originated from a decoupling between the water vapor sources needed to produce convection triggers and surface water vapor fluxes.
The coupling between soil moisture, land surface fluxes, and the initiation of convection, which may lead to convective rainfall, remains an open research problem and has attracted much recent inquiry (Dirmeyer et al. 2006; Kim and Wang 2007; Kochendorfer and Ramirez 2005; Koster et al. 2004; Koster and Suarez 2004; Lawrence and Slingo 2005; Mahanama and Koster 2005; Santanello et al. 2007; Wang et al. 2007; Fennessy and Shukla 1999). The reason this problem remains vexing stems from the fact that the coupling between soil moisture and convective rainfall involves a large number of interacting processes occurring within the soil–plant–atmosphere system that vary over a wide range of space and time scales.
Belowground and surface processes involve the dynamics of water movement from the soil into the atmosphere (rooting system, plant hydrodynamics, and stomatal regulation dictating water movement from the roots and out of the stomata as water vapor after the phase transition), the canopy aerodynamics (affecting the transport of heat and water vapor from the canopy into the mixed layer), and the partitioning of net radiation (Rn) into latent, sensible, and soil heat fluxes [thereby influencing skin temperature and directly affecting the dynamics of mean air temperature and water vapor concentration in the atmospheric boundary layer (ABL)]. On the other hand, the ABL, with its unique coexistence of mechanically and thermally generated turbulence, acts as an integrator of these surface processes with larger length and slowly evolving synoptic-scale processes affecting entrainment. The dynamics of these land surface fluxes and soil–plant–atmosphere state variables control the simultaneous growth of the convective boundary layer and lifting condensation level (LCL) and thus their crossing—a necessary but not sufficient condition for triggering convection. Even in the most idealized cases, any exploration of the feedback mechanisms between soil moisture and convection triggers must account, at minimum, for all these pathways at the appropriate scales. To do so, small-scale processes, such as root water uptake, must be spatially upscaled, and fast processes, such as entrainment of turbulent fluxes, must be temporally averaged in such a way that the dynamics of each process is still preserved and the coupling between them can be properly investigated.
Numerous studies have already been conducted to explore, numerically and/or experimentally, the interaction among a subset of these processes across a wide range of spatial scales. It is not the intent of this work to review all of them but, broadly speaking, many of these studies focused on aboveground pathways with simplified accounting for belowground soil moisture dynamics (Atlas et al. 1993; D’Odorico and Porporato 2004; Findell and Eltahir 2003b; Giorgi et al. 1996; Pan et al. 1996; Trenberth and Guillemot 1996; Clark and Arritt 1995; daRocha et al. 1996). Indeed, it is possible to establish relationships between canopy transfer coefficients and bulk soil moisture without resolving the detailed dynamics of soil water movement. However, these relationships are generally reasonable when daily time scales are of interest and subdaily scale soil–plant dynamics, such as diurnal hysteresis in bulk canopy conductance (Tuzet et al. 2003) or hydraulic redistribution (i.e., reversed water movement from the roots into the soil), do not affect the dynamics of ABL growth. With regard to the latter example, Lee et al. (2005) incorporated a simple parameterization of root hydraulic redistribution in an atmospheric general circulation model and found that it can contribute an approximately 40% increase in transpiration over a 3-month dry season period. With regard to the former example, a number of models have already documented the role of soil–root hydrodynamics on water vapor fluxes and soil moisture redistribution. For example, the model of Amenu and Kumar (2007) for root water uptake includes soil moisture dynamics and explicitly resolves the root hydraulic system. They did not consider the important effects of “microscale” soil moisture dynamics adjacent to the rooting system (or other effects on ABL dynamics and its concomitant atmospheric state variables). In fact, microscale soil moisture dynamics is complicated by radial water flow from neighboring soil elements into the rooting system occurring on scales of millimeters (Tuzet et al. 2003), vertical redistribution of soil moisture occurring on scales of meters (“macroscale”), and other “nonlocal” processes, such as hydraulic redistribution (Dawson 1993; Emerman and Dawson 1996; Williams et al. 1993). Thus, the inclusion of such prognostic equations for modeling soil–moisture dynamics may reveal some feedback mechanisms between soil moisture and convective triggers that have not been historically explored.
For the purpose of quantifying such feedbacks between soil moisture and convective rainfall potential, Findell and Eltahir (2003a) developed a framework based on the low-level temperature and humidity structure of the early morning atmosphere [convective triggering potential (CTP) and low-level humidity index characterization (HIlow)]. They used the CTP-HIlow to catalog possible feedbacks of soil moisture on summertime convective precipitation in the continental United States and classified the continental United States into regions according to the possible feedbacks introduced by soil moisture (Findell and Eltahir 2003b). These models identified the southeastern (SE) United States as a region with likely positive feedback, meaning that subsequent convective precipitation events are highly correlated to antecedent wet soil moisture states and to high evapotranspiration rates during the late morning and early afternoon.
On a more local scale, the ensemble of direct surface measurements of antecedent soil moisture content and air relative humidity just prior to a convective rainfall event may offer clues on these feedbacks (Juang et al. 2007b). A methodology to conditionally sample precipitation and identify whether it is convective using a long-term point rainfall record can be developed if a concomitant time series of measured surface sensible heat flux and ancillary meteorological variables are available (Juang et al. 2007a,b). In Juang et al. (2007b), the conditional sampling scheme makes use of a simple ABL slab model driven by measured surface fluxes and surface environmental conditions to estimate ABL height (hBL) and lifting condensation level (hLCL). Convective precipitation was identified when the convectively modeled hBL intersects the modeled hLCL just prior to a recorded rainfall event. This identification scheme must be viewed as a simplification because convection and convective cloud formation are triggered when a fluid parcel reaches the level of free convection (LFC), which is higher than the LCL. Hence, for a fluid parcel originating within the boundary layer to reach the LFC, the parcel must carry enough energy to overcome the negative buoyant force between hLCL and LFC, which is not guaranteed when the hBL intersects the hLCL. It is also known that convection can be triggered by other vertical velocity perturbation, such as surface physiography, atmosphere structure, or surface heterogeneity (Rogers and Fritsch 1996). However, Juang et al. (2007b) were able to successfully describe convective precipitation timing in 92% of the rainfall events identified as convective, even when using an averaged early morning temperature lapse rate from a nearby airport station for the entire growing season. Hence, identifying convection triggers with such an intersection is a reasonable surrogate for the precise triggering mechanism. The advantage of using this surrogate is that hLCL, contrary to LFC, is a local ABL property that can be derived independently from the conditions above it.
Although the analysis in Juang et al. (2007b) supported the conclusions of Findell and Eltahir (2003b) in that most convective rainfall events were consistent with their classification, they did identify several convective precipitation events under conditions of low soil moisture and low atmospheric humidity. Juang et al. (2007b) postulated that these soil moisture and relative humidity states be interpreted as an indication that a negative feedback between soil moisture and subsequent precipitation may exist.
Motivated by these analyses, the primary objective of this work is to explore via model calculations the conditions that can sustain such a negative feedback (if it exists). To achieve this objective, two tasks must be completed: 1) couple an ABL dynamics model with a detailed soil–plant model that resolves both, the microscopic and macroscopic soil moisture dynamics; and 2) use this model to explore numerically the possible feedback mechanisms between antecedent soil moisture and convective precipitation precursor for similar synoptic-scale conditions and atmospheric forcing. The soil–plant model should be able to simulate root water uptake mechanistically and capture vegetation responses to water stress, including soil moisture modulations of the Bowen ratio (β), in a manner consistent with what is known about the soil–plant hydraulic system. The ABL dynamics model should be sufficiently realistic to provide reasonable estimates of entrainment fluxes, which brings information from synoptic-scale conditions to ABL state variables and consequently to hBL and hLCL. The coupled model will be applied to the pine plantation described in Juang et al. (2007b) and simulation results will be used to investigate possible feedback mechanisms between soil moisture on convective precipitation triggers. Although Juang et al. (2007b) included three representative cover types for the region in their analysis, this exercise will focus on the pine plantation because almost all parameters needed for model simulations (see Table 1) were independently measured, thereby enhancing the reliability of the model results. As will be shown later, model outputs follow similar trends in the pairs of soil moisture–relative humidity, triggering convection as in Juang et al. (2007b), thereby providing robustness to the conclusions. The same criteria as in Juang et al. (2007b) for the trigger of convection will be considered here. The feedback will be evaluated from the behavior of the hBL and hLCL as the soil moisture becomes a major limitation to transpiration.
2. Model description
a. Soil–plant model
The soil–plant model used here is similar to the one presented in Siqueira et al. (2008). However, for completeness, a brief description is provided. Figure 1b presents a schematic diagram of the model framework.
In this soil–plant model, root water uptake is hydraulically controlled and is a function of the difference between local root water potential and local soil water potential at the root–soil interface. Water movement within the soil, needed to estimate water potential in the root vicinity, is governed by the Richards equation (Campbell 1985; Richards 1931). The main novelty of the Siqueira et al. (2008) model is that instead of solving the entire highly nonlinear 3D Richards equation for complex geometry imposed by root morphology (hence, extremely numerically demanding), the model assumes horizontal homogeneity in root length distribution and decomposes the 3D Richards equation in two 1D directionally distinct coupled components: 1) a radial direction component for computing water flow toward (or away from during hydraulic redistribution events) the rootlet and 2) an estimation of bulk (layer averaged) vertical motion of water. These coupled components can be expressed as
where r (m) and z (m) are radial and vertical coordinates, with r being the distance from the root–soil interface in the vicinity of a rootlet and z the soil depth, t (s) is time, θ (m3 m−3) is the soil water content, ϕ (m2 s−1) is the transformed soil “matric” water potential [after applying the Kirchhoff integral transformation (Campbell 1985; Redinger et al. 1984)], K (m s−1) is hydraulic conductivity, qz (m s−1) is water flow in the vertical direction, and Es (s−1) accounts for soil water evaporation. In the equations, angle brackets represent a layer-averaged value and parentheses and square brackets following dependent variables indicate dependence (function of) relationship. Note that the radial component refers to distance of a given point in the soil to the closest rootlets (on the order of millimeter) in each layer and not the distance from a single tree (on the order of meters).
Water retention and hydraulic conductivity functions required to close the system are given by
where ψ (m) is soil water potential; θs (m3 m−3), ψe (m) and Ks (m s−1) are the saturation water content, air entry water potential, and saturated hydraulic conductivity, respectively; and b is an empirical parameter that varies with soil type (Campbell 1985).
For consistency, the layer-averaged variables should be estimated using functional relationships given by Eq. (2), with 〈θ〉 being the mean value of the soil moisture at each z location being expressed as
where λ (m m−3) is root length density assumed to be uniformly distributed horizontally, rr (m) is root radius and R (m) is the radial domain size, which is the average (halfway) distance between rootlets for each layer.
When root water uptake qr (m s−1) is hydraulically controlled, it is given by
where Kr (s−1) is a root membrane permeability, and ψr (m) is the root pressure referenced to ground level. Hence, boundary conditions on the radial flow equation are a zero flux at R (because of symmetry) and Eq. (4) at the root–soil interface.
Transpiration TR (m s−1) is assumed to be identical to the sap flow (SF) per unit of ground area (no capacitance) and is given by
where ψυ (m) is the leaf pressure, χ (s) is root-to-shoot plant hydraulic resistance and ZR (m) is rooting system depth. In the model calculations here, ψr is the main link between the soil and the leaf system. It is assumed to be constant with depth (adjusted for hydrostatic contributions) but evolving in time.
Stomatal response to a drying soil is modeled with a logistic function given by (Tuzet et al. 2003)
where gs (m s−1) is stomatal conductance, and g0 (m s−1) and gmax (m s−1) are residual and maximum stomatal conductance, respectively. The fψ is a reduction function with empirically determined sensitivity parameter sf (m−1) and reference potential ψf (m).
The partial differential equations (1a and 1b) are discretized using a control volume approach with central differencing scheme for spatial derivatives and implicit scheme for time derivatives and solved using an iterative Newton–Raphson method. Only Eq. (1b) is solved for root-free zone with the addition of an evaporation term. To account for vertical gradients in soil properties, the method of Ross and Bristow (1990) was adapted to the discretization technique used here.
Evaporation required for solving these equations [(1a) in the root zone and (1b) in the root-free zone] is modeled using a simple water vapor diffusion equation with fractional relative humidity as the driving gradient. Thus, evaporation also depends on soil temperature Ts, which is computed by solving the heat flow equation (not shown here). Additional equations necessary to “close” the problem are provided by the energy balances for leaf and canopy air and a radiation partitioning model. In contrast with Siqueira et al. (2008) in which the model was driven with atmospheric forcing variables measured just above the canopy, the atmospheric inputs here are synoptic-scale variables prescribed above the ABL.
b. ABL model
Modeling the surface layer structure and coupling this surface layer to the convective boundary layer aloft is necessary to link the atmospheric forcing variables measured just above the canopy to synoptic-scale variables. A logarithmic profile for the time-averaged state variables and constant flux within the surface layer were assumed, with surface layer height (hs) set to be the minimum of the absolute value of the Obhukov height (L) and a preset fraction of the hBL. The minimum avoids numerical instabilities at very high values of L during transition phases (i.e., stable-to-unstable or unstable-to-stable atmospheric stability states). Monin–Obhukov similarity theory was used to derive the transfer coefficients needed to compute the surface fluxes (Garratt 1994).
As earlier noted, a necessary but not sufficient condition for convective cloud formation is that the hBL intersects the hLCL as shown in Fig. 1. To explore under what soil moisture conditions this intersection occurs, it is necessary to describe the ABL dynamics with reasonable accuracy. Multilayer models that resolve the vertical structure of the turbulent field, such as large eddy simulation (LES) or even higher-order Reynolds-averaged Navier–Stokes (RANS) equations are computationally too demanding for such long-term applications (e.g., dry-down periods of several weeks are required in such simulations). Assuming that (1) the turbulence within the boundary layer is intense enough to sustain well-mixed conditions, (2) the atmospheric variables are horizontally homogeneous (no net large-scale advection), (3) the heat capacity of the surface layer is negligible compared to whole boundary layer, and (4) the inversion that caps the boundary layer is sufficiently thin to be approximated by a step change, the simplified equations for heat (both sensible and latent) transfer within the ABL can be expressed as (Driedonks 1982; Tennekes and Driedonks 1981)
where Θm (K) is the time-averaged mixed-layer potential temperature, ΔΘ (K) is the discontinuity jump at the inversion base, γΘ (K m−1) is the potential temperature lapse rate at hBL, w (m s−1) is the vertical velocity component, and primes represent turbulent excursions from mean values (typically half hourly to hourly). For the fluxes, the overbar represents time averaging and subscripts h and a indicate variables at hBL and at the surface, respectively, such that a (K m s−1) (= Ha/ρ Cp) is the surface sensible heat fluxes with ρ (kg m−3) and Cp (J kg−1 K−1) being air density and specific heat at constant pressure, respectively; and h (K m s−1) is the entrainment sensible heat flux. Here, LWnet (W m−2) is the net longwave (LW) radiation exchange of the mixed layer. Longwave exchange is relevant for long-term adjustments of the boundary layer state (Kim and Entekhabi 1998), and it may also have an important cooling effect in the late afternoon. The longwave radiation budget model applied here is similar to Brubaker and Entekhabi (1995) but with emissivity functions given by Garratt and Brost (1981). With the exclusion of LWnet, identical equations can be written for water vapor mixing ratio [Q (kg kg−1)] replacing Θ by Q (Driedonks 1982), with a (m s−1) (=LEa/ρ lυ) being surface water vapor flux, where LEa (W m−2) and lυ (J kg−1) are surface latent heat and latent heat of vaporization, respectively; h is the entrainment latent heat flux; and γq being the lapse rate for water vapor.
The equations above still require closure approximations, as usual in turbulence modeling. Assuming that the budget equation of ABL turbulent kinetic energy (TKE) can be represented by a balance between the TKE rate of change, TKE flux convergence, and buoyancy flux at the inversion base (Tennekes and Driedonks 1981) along with parameterizations for the rate of change term given in Tennekes (1973) and flux convergence term by Zilitinkevich (1975), an equation for hBL can be written as
where h (m2 s−3) is the entrainment buoyancy flux (includes heat and water vapor); σw, w* [= (hBLs)1/3], and u* (m s−1) are turbulent, convective, and friction velocity scales, respectively; T and Tυ (K) are temperature and virtual temperature, respectively; g (m s−2) is gravitational acceleration; and cF, cT, and A are similarity parameters (see Table 1 for values).
Implicit in this model formulation is the canonical shape of the state and flux profiles within the ABL. Figure 2 illustrates the model profiles along with an ensemble-averaged potential temperature and water vapor mixing ratio profiles from late afternoon sounding data (when a clear boundary layer was identified) at a nearby airport station [Piedmont Triad International Airport (GSO); see section 3a for details] normalized by hBL and referenced to Θm or Qm. The assumed shapes of the state variable profiles are consistent with the ensemble shape values derived from the sounding data. Although large variability does exist, primarily as a result of the jump (or its assumed thickness), the well-mixed condition seems to be a robust approximation. Additionally, slab models, such as the one used here, have been successfully used in several recent numerical and experimental studies, providing some confidence in their realism (Cleugh et al. 2004; Conzemius and Fedorovich 2007; Garcia et al. 2002; Kim et al. 2006; Lewis 2007; Pino et al. 2006a,b; Villani et al. 2005).
Because the application here requires continuous simulation, the nocturnal boundary layer (NBL) must also be modeled. Conceptually, the NBL is a stable layer that grows after the unstable-to-stable atmospheric stability transition occurred, typically around sunset, underneath the now termed residual layer [(RL) formerly convective ABL; see Garratt (1994), p. 146 for illustration]. This layer suppresses the energy supply for TKE generation, which will decay in the RL above the NBL and, consequently, invalidate the underlying assumptions implied in Eqs. (7)–(10). A simple NBL model is used here because there is no particular interest in capturing all the detailed dynamics of the NBL, given that it will simply provide a physically sound near-surface potential temperature lapse rate that the ABL will need to overcome in the early morning hours. During these hours, this potential temperature lapse rate may have a significant control on ABL evolution and may be a function of the soil water state (Findell and Eltahir 2003a; Santanello et al. 2005).
A stationary NBL model prescribes an NBL height (he) immediately after the unstable-to-stable atmospheric stability transition occurs (Garratt 1994). This height remains constant in time thereafter as long as nocturnal conditions prevail. The mean potential temperature profile is assumed to be linearly distributed within the NBL, and no exchange between the NBL and the RL aloft is allowed. The NBL energy balance along with geometric considerations provides the additional equation needed to solve the soil–plant–atmosphere system during nighttime. Around sunrise, at the stable-to-unstable atmospheric stability transition, the hBL is initialized to a prescribed value and the mean potential temperature profile is initialized with a lapse rate estimated from the NBL model up to he. The synoptic-scale lapse rate is used above the he. The water vapor that accumulated during the night as a result of evaporation (and some transpiration) is allowed to be flushed out at the stable-to-unstable atmospheric stability transition, thereby conserving the mass of water in the system. The time evolution of the hLCL was computed using the modeled ABL water vapor mixing ratio and potential temperature using standard formulations.
The model described above was parameterized for the Duke Forest Pine Plantation, located near Durham, North Carolina. Briefly, the site is composed of a 23-yr-old (in 2006) Loblolly pine stand planted over a Enon series soil with a clayey loam texture in the upper 0.3 m and a clay soil from 0.3 to 0.7 m where the bedrock typically resides (Oren et al. 1998; Stoy et al. 2007). Energy and mass fluxes are currently measured with an eddy covariance system composed of a CSAT-3 triaxial sonic anemometer (Campbell Scientific) and a Licor 7500 infrared gas analyzer (Licor), along with other environmental variables. Soil moisture was measured by a network of 24 CSI 615 vertically arrayed rods (Campbell Scientific) in the top 30 cm (Katul et al. 2007). When compared to sap flow measurements, it was shown that the top 30 cm explained much of the root water uptake activity (>90%). The details of the site characteristics and measurements are described elsewhere (Oren et al. 1998; Siqueira et al. 2006; Stoy et al. 2007). The soil–plant parameters pertinent to the model calculations are shown in Table 1. Hereafter, unless otherwise stated, the term soil moisture refers to the root-zone depth-averaged soil moisture in the top 30 cm (〈〈θ〉〉).
The results presented are model outputs for a 30-day dry-down period with repeated daily cycle of climatic forcing initialized assuming well-watered soil conditions. To get reasonable temperature distribution for soil, canopy, and ABL, a spin-up run, in which no soil moisture dynamics were considered, was performed prior to the main simulation period. The spin-up run was for a period of time sufficiently long such that the assumed initial temperature memory was lost.
a. The forcing
Climatic forcing including incoming free atmosphere (FA) shortwave (SW) and LW radiation and u* are described as repeated diurnal cycles (see Fig. 3). These forcing variables are intended to represent typical drivers and are derived from ensembles averaging the summertime (June–August) surface flux data each half hour at the site. Close to the surface, u* is the outcome of the imbalance between synoptic-scale pressure and Coriolis forces. Here, u* is prescribed in lieu of these pressure fields and the geostrophic winds because it is more relevant to surface flux estimation and is consistent with the pine plantation roughness. However, for the synoptic-scale potential temperature and water vapor lapse rates, sounding measurements at the nearby Greensboro, North Carolina, airport (GSO; 36°05′N, 79°57′W) are used. The GSO is located 79 km west of the site and 270 m above sea level (Juang et al. 2007b). The sounding data are maintained by the Department of Atmospheric Science at the University of Wyoming and collected at 0000 and 1200 UTC [0700 and 1900 local time (LT), respectively] everyday. The representative lapse rates for both scalars are estimated by averaging early morning profiles from 300 to 3000 m for the growing season of 2006 (June–August). Their numerical values are 5 × 10−3 K m−1 for γΘ and −1.45 × 10−6 kg kg−1 m−1 for γQ. Repeated climatic forcing of radiation and friction velocity for a period of 30 days does not reflect realistic atmospheric forcing conditions given their largely stochastic nature (e.g., in the wind fields). However, the use of such idealized forcing is to filter out the transient ABL dynamics caused by variability of the forcing from the ABL dynamics as a result of soil water limitations. A more realistic stochastic forcing would mask these effects and are beyond the scope of this study.
b. Soil–plant hydraulics
To account for soil moisture controls on land surface fluxes, logistic function parameters [sf and ψf; see Eq. (6)] are derived from xylem vulnerability curves, which are developed from the relationships between xylem conductivity and water potential for different xylem tissues (roots, trunk, and branches; Tyree and Sperry 1989) and also measured at the pine plantation site and reported elsewhere (Hacke et al. 2000). The relationship between xylem conductivity and water potential emerges because at low (highly negative) water potential, xylem cavitation occurs, thereby reducing its conductivity. To obtain the parameter values (see Table 1), it was assumed that stomata close to avoid xylem runaway cavitation. Figure 4a presents the logistic function variation with decreasing soil water potential (i.e., becomes more negative). The effects of such a leaf-based logistic function on the entire canopy conductance are explored using the model results for the dry-down numerical experiment. In Fig. 4b, the midday (1000–1500 local time) leaf-specific canopy conductance (Gcl) from the model is compared to sap flow data reported from a short-term intensive experiment conducted at the site in the growing season of 1994 (Oren et al. 1998). Root membrane permeability (Kr) controls the amount of hydraulically lifted water and somewhat the onset of water stress and the slope of the response of canopy conductance to soil moisture (Siqueira et al. 2008). Because no direct measurements are available for Kr, its value (see Table 1) was established to optimize soil moisture dynamics. This is the only “fitted” parameter from the list of parameters reported in Table 1. Although the model slightly overestimates the soil moisture at the onset of water stress, it reproduces the slope of the stomatal closure and provides some confidence in the realism of the model to simulate vegetation water stress from soil–plant hydrodynamics principles and independent soil and plant hydraulic parameters.
Figure 5 shows the soil moisture dynamics for the 30-day dry-down simulation. After an initial phase of high drainage, lasting until field capacity is achieved (∼3 days), the diurnal cycle in soil moisture becomes apparent. Note the drying during the day and the recharge during the night as a result of hydraulic lift (explicitly resolved in the model). When soil moisture is largely depleted, hydraulic lift is gradually shut down and transpiration is inhibited (resulting in almost steady-state soil moisture, whose value is close to the minimum recorded value at the site over an 8-yr period as shown in Katul et al. 2007).
c. Surface energy fluxes
In Fig. 6, measured and modeled surface energy fluxes for well-watered soil moisture conditions are compared for one particular day in 2006, with similar climatic forcing as those assumed in Fig. 3 (i.e., the forcing during this day are representative of the ensemble averages). Model results are shown for a steady-state daily cycle after a spin-up run with repeated daily forcing. The rapid increase in the modeled latent heat just after sunrise represents a rapid “flush” of mean water vapor concentration trapped in the NBL when the atmospheric stability switches from stable to unstable (and hence is not related to transpiration but to storage fluxes). Figure 6 suggests that model environmental conditions favor lower β’s when compared to the measurements. As a matter of fact, root-zone soil moisture values on this particular day were on the order of 0.25 m3 m−3. Although higher than the critical value (≈ 0.2), this soil moisture could affect surface fluxes, producing more sensible heat. Surprisingly, the NBL model also seems to capture reasonably well the nocturnal sensible heat fluxes despite the crude parameterizations adopted here.
The components of the energy balance are shown for both model and measurements in Fig. 7, where the sum of surface latent (LEa) and sensible (Ha) heat fluxes is plotted against net radiation Rn. Note that the surface energy budget derived from the model is conserved on a daily basis. However, the components of the energy fluxes were measured independently and an imbalance of up to 40 W m−2 for a daily cycle was reported for this site (Stoy et al. 2006), which is consistent with the maximum differences shown in Fig. 7 (occurring at the higher end of the energy budget).
Because model results in Figs. 6 and 7 represent a steady-state daily cycle with repeated forcing while the measurements are the result of “stochastic” natural forcing, one-to-one comparisons were not possible. Hence, the comparisons in Figs. 6 and 7 cannot be considered a rigorous model validation, but the agreement level in surface fluxes and energy balance is a good indication that the model is appropriate for the application intended here.
d. Boundary layer dynamics
Figure 8a shows how soil–plant hydraulics regulates daytime β, which then affects the growth of the ABL. Because sensible heat is more effective than latent heat in buoyant production of TKE, more sensible heat allows the ABL to grow at a faster rate and become deeper. This is demonstrated in Fig. 8b, which shows the daily maximum hBL as a function of daily averaged soil moisture. Also in Fig. 8b, the maximum, minimum, and averaged hBL derived from the evening sounding data at GSO airport are reported. Model predictions are well within these recorded values except for extremely strong water stress conditions, which did not occur in 2006. As a matter of fact, although the very extreme (〈〈θ〉〉<0.1) events are rare in occurrence, the transitions states are actually close to the mode of soil moisture probability density function (PDF), as evidenced by the inset in Fig. 8, which presents the growing-season soil moisture PDF for the years 2000–04 (Juang et al. 2007b). This is not surprising given the fact that, at these stages, soil moisture exerts strong controls on transpiration and has a negative feedback on its own depletion. This is a good indication that, despite of all the simplifying assumptions adopted in the model, the sharp transition imposed by the soil–plant hydraulics is realistically reproduced.
The mean near-surface potential temperature lapse rate computed by the model through the NBL parameterization ranges from 0.019 K m−1 (under water stress) to 0.032 K m−1 (for well water condition, which is consistent with Santanello et al. 2007), possibly as a result of warmer surface temperature under water-limited conditions. The ensemble-averaged clear-sky early-morning sounding lapse rate for the first 300 m was 0.0221 ± 0.0101 K m−1, which suggests that the simple NBL parameterization provides reasonable near-surface potential temperature lapse rates.
Because the model appears to capture all the key components responsible for coupling the soil–plant hydrodynamics with evolving ABL, the feedbacks of soil moisture dynamics on convection triggers is explored next. Figure 9 presents the hBL evolution for two different days, along with the computed hLCL derived from boundary-layer state variables. As mentioned before, a necessary but not sufficient condition for convection, and consequently convective precipitation, to occur is that the hLCL intersects hBL.
For illustration, the two days selected represent the end-member extremes of soil moisture states, with a wet state shown in Fig. 9a and a dry state shown in Fig. 9b. During early morning hours, the hBL increases at a slower rate (see Fig. 9a) because it has to expand against a steeper potential temperature lapse rate that developed during the night in the NBL. Again, recall that although the NBL height is constant during the night, the temperature lapse rate evolves in time because of the radiation balance. After overcoming the stable inversion layer, hBL grows faster given the higher surface sensible heat flux available. Later in the day, the buoyancy flux loses strength and the ABL starts to flatten. On the other hand, hLCL increases more rapidly in the early morning, slows down in the early afternoon, but then speeds up in the late afternoon. Interestingly, in the early evening the hLCL decreases as a result of longwave cooling. However, model calculations here suggest that this decrease in hLCL is not sufficient to produce favorable conditions for a convection trigger.
A negative feedback of soil moisture on convection, and consequently convective rainfall, meaning a dry soil moisture state increases the likelihood of precipitation and subsequently shifting toward wet soil state, can be characterized by a higher convergence rate of the hBL to the hLCL when compared to the convergence rate for no water stress conditions. Convergence rate here refers to the decrease in height difference between hLCL and hBL as time progresses. As expected, model results show that the ABL evolution under water limitation overtakes the initial stable boundary layer faster as a consequence of the higher sensible heat and less steep near-surface potential temperature lapse rate (see Fig. 9 insets), and grows at a higher rate thereafter when compared to their high soil moisture counterpart (for the same forcing). However, for low soil moisture states, the hLCL also grows faster with no slowing down in the mid-afternoon. Additionally, because of the differences in soil, surface, and ABL temperature developments under water stress, there is no longwave-induced decrease in hLCL in the evening. Hence, no negative feedback from soil moisture onto the convection trigger is likely to occur from such “endogenous” factors. Although one day of water-stressed ABL evolution is presented in Fig. 9 (for illustration), the trend described in this figure is consistent for the entire period after the onset of water stress (not shown).
Another manifestation of the negative feedback would be if conditions in the free atmosphere change—say, because of a larger-scale or synoptic process—such that more water vapor becomes available above the ABL that can be “harvested” by the rapidly growing ABL. Hence, under such circumstances, a faster growing ABL may enhance the convection trigger. To explore this possibility, additional calculations were performed assuming different free atmospheric water vapor mixing ratio (QFA) values. The goal of these model calculations was to find the minimum QFA that can produce conditions for a convection trigger during daytime conditions. A trial-and-error procedure was utilized so that for each day, the QFA was gradually increased until first crossing of hBL and hLCL occurred near the end of the day. Figure 10 shows the model results for the ABL evolution for the same days as in Fig. 9, but now with hBL and hLCL intersecting at the end of the day (i.e., around the time when the transition from unstable to stable atmospheric stability conditions occurs).
Note that the QFA necessary for this trigger need not be the same for different days. Figure 11 presents the QFA necessary to generate convection trigger conditions at the end of the day as a function of soil moisture. Also shown are the average and the maximum water vapor mixing ratios measured from the soundings. The calculated QFA values necessary to trigger convection are within the plausible values measured by the soundings. Clearly, soil moisture limitations require larger amounts of QFA to trigger convection. More interesting is that the model results here suggest that there are three stages of soil moisture interaction with QFA needed to produce a convection trigger. At high soil moisture (stage I, Fig. 11), there is no feedback because surface fluxes are almost decorrelated from soil water. At intermediate soil water states (stage II, Fig. 11), there is a sharp transition from no soil moisture control to strong soil moisture control on QFA, characterizing a threshold, because of the strong flux-partition sensitivity to soil moisture stemming from Gcl. After this short (in terms of soil moisture) transition period, the QFA needed to trigger convection becomes independent of soil moisture (stage III, Fig. 11). Again, this behavior cannot be classified as a negative feedback because the source of the free-atmospheric water vapor is external to the system and does not originate from the soil.
In a previous study, convective precipitation conditions were identified for low soil moisture and low near-surface air relative humidity (RH) conditions (Juang et al. 2007b). As a consequence, the existence of a negative feedback of water stress on convective rainfall trigger was postulated. However, the model calculations here indicate that such a negative feedback is unlikely. To explain this apparent inconsistency between this previous data analysis and the model results here, Fig. 12 shows the data (for convective rainfall) of Juang et al. (2007b), along with the modeled near-surface air relative humidity and soil moisture just prior to the convection trigger. The model result presented in Fig. 12 refers to the simulations with modified QFA to induce a convection trigger. The model exhibits similar trends as the surface data, in which near-surface air humidity decreases for low soil water. This suggests that the moisture needed to generate the trigger under soil water limitations is mostly provided by external sources, which appears decoupled from the immediate surface conditions. This may be a plausible explanation as to why a wider scatter in the near-surface air humidity populates this soil moisture regime in their data.
Figure 13 explores this point further by showing the sources of water vapor likely to trigger convection. This figure shows the total amount of water vapor (per unit of ground area) in the ABL column (W) required to produce a convection trigger as a function of soil moisture (Fig. 13a) and, of this amount, the percentage that is originating from the surface latent heat fluxes versus the one originating from the free atmosphere (Fig. 13b). It is clear that the amount of water vapor required during low soil water states is higher as a result of a larger air volume residing in the deep ABL. However, with no soil water stress, about 30% of the water available to trigger convection originates from the surface through evapotranspiration. As the soil dries, two immediate consequences arise: 1) evapotranspiration decreases and 2) the rate of entrainment increases because it is proportional to dhBL/dt, which is higher as a result of the increased sensible heat flux. Under these circumstances, the convection trigger in the model becomes mostly dependent on the available atmospheric water vapor above ABL. In more realistic 3D situations, the source of moisture is often a result of horizontal convergence, also known to contribute to vertical velocity perturbation needed to overcome negative buoyancy force above LCL (Rogers and Fritsch 1996). The simple 1D ABL model used here is unable to account for these processes. However, this fact need not alter the basic finding here that additional moisture sources are needed to produce convection trigger. Moisture originating from horizontal convergence is still external to the system.
In this study, all vegetation structure variables, such as leaf area index (LAI), were held constant. It is possible that, depending on the timing of the drought, there might be some structural response to water stress, such as bud break delays. Leaf phenology would be the structural process most prone to respond fast enough to affect the feedbacks discussed here. However, for the months in which the climatic forcing were estimated (June–August), the LAI at the pine site was close to its maximum with little change (Siqueira et al. 2006). Under these circumstances, it is not likely that vegetation structural response would have any influence on feedbacks (at least for the 30-day time scale). Yet, it is conceivable that at longer time scales, a sustained drought could trigger structural changes that lead to different feedback mechanisms, which are outside the scope of this work. We note that field evidence of such structural changes has been reported. For example, McCarthy et al. (2007) identified an increased leaf area loss for this site after a 3-month drought period during the summer of 2002.
A soil–plant hydrodynamics model was coupled to a simplified ABL 1D heat and moisture budget to explore the feedback mechanisms between soil moisture and convective rainfall. The focus was on soil moisture controls of the hBL and hLCL dynamics and their intersection, which was used here as a surrogate for a convection trigger. The main novelty of the proposed model is that in the calculation of ABL evolution, micro- and macroscale [(O(1 mm) and O(1 m), respectively] processes affecting soil moisture dynamics are resolved.
The model was parameterized for a pine plantation, and a 30-day dry-down simulation was used to investigate possible feedbacks between soil moisture and convective rainfall triggers. The combined model was forced with cyclical incident short- and longwave radiation and friction velocity, and plausible lapse rates of potential temperature and water vapor. These forcing variables were the same through out the dry- down simulation, enabling direct assessments of the feedbacks originating from soil moisture depletion. Even though no one-to-one comparisons between data and model results were made as a result of the cyclic forcing imposed on the simulation, the model captures surface flux partitioning and water stress controls on the Bowen ratio. Furthermore, the model provides realistic ABL dynamics and was deemed appropriate for the investigation of the soil moisture–convection trigger feedback mechanisms.
According to model simulation results, neither endogenous nor exogenous factors seem to induce negative feedbacks between soil moisture and convection triggers. On the other hand, drastic changes in water vapor sources to the ABL are needed for triggering convection when soil moisture is depleted.
More broadly, the implication of such findings is that the water state of an ecosystem becomes tightly coupled to the water state in the free atmosphere as soil moisture limitations establish themselves (assuming no horizontal net advection of water or heat). If the arrival time between rain events is such that the ecosystem exchange becomes soil moisture controlled, especially if soil moisture drops below the soil–plant hydraulic threshold in relation to QFA to produce rainfall trigger, water recycling could be compromised and the restoration of soil water must rely on external water vapor input.
The coupled model developed and used in this study reproduces reasonably well the dynamics of the soil–plant–atmosphere system at a local scale. As such, the feedbacks considered here are related to the ability of the soil–plant–atmosphere system to self-regulate its own water status. This issue is relevant under future climate change scenarios in which “external” water sources may be compromised as a result of large-scale weather pattern changes. Furthermore, regional and global scale models may benefit from studies like the one presented here for mapping future “hot spots” for regions in which self-regulation of water status is dependent on external water vapor sources. This idea mirrors the work of Koster et al. (2004), who presented “hot spots maps” for regions in which parameterization of initial soil moisture state influence the feedback between the land surface and the atmosphere.
This study was supported by the U.S. Department of Energy (DOE) through the Office of Biological and Environmental Research (BER) Terrestrial Carbon Processes (TCP) program (Grants 10509-0152, DE-FG02-00ER53015, and DE-FG02-95ER62083), the National Science Foundation (Grants NSF-EAR 0628342 and NSF-EAR 0635787), and the Binational Agricultural Research and Development (BARD) fund (Grant IS-3861-06).
Corresponding author address: Mario Siqueira, Nicholas School of the Environment and Earth Sciences, Box 90328, Duke University, Durham, NC 27708-0328. Email: email@example.com