Modeling of vegetation dynamics in hydrological models for the assessment of the effects of climate change on evapotranspiration and groundwater recharge

Vegetation affects water balance of the land surface by e.g. storage of precipitation water in the canopy and soil water extraction by transpiration. Therefore, it is essential to consider the role of vegetation in affecting water balance by taking into account the temporal dynamics of e.g. leaf area index, rooting depth and stomatal conductance in hydrological models. However until now, most conceptual hydrological models do not treat vegetation as a dynamic component. This paper presents an analysis of the effects of the application of two different complex vegetation models combined with a hydrological model on the model outputs evapotranspiration and groundwater recharge. Both model combinations were used for the assessment of the effects of climate change on water balance in a mesoscale catchment loctated in the Northeastern German Lowlands. One vegetation model assumes a static vegetation development independent from environmental conditions. The other vegetation model calculates dynamic development of vegetation based on photosynthesis, respiration, allocation, and phenology. The analysis of the results obtained from both model combinations indicated the importance of taking into account vegetation dynamics in hydrological models especially if such models are used for the assessment of the impacts of climate change on water balance components.


Introduction
Future climate changes such as an increase in temperature and a decrease in summer precipitation may lead to higher evapotranspiration rates, lower groundwater recharge rates and to a prolongation of drought periods in summer, which Correspondence to: M. Wegehenkel (mwegehenkel@zalf.de)might affect water balance, water availability, and also regional distribution as well as type of vegetation (e.g.Strzepek and Yates, 1997;Hatterman et al., 2008).Such assessments of the effects of climate change on water balance were carried out by e.g. using hydrological models and meteorological data from climate change scenarios (e.g.Hattermann et al., 2008).
Vegetation affects climate and water balance of land surfaces by evapotranspiration and interception.Mainly two physiological properties, leaf area index (LAI) and stomatal conductance are the basis of evapotranspiration and interception parameterizations in complex physically based hydrological models such as TOPOG (Silberstein et al., 1999).However, conceptual hydrological models e.g.HEC-HMS (USACE, 2000) or HBV (Lindström et al., 1997) do not parameterize vegetation as a dynamic component.In these models, the seasonal evolution of vegetation properties such as LAI or rooting depth is prescribed, and e.g.monthly values of these vegetation properties are kept constant year after year and no impact of climate change on vegetation development is taken into account.However, vegetation growth processes such as photosynthesis, respiration, allocation, and phenology, which are all strongly dependent on environmental conditions, make vegetation a dynamic component.Therefore, a dynamic vegetation model should be included in hydrological models especially if these models are used for the assessment of the effects of climate change on water balance of catchments (e.g.Arora, 2002;Gerten et al., 2004;Montaldo et al., 2005;Rodriguez-Iturbe and Porporato, 2005).
In our case study, we analyzed the results of the application of two different complex vegetation models coupled with a hydrological model regarding evapotranspiration and groundwater recharge using data from a climate change scenario for the time period from 1950-2100.This analysis was carried out in a mesoscale catchment located in Northeastern Germany.One vegetation model assumes a static vegetation development independent from environmental conditions.The other vegetation model calculates dynamic development of vegetation based on photosynthesis, respiration, allocation, and phenology.

Test area
The location of the Ucker catchment with an area of 2415 km 2 is shown in Fig. 1.The actual land cover consists of 5% settlements, 4% water bodies, 13% meadows, 3% wetlands, 25% forests and 50% arable land.The areas located in the north of the Ucker catchment are covered with sandy soils with high infiltration and low soil water storage capacity.In the southern and western part, loamy soils are located.In the flat river plains, wetlands and peat soils with high groundwater tables dominate.
For the model application, we used time series covering the period from 1950-2100 with daily rates of precipitation, minimum and maximum air temperature, saturation deficit of air, wind speed, and global radiation.These time series were obtained from the Remo-UBA grid data set generated by the Max Planck Institute of Meteorology (MPI-M) Hamburg in cooperation with the Federal Environment Agency (UBA) (Jacob et al., 2008).We used the dataset based on the CO 2 -Emission-scenario A1B.The A1B scenario describes a future world of very rapid economic growth, global population that peaks in mid-century and declines thereafter, and the rapid introduction of new and more efficient technologies.The A1B scenario assumes a balance across fossil and non-fossil energy sources with an increase in CO 2 -Emission up to the mid-century in 2050 and a decrease afterwards.The other spatial data used for the simulations consists of a land cover map, a soil map, a river net map, a digital elevation model with a grid size of 50 m×50 m, and a subbasins map.More information about the catchment and the spatial data can be obtained from Wegehenkel et al. (2006).

The Theseus model
In our study, we used the hydrological simulation model THESEUS, which is coupled to a Geographical Information System (= GIS) for the spatially distributed calculation of catchment water balance (Wegehenkel, 2002).We combined two different complex vegetation modules with this hydrological model.
In both vegetation models, potential grass reference evapotranspiration ET p is determined by a modified PENMANapproach (Supit et al., 1994), stated as follows: where is the slope of the saturation vapor pressure curve in kPa • C −1 ; R na is the net radiation defined as evapotranspiration equivalent in mm d −1 ; γ is the psychrometer constant =0.65 in kPa • C −1 ; e s is the saturated vapor pressure and e a is the actual vapor pressure, both in kPa; f is an empirical constant =1.0 and c is a empirical coefficient calculated from the difference between daily maximum and minimum air temperature; u(2) is the mean wind speed at 2 m height in m s −1 .Specific potential evapotranspiration of each vegetation or crop type ETPOT in mm d −1 is given by where F (t) is a conversion factor.The dynamic vegetation model was obtained from the WOFOST6.0 model and simulates crop growth based on ecophysiological processes (Supit et al., 1994).The physiological properties of each crop are defined by a parameter file consisting of data e.g.conversion factor F (t) and temperature sum from sowing to emergence.Based on this parameter set, the dynamic model calculates phenological development, CO 2 -assimilation, transpiration, growth and maintenance respiration, distribution of assimilates on stem, leaf, fruit and root as well as dry matter formation (Supit et al., 1994).The crop growth is only limited by soil water availability and temperature stress (Supit et al., 1994).This model is restricted to agricultural crops.In former studies, this model showed a sufficient agreement between simulated and measured above ground biomass and yield (e.g.Wegehenkel et al., 2004;Marletto et al., 2005).
The static vegetation model is based on simple empirical two-dimensional table functions with linear interpolation between neighbouring table values for rooting depth, plant height, conversion factor F (t) and soil coverage of a given crop or vegetation type for the calculation of transpiration, interception and evaporation (Koitzsch and Günther, 1990).These table functions consist of a certain number of data pairs (crop parameter; Julian day of the year), which correspond to significant phenological stages of the crop or the vegetation type.Due to a heat accumulation approach, the starting point of these table functions can be adjusted to the actual meteorological boundary conditions.This model includes agricultural crops, grassland and forests.
In both vegetation models, calculation of actual evapotranspiration ET r in mm d −1 and the subdivision of ET r into transpiration and evaporation is based on a density function for potential water extraction from the soil layers combined with a reduction factor dependent on the soil layer (Koitzsch and Günther, 1990) where SCD is the relative crop soil cover (0-1); r i is a soil layer dependent reduction factor; g i is the fraction of potential water extraction by transpiration; and h i is the fraction of potential water extraction by evaporation.The reduction factor (0≤r i ≤1) is a function of available soil water in the layer.
If the actual soil water content in one layer is less than a critical layer dependent soil water content, r i decreases linearly to zero.For transpiration, r i is zero at the wilting point.For evaporation, r i is zero at air dryness.The fractions g i and h i in Eq. ( 3) are derived from the functions g (x) and h (x).
In the case of transpiration, the function g (x) describes the distribution of daily water extraction between the soil surface (x = 0) and the actual rooting depth (L in dm, x = L).In the case of evaporation, the function h (x) describes the water extraction between soil surface and the maximum water extraction depth (E in dm, x = E) (Koitzsch and Günther, 1990) where x is the difference between soil surface and the corresponding soil layer i in dm and c is a parameter for the position of water extraction in the soil profile.The higher the value of c, the higher is the root water uptake in the upper soil layers.Here, it is assumed that c = 10 is for vegetation cover in the case of transpiration and c = 20 as well as the maximum water extraction depth E = 2 dm for bare soil in the case of evaporation (Koitzsch and Günther, 1990).Both models simulated interception using a simple single linear storage approach based on the maximum interception storage capacity depending on the type of vegetation.
In the soil water balance model, infiltration of precipitation or throughfall into the topsoil layer and runoff generation are calculated by a modified approach according to Holtan (1961).Actual water content and water fluxes of each soil layer are calculated using a multilayer capacity approach with nonlinear storage routing according to Glugla (1969).The water flux across the lower boundary of the soil profile is defined as groundwater recharge.The soil type of each hydrotope is described by a corresponding representative 1-D vertical soil profile down to 2 m depth.Based on texture and bulk density of each soil layer, obtained from these profiles, parameters like field capacity, and wilting point were evaluated according to AG Boden (2005).In the presence of groundwater in the soil profile, capillary rise is calculated according to AG Boden (2005) depending on soil texture, bulk density, distance from the soil layer to groundwater table, and soil water content in the layer.More information can be obtained from Wegehenkel (2002).

Modelling procedures
For the Ucker catchment, we used the 35 corresponding Remo-UBA-Grids covering the catchment (Fig. 1).Each grid consists of a time series with daily rates of precipitation, minimum and maximum air temperature, saturation deficit of air, wind speed and global radiation, split in a reference period from 1950-2000 and the climate change period from 2000-2100.In our study, we applied the hydrological model with the two different vegetation models separately for each grid without a spatially distributed calculation.The simulation runs for each grid were carried out with a daily time step starting at 1 January 1950 and ending at 31 December 2100.The initial soil water contents were set to field capacity.For each grid, we estimated the dominating soil type and land cover by GIS-procedures.Due to the actual restrictions in our dynamic vegetation model, we focused on regions and fields with agricultural land use in the catchment in our study.No additional correction of the precipitation of each grid was carried out.Agricultural land use in the Ucker catchment covers an area of 1101 km 2 and 47% of this area is located on a loamy soil with a clay content within 11-20% and a silt content within 15-21%.Therefore, we focused in our analysis on this soil type combined with agricultural land use.According to the FAO-classification from the year 1974, the soil is an Orthic Luvisol.The other soil physical properties used for the simulation runs can be obtained from Table 1.The crop rotation scheme for this soil type was estimated using the German soil quality index and management procedures such as organic farming or conventional farming (Kersebaum et al., 2003).In the case of the soil type Orthic Luvisol, the German soil quality index was within 35-45 and we assumed a conventional farming management.The corresponding crop rotation scheme consists of Maize -Winter wheat -Winter barley -Grass (Lolium Perenne) -Winter rye -Winter wheat -Winter barley -Grass (Lolium Perenne) -Sugar beet -Grass (Lolium Perenne) -Maize.This crop rotation scheme covers a time period of 8 years.More details about this estimation procedure of crop rotations can be obtained from Kersebaum et al. (2003).

Results and discussion
In the following analysis, we compared the reference period from 1950-2000 only with the climate change period from 2050-2100 to enable a better overview of the effects of climate change.As an example, the effects of climate change on precipitation and potential evapotranspiration ET p for one grid located in the centre of the catchment (location see Fig. 1) is illustrated in Fig. 2. Both models used the same approach according to Eq. ( 1) for calculating ET p .Cumulative precipitation for the climate change period decreased in an order of magnitude of 19% of the reference period and cumulative ET p as a general indicator for the climatological conditions for the Ucker catchment increased in an order of magnitude of 10% of the reference period (Fig. 2).For the total Ucker catchment, the decrease in precipitation was within 11-21% and the increase for ET p within 15-18% for the corresponding Remo-UBAgrids (Figs. 1-2).
For a more detailed analysis, we selected randomly two shorter time periods, each period covering 16 years with two complete crop rotations for the selected grid (Fig. 1).The selected time period from 1950-1966 presented the reference period from 1950-2000 and the selected time span from 2078-2094 corresponded to the climate change period from 2050-2100.and groundwater recharge (= C.GWR) in mm for the selected REMO-UBA-grid, calculated with the dynamic vegetation model for the selected reference period from 1950-1966(= Dyn. 1950-1966) ) and the selected climate change period from 2078-2094 (= Dyn.2078-2094) and simulated with the static vegetation model (= Stat. 1950-1966, Stat. 2078-2094).
For the selected climate change period from 2078-2094, simulated real evapotranspiration ET r showed a decrease in an order of magnitude of 22-23% of the selected reference period from 1950-1966 for both vegetation models (Fig. 3, upper graph).Due to the decrease in precipitation, cumulative groundwater recharge GWR in the selected climate change period simulated by the dynamic vegetation model showed a decrease in an order of magnitude of 9% of the selected reference period.In contrast to that, the static vegetation model calculated a higher decrease of cumulative GWR corresponding to 22% of the selected reference period (Fig. 3, lower graph).1950-1966(= Dyn. 1950-1966) ) and the selected climate change period from 2078-2094 (= Dyn.2078-2094) and simulated with the static vegetation model (= Stat. 1950-1966, Stat. 2078-2094).
These discrepancies in the results for cumulative GWR are due to the conceptual differences between both vegetation models regarding the calculation of the effects of soil water availability on crop growth.In the dynamic crop growth model, the feedback between transpiration, soil water availability and crop growth is simulated.However, in the static vegetation model, this feedback process is ignored.The ratio of ET r /ET p is a common indicator for crop water stress.In Fig. 4, simulated cumulative ET r /ET p for the selected reference period from 1950-1966 is compared with that for the selected climate change period from 2078-2094.The closer the ratio of ET r /ET p is to zero, the higher the water stress for the corresponding crop.Therefore, a decrease in the slope of cumulative ET r /ET p means an increase of crop water stress and, an increase of the slope a decrease of crop water stress (Fig. 4).Both vegetation models calculated an increase of crop water stress for the selected climate change period in an order of magnitude of 12-20% of the selected reference period (Fig. 4).However in comparison with the static vegetation model, the dynamic vegetation model simulated a higher increase of crop water stress with an amount of 20% (Fig. 4).
Soil covering degree SCD depends on the amount of above ground biomass and leaf area index LAI.For the selected climate change period, the dynamic vegetation model simulated a strong reduction of SCD for some crops in comparison with the selected reference period due to a corresponding increase of crop water stress .In contrast to that, the static vegetation model showed no significant impact of increased crop water stress on the time series of SCD (Fig. 5).
In both vegetation models, the conversion factor F (t) and SCD determine the amount and the division of ET r in transpiration and evaporation according to Eqs. ( 2) and (3) as well as the amount of interception.Due to the differences   (= Stat. 1950-1966, Stat. 2078-2094).
in the simulated time series of SCD and in the values of F (t) between both approaches, the static vegetation model calculated higher cumulative transpiration and interception as well as lower cumulative evaporation in comparison with the dynamic vegetation model for both, the selected climate change period and the selected control time period (Figs.5-6).Cumulative transpiration for the selected climate change period simulated by the static vegetation model decreased in an order of magnitude of 17% of the selected reference period (Fig. 6).The reduction of cumulative transpiration for the selected climate change period simulated by the dynamic vegetation model corresponded to 33% of the selected reference period (Fig. 6).Cumulative evaporation for the selected climate change period calculated by the static vegetation model decreases in an order of magnitude of 25% of the selected reference period, the dynamic vegetation model simulated a decrease of 8% (Fig. 6).Cumulative interception simulated by the static vegetation model for the selected climate change period decreased in an order of magnitude of 33% of the selected reference period, cumulative interception calculated by the dynamic vegetation model was reduced in an order of magnitude of 41% (Fig. 6).Due to lower rates of interception and corresponding higher amount of throughfall, the dynamic vegetation model simulated higher groundwater recharge rates for both periods in comparison with the static model (Figs. 3,(5)(6).

Conclusions
In our study, the comparison of the results obtained from the static vegetation model with those obtained from the dynamic vegetation model indicated the importance of taking into account vegetation dynamics and feedback between vegetation cover and soil moisture by transpiration in hydrological models especially if such models are used for the assessment of the effects of climate change on water balance components such as evapotranspiration and groundwater recharge.
Actually, there exist ecohydrological models with dynamic vegetation components such as the model SWIM (Krysanova et al., 1998) or the model SWAT (Arnold and Fohrer, 2000).In different studies, these models simulated vegetation dynamics with an acceptable precision (e.g.Post et al., 2007;Arabi et al., 2008).However in both models, the quality of vegetation modeling depends on the availability of a more or less large amount of crop specific input parameters.In addition, both models require a precise and detailed map of the actual regional distribution of agricultural crops such as winter cereals or sugar beet, of grassland and of forest types e.g.Oak-Beech forests or Scots-Pine including the age of the forest stands.Furthermore, crop specific input parameters required for the application of crop growth models are often difficult to identify, and in some cases, they have to be calibrated for the local conditions in the catchment.Information about agricultural and forest management practices such as fertilization in the catchment is also required.Such information is often not available in mesoscale catchments.This is the main problem of applying such vegetation dynamic simulation approaches integrated in models like SWIM or SWAT in mesoscale catchments.This leads to the conclusion, that in the case of limited data availability regarding e.g.areal crop distribution and agricultural management practices at the catchment scale, the incorporation of vegetation dynamics in hydrological models might need a more simpler model taking into account only some key ecosystem processes including e.g.vegetation establishment, growth, mortality, photosynthesis, phenology and physiognomy for all vegetation types with a low demand of input parameters (e.g.Arora, 2002;Gerten et al., 2004).However, such a simple approach is only driven by water and climate conditions and cannot simulate additional impacts on crop growth such as plant disease, increase in CO 2 and fertilization.

Fig. 1 .
Fig. 1.Location of the Ucker catchment in germany, Remo-UBAgrids, weather stations, regional distribution of soil type Orthic Luvisol and the selected grid for the analysis of the effects of climate change.

Fig. 2 .
Fig. 2. Comparison of cumulative precipitation (= Prc) and cumulative potential evapotranspiration (= ET p ) in mm for the reference period from 1950-2000 with those for the climate change period from 2050-2100, calculated for the selected REMO-UBA-grid.

Fig. 5 .
Fig. 5. Comparison of the time series of soil covering degree (= SCD) in % for the selected REMO-UBA-grid, calculated with the dynamic vegetation model for the selected reference period from 1950-1966 and for the selected climate change period from 2078-2094 (= SCD Dynamic, upper graph) and simulated with the static vegetation model (= SCD Static, lower graph), 16 year crop rotation, M= Maize is the starting point of each crop rotation.

Table 1 .
Selected soil physical properties of the OrthicLuvisol, Ucker catchment, field capacity FC and wilting point WP according to AG Boden (2005).