450 likes | 717 Vues
The Structural Design and Operational Behavior of a Specific SVAT * Model. The particular land surface model (LSM) examined here is the “Mosaic LSM”. Although this model has some unique features, its description should nevertheless give a sense for how typical SVAT models work. References:
E N D
The Structural Design and Operational Behavior of a Specific SVAT* Model The particular land surface model (LSM) examined here is the “Mosaic LSM”. Although this model has some unique features, its description should nevertheless give a sense for how typical SVAT models work. References: Koster, R., and M. Suarez, Modeling the land surface boundary in climate models as a composite of independent vegetation stands, J. Geophys. Res.,97, 2697-2715, 1992. Koster, R. and M. Suarez, Water and Energy Balance Calculations in the Mosaic LSM, NASA Tech. Memo. 104606, Vol. 9., 1996. *SVAT stands for “soil-vegetation-atmosphere transfer”. SVAT models include SiB and BATS.
MOSAIC LSM: OVERALL STRUCTURE Mosaic Strategy: Using vegetation maps, the heterogeneous vegetation cover within a grid cell is subdivided into a “mosaic” of “tiles”. Separate energy and water budgets are computed over each (relatively homogeneous) tile. The GCM atmosphere responds to the areally-weighted fluxes. Bare soil (9%) Deciduous Trees (35%) Needleleaf Trees (24%) Grassland (32%) TYPICAL TILE BREAKDOWN FOR A GCM LAND SURFACE GRID CELL
Percent coverage of vegetation type within grid cell
Structure of a Mosaic LSM tile: Water Balance precipitation evaporation + transpiration INTERCEPTION RESERVOIR throughfall surface runoff infiltration SURFACE LAYER ROOT ZONE LAYER soil moisture diffusion RECHARGE LAYER drainage
Structure of a Mosaic LSM tile: resistance network Sensible heat network Evaporation network Not shown on this diagram is the “zero resistance” associated with evaporation from the canopy interception reservoir.
Operations performed at each time step: time step n Turbulence subroutine computes Eo, Ho (and tenden- cies) from Tsold and eaold Update Seasonally- Varying Parameters Reflectances computed => Net shortwave, PAR flux LSM computes energy and water balances preprocessing LSM call
Seasonally-varying parameters include: “greenness fraction”: the fraction of vegetation leaves that are alive LAI: the leaf area index roughness length and other boundary layer parameters root length We prescribe values to these parameters, using one of two approaches: 1. Assign values based on vegetation (or soil) type and time of year. (This is a necessary approach for many parameters.) 2. Assign geographically (and seasonally) varying parameter values from maps derived, e.g., from remote sensing data. LAI from tables (=f(veg type)) LAI from satellite data
er Eo ra-old ea-old rc-eff es(Ts-old) Preliminary Turbulent Flux Calculation Without considering any energy or water balance requirement, we can compute“preliminary” values of evaporation (Eo) and sensible heat flux (Ho) based on the values of surface temperature (Ts-old) and canopy air vapor pressure (ea-old) determined in the previous time step. A function of Ts-old, ea-old, Tr, er, roughness, etc. Tr ra-old The surface temperature, Ts, is assumed to apply to the canopy air, as well. Ho Ts-old ea is the vapor pressure in the canopy air. We treat ea as a prognostic variable, keeping track of its value between time steps.
.622r ea-old - er Eo = ps ra-old Ts-old - Tr Ho = ra H H ra E E ra ra Cpr = = = = ra-old ea T T T ea ea T ea ea-old,Ts-old ea-old,Ts-old ea-old,Ts-old ea-old,Ts-old ea-old,Ts-old ea-old,Ts-old ea-old,Ts-old ea-old,Ts-old -.622r -.622r -.622r Ts-old - Tr ea-old - er ea-old - er Ts-old - Tr ps ps ps ra-old ra-old2 ra-old2 ra-old2 ra-old2 -Cpr -Cpr Under this framework, compute Eo, Ho, and their tendencies. As will be seen later, these will be necessary for the energy balance calculation. + Cpr + ra-old
Preliminary Reflectance Calculation Reflectances (and thus net shortwave radiation and photosynthetically active radiation [PAR]) are assumed not to be affected by the energy and water balance calculations, which means we can compute them ahead of time. Shortwave radiation is divided into four components: a) Visible direct radiation b) Visible diffuse radiation c) Near-infrared direct radiation d) Near-infrared diffuse radiation. We calculate a reflectance for each component using a simple empirical formula that approxi- mates the results of the full two-stream calculation. For full details, see NASA Tech. Memo. #104538 (1991). Note: in the Mosaic LSM, albedo is not a function of surface water content.
Now that the preliminary calculations are done, it’s time to call the Mosaic LSM itself. The parameter list includes inputs, updates, and outputs: INPUTS: - Vegetation type and time step length - GCM “weather”: rainfall, wind speed, vapor pressure in air, etc.) - Seasonally-varying parameters - Eo, Ho, and tendencies - Radiation quantities (with shortwave reduced by pre-calculated albedo) UPDATES (i.e., prognostic variables): - Surface/canopy temperature Ts - Deep soil temperature Td - Canopy vapor pressure ea - Water contents of three soil layers - Water content of interception reservoir - Snow water equivalent OUTPUTS/DIAGNOSTICS: - Evaporation E and sensible heat flux H - Surface runoff and drainage out of the column - Anything else that might me of interest
Dry fraction -1 1 1 ea rdry-eff ( ) + = rc rd + rsurf ea subcanopy aerodynamic resistance rd rc rdry-eff surface resistance es(Ts) es(Ts) canopy resistance (a function of environmental stress) rsurf es(Ts) INSIDE THE MOSAIC LSM First step: Compute rc-eff, a single surface resistance that accounts for all evaporation pathways (transpiration, bare soil evaporation, interception loss, snow evaporation). Assuming no snow, we assume that the tile is covered by a “wet fraction” (over which the interception reservoir is full) and a “dry fraction” (over which it is empty). Note similarity to electrical resistance network calculation Wet fraction + Dry fraction We find the single effective resistance reff (for the entire surface) that would give the same evaporation as the dry area evap- oration (computed with rdry-eff) added to potential evaporation from the wet area.
H E H E T ea ea T ea-old,Ts-old ea-old,Ts-old ea-old,Ts-old ea-old,Ts-old The energy balance calculation has two unknowns: DTsand Dea. It thus needs two equations. The first one has been seen before: Sw + Lw= Sw + Lw + H + lE + G basic energy balance In the Mosaic LSM, we assume: Snowmelt (M) is a special case, to be treated later. Emissivity =1, so that the upward longwave radiation = sT4 The ground heating, G, is composed of two terms: heating of the surface system (CpDTs/Dt) and a heat flux into the deep soil (assumed proportional to Ts - Td). Ho, Eo, and their “tendencies” are provided by the GCM as described above, so that we can assume H = Ho + DTs Dea + Dea DTs + E = Eo +
Assuming Ts = Ts-old + DTs and ea = ea-old +Dea, we can show that the energy balance equation reduces to: Qo = [ 4sTs-old3 + dH/dT + ldE/dT + Cp/Dt + b] DTs + [dH/dea + ldE/dea]Dea, where Qo = Sw - Sw + Lw - sTs-old4 - Ho - lEo - b(Ts-old -Td), and b = deep soil heat flux proportionality constant. Equation #1 How do we get the second equation? er Assume that evaporative flux from canopy air to reference level... E ra ea ...equals the flux from the saturated surface (within stomates, etc.) to the canopy air. (That is, the canopy air can’t “build up” moisture.) rc-eff E es(Ts)
In other words, Eo + (dE/dT) DTs + (dE/dea) Dea = (0.622r/ps) (es(Ts) - ea) / reff Flux from canopy air to reference level Flux from surface to canopy air Expanding, and neglecting 2nd order terms, gives Eo - (0.622r/ps) (es(Ts-old) - ea-old) / reff-old = (1/reff-old) [(0.622r/ps) des/dT - (dE/dT) reff-old - Eo(dreff/dT) ] DTs + (1/ reff-old ) [(-0.622r/ps) - (dE/dea) reff-old - Eo(dreff/dea) ] Dea Equation #2 The two equations are solved for DTs and Dea. Afterward, snowmelt is accounted for, if necessary. Note that the equations simplify in cases of dewfall or snow evaporation, for which we assume reff = 0.
Next: compute water fluxes (i.e., solve the various water balance equations). 1. Evaporation. With DTs and Dea computed, the evaporation rate for the time step is known. We remove this moisture from the interception reservoir, the surface soil layer, and the root zone soil layer in amounts consistent with the resistances. 2. Moisture transport between soil layers. We use a discretized version of Darcy’s law for unsaturated flow. (See water balance lecture. Some features differ; e.g., we use an “upstream” hydraulic conductivity.) Moisture flow from surface layer to root zone layer accounts slightly for subgrid heterogeneity. 3. Assign precipitation water to reservoirs. Assume a uniform precipitation depth within a prescribed fractional wetted area, and allow a fraction of this storm area to consist of previously wetted leaves. Surface runoff and infiltration are computed from resulting throughfall.
The Mosaic LSM, like any other SVAT model, has a drawback -- it requires “realistic” values for numerous parameters: Soil: Layer capacities Porosity Saturated soil matric potential Saturated soil hydraulic conductivity Soil pore size distribution index Bedrock slope Surface heat capacity Vegetation: Surface “type” (one of 10 generalized types) Leaf area index Greenness fraction Roughness height Vegetation height Unstressed canopy resistance parameters (6) Vapor pressure deficit stress parameter Temperature stress parameters (3) Leaf water potential stress parameters (5) Subcanoy aerodynamic resistance parameters (2) Other: Storm fractional area The vegetation type assigned to the tile defines the values used for most of these parameters Note that some of the parameters cannot be directly measured
The Mosaic LSM’s structure allows the breakdown of total evaporation by vegetation tile... … and the breakdown of each tile’s evaporation by component.
GCM P, radiation, Tair, etc. E, H, upward longwave LSM How do we evaluate the performance of such an LSM? “Online” approach: test GCM output against observations. Advantage: The coupling effects can be studied, and various sensitivity tests can be performed. Disadvantage: The model forcing (precipitation, radiation, etc.) can be wrong, so validating the land surface model can be very difficult. (“Garbage in -- Garbage out”) Example from GISS GCM/LSM: The Amazon river is poorly simulated, but we can’t tell if this is due to a bad LSM or poor precipitation from the GCM.
Better approach: Offline forcing (one-way coupling) Forcing Data Advantage: Land surface model can be driven with realistic atmospheric forcing, so that the impact of the LSM’s formulations on the surface fluxes can be isolated. Disadvantage: Deficient behavior of the LSM may seem small in offline tests but may grow (through feedback) in a coupled system. Thus, offline tests can’t get at all of the important aspects of a land surface model’s behavior. Output File P, radiation, Tair, etc. E, H, Rlw , diagnostics LSM PILPS model intercomparisons (to be discussed in a later lecture) have largely focused on such offline evaluations.
model obs The resulting mean seasonal cycle of runoff was compared to observa- tions. In this particular test, the Mosaic LSM did quite well. Mosaic LSM’s behavior in PILPS 2c (a study based in the Red-Arkansas River Basin). Forcing data covering several years for each of 61 1o X 1o grid cells in the Red Arkansas Basin were provided to participants.
Parameter values make a difference! In a recent study, it was found that the apparently poor behavior of the Mosaic LSM in an offline study using Oklahoma measurements was associated in part with an inaccurate setting of the ground heat capacity. Ground heating rates improve when the right heat capacity is used. Ground heating rates for mosaic are way off, throwing off the other fluxes. Robock et al., JGR, 108, D22, 8846, doi:10.1029/2002JD003245, 2003.
Coupled System Analysis Analysis, using Mosaic LSM, of what makes a SVAT model act differently from the standard Bucket model...
Sensitivity test: Addition of vapor pressure deficit stress Precipitation differences: with VPD stress minus w/o VPD stress Inclusion of vapor pressure deficit stress leads to large decreases in rainfall in some regions. Why? “Stomatal suicide” -- a serious positive feedback in the coupled system: Reduced humidity leads to higher VPD stress Higher VPD stress leads to reduced evaporation Reduced evaporation leads to precipitation Reduced evaporation leads to reduced humidity
Sensitivity tests: Removal of temperature stress; removal of interception loss mechanism Precipitation (top) and evaporation (bottom) differences: “with interception loss allowed” minus “w/o interception loss allowed” Precipitation differences: “with temperature stress” minus “w/o temperature stress”
Precipitation (top) and evaporation (bottom) differences: “standard formulation” minus “bucket-type formulation” Sensitivity tests: Coupling strategy A simulation was performed with a “pseudo-bucket”, one that used a “bucket-style” coupling to the atmosphere but was carefully controlled to reproduce the Mosaic LSM’s long-term surface energy budget in offline simulations. In the plots, large differences are seen in simulated evaporation and precipitation rates. These differences result strictly from feedbacks between the land and the atmosphere.
Main conclusions from coupled sensitivity analysis (Koster and Suarez, Advances in Water Resources, 17, 61-78, 1994.) 1. Of the environmental stresses that increase canopy resistance, -- temperature stress is not significant -- vapor pressure deficit stress is significant, partly due to feedback. 2. Of the main differences between the two model types, the presence of the interception reservoir in the SVAT model has the largest effect on evaporation rates. 3. The incorporation of a bucket model structure appears to have an effect on precipitation rates in the tropics and subtropics, perhaps due to the damping of diunal and synoptic-scale variability in land surface control. The differences, in any case, reflect land-atmosphere feedback. (In a later study, with the same LSM but a modified AGCM, the impact on the general circulation was found to be reduced.)
COMPUTER LAB: RUNNING A LAND SURFACE MODEL (A Retropective!) • This model is designed to simulate a tropical forest’s response to prescribed atmospheric forcing over a repeated full seasonal cycle. The relevant files are: • Model: gm_model.f (Includes driver; written in FORTRAN.) • Forcing file: TRF.DAT.30 (Includes rainfall rates, radiation forcing, etc., at a 30 minute time step over a full annual cycle. Model automatically interpolates to a 5 minute time step.) • Initialization file: input/lsm_input.dat (Includes parameter values to change for class experiments.) • How to run the model: • 1. Create input and output directories below the current directory. (This assumes a UNIX system.) • 2. Place lsm_input.dat in the input directory. • 3. Find a directory that can comfortably hold trf.dat.30.diur (1.4 Mb) • 4. Compile the program gm_model.f • 5. Modify the model parameters in lsm_input.dat as appropriate. • 6. Run the program. • 7. Four output files will be produced in the output directory: • mosaic.trf.mon.xxxx (4.5 Kb) • mosaic.trf.dat.xxxx (388 Kb) • mosaic.trf.tra.xxxx (12.9 Kb for 3-year run) • mosaic.trf.123.xxxx (291 Kb) • where xxxx is the label for the particular experiment. • 8. For new experiments, start at instruction 5.
INPUT FILE:/land/koster/pilps/TRF.DAT.30 This is the forcing data: modify path as necessary. VEGETATION IDENTIFIER:trf Leave as is EXPERIMENT IDENTIFIER:gp7 By changing this according to your own system of codes, you control the labeling of the output files of different experiments. TIME STEPS T.S. LENGTH DIAGS 1ST FORCING ALAT 534529 300. 2880 0 -3. 534529 = (365x3 + 31) x 24 x 12 + 1 = # of time steps in 3 years + 1 January + 1 time step. 300 = number of seconds in the 5 minute time step. DIAGS, 1ST FORCING, ALAT do not need to be changed. NUMBER OF TILES: 1 TYPE FRACTION 1 1.0 Type 1 = tropical forest Fraction = 1 means a homogeneous cover INITIALIZATION: TC TD TA TM 300.0 300.0 300.0 300.0 TC = Initial canopy temperature TD = Initial deep soil temperature TA = Initial near-surface atmospheric temperature TM = Initial assumed first forcing temperature WWW(1) WWW(2) WWW(3) CAPAC SNOW 0.5000 0.5000 0.5000 0.5 0. WWW(i) = Initial degree of saturation in soil layer i CAPAC = Initial fraction of interception reservoir filled SNOW = Initial snow amount
EXPERIMENT 1 HEAT CAPACITY WATER CAPACITY FACTOR TURBULENCE FLAG 70000. 1. 0 Heat capacity is in J/oK. If water capacity factor is 0.5, then the default capacity is halved; if it is 2, then the default capacity is doubled, etc. Turbulence flag: you won’t need this. EXPERIMENT 2 INTERCEPTION PARAMETER PRECIP. FACTOR 1. 1. Interception parameter: you won’t need this. Precip. factor: factor by which to multiply all precipitation forcing. EXPERIMENT 3 ALBFIX RGHFIX STOFIX 0 0 0 ALBFIX: If this is 1, you are using tropical forest albedo. RGHFIX: If this is 1, you are using tropical forest roughness heights STOFIX: If this is 1, you are using tropical forest water holding capacities. EXPERIMENT 4 FRAC. WET PRCP CORRELATION 0.3 0. FRAC. WET: The assumed fractional coverage of a storm; equivalent here to the assumed probability that a rainfall event will be applied to the land surface model. PRCP CORRELATION: Imposed time-step-to-time-step autocorrelation of precipitation events.
EXPERIMENT 1: CHANGE IN MODEL PARAMETERS • Background: • The heat capacity of the soil surface has an important effect on the land surface model’s surface energy budget calculations. Presumably, the higher the heat capacity, the more slowly the surface temperature will change under a given forcing, leading to a smaller amplitude of the diurnal temperature cycle. This could have profound effects on the annual energy balance. • The water holding capacity of the soil has an important effect on the annual water balance and thus on the annual energy balance. A larger water holding capacity, for example, means that high precipitation rates in the spring can more easily lead to high evaporation rates during a subsequent dry summer. • Possible experiments: • .Modify the heat capacity. You may have to modify it by an order of magnitude or so to see significant effect on the energy budget terms. • .Modify the water capacity factor. For starters, try 0.5 and 2. • Questions to answer (choose 1) • 1. How does varying the heat capacity affect the diurnal energy balance, in particular the amplitude of the diurnal temperature cycle? How large does the change have to be to see an effect? Is the effect in the expected direction? • 2. How does varying the heat capacity affect the annual energy balance? • 3. How does varying the water holding capacity affect the diurnal and annual energy and water budgets? Does a higher capacity imply a larger annual evaporation?
time of day time of day (Recall from earlier lecture:) The choice of the heat capacity can have a major impact on the surface energy balance. Low heat capacity case High heat capacity case -- Heat capacity might, for example, be chosen so that it represents the depth to which the diurnal temperature wave is felt in the soil. -- Note that heat capacity increases with water content. Incorporating this effect correctly can complicate your energy balance calculations.
Effect of Heat Capacity Change Diurnal Cycle of Surface Temperature: Black: Cp = 70000 J/kg (control) Purple: Cp = 700000 J/kg Yellow: Cp = 7000 J/kg
Effect of Heat Capacity Change Annual Cycle of Surface Temperature: Black: Cp = 70000 J/kg (control) Purple: Cp = 700000 J/kg Yellow: Cp = 7000 J/kg
Effect of Water Holding Capacity Recall from earlier lecture:
Effect of Water Holding Capacity Change control value x1 control value x2 control value x5 control value x 1/2
EXPERIMENT 2: CHANGE IN MODEL INITIALIZATION Background: All models require a “spin-up” period to remove the effects of initialization. In other words, the initial conditions imposed in a model may be inconsistent with the preferred model state, and this inconsistency may lead to energy and water budget terms that are unrealistic – they reflect the inappropriate initial conditions imposed rather than the model parameterizations or the atmospheric forcing. The length of the spin-up period is a function of the model (in particular its heat and moisture capacities) and the forcing. Possible experiments: Initialize the soil moisture reservoirs to complete saturation: set WWW(1), WWW(2), and WWW(3) to 1. Initialize the soil moisture reservoirs to be completely dry: set WWW(1), WWW(2), and WWW(3) to 0.0001. Initialize the soil moisture reservoirs to be completely dry, and double the water holding capacity: set WWW(1), WWW(2), and WWW(3) to 0.0001, and set the “water capacity factor” (from experiment 1) to 2. Complete drydown. Set WWW(1), WWW(2), and WWW(3) to 1, and set the “precip. factor” to 0. (This turns off all precipitation.) Note: for these experiments, you may want to increase the number of time steps. (You won’t know if you need to until you run them.) If n is the number of years you want the model to run, set the # of time steps to [(365*n)+31)]*24*12+1. Questions to answer (Choose 1): 1. How does the transient model response differ in the drydown and wet-up simulations (1 & 2)? 2. How does doubling the water holding capacity affect the wet-up period? 3. How long does complete drydown take (simulation 4)? Is equilibrium ever really achieved? Can you define a time scale for the drydown?
First 5 years First year, including equilibrium cycle
EXPERIMENT 3: CHANGE IN MODEL BOUNDARY CONDITIONS • Background: • GCM deforestation experiments have examined how replacing the Amazon’s forest with grassland can affect the regional climate. In a land surface model, forest and grassland are distinguished from each other only by the values used for various parameters. The experiments below examine “deforestation” in an offline environment. (Of course, deforestation effects in a fully coupled GCM environment may be different.) • Possible experiments: • .Perform a control simulation, using TYPE =1 (tropical forest). • .Replace the tropical forest with grassland: set TYPE=4. • .Replace the tropical forest with grassland, but maintain tropical forest albedo: set TYPE=4 and ALBFIX=1. • .Replace the tropical forest with grassland, but maintain tropical forest roughness: set TYPE=4 and RGHFIX=1. • .Replace the tropical forest with grassland, but maintain tropical forest water holding capacity: set TYPE=4 and STOFIX=1. • .Replace the tropical forest with grassland, but maintain tropical forest albedo, surface roughness, and water holding capacity: set TYPE=4, ALBFIX=1, RGHFIX=1, and STOFIX=1. • Questions to answer (choose 1) • 1. What is the effect of deforestation on the annual energy and water budget? What effect does it have on diurnal cycles? • 2. How do albedo change, roughness change, and storage change contribute to the tropical forest / grassland differences? Which effect is largest? • 3. Are the impacts of albedo change, roughness change, and storage change linear? E.g., do the changes induced by these three parameters alone add up to the changes seen in simulation 6?
Main effect of “deforestation”: Black = forest, purple = grassland Black = forest purple = grassland yellow = grassland with forest albedo Black = forest purple = grassland yellow = grassland with forest storage
EXPERIMENT 4: CHANGE IN MODEL FORCING • Background: • The precipitation forcing, which comes from a GCM, need not be assumed to fall uniformly within the GCM’s grid cell area. If the typical areal storm coverage is, say, only half the grid cell’s area, then one can consider an alternative interpretation: that whenever the GCM provides precipitation for a grid cell, the probability that it occurs at a given point within the cell is ½, and when it does occur there, the GCM’s precipitation intensity is doubled. A further consideration is the temporal autocorrelation of storm events, i.e., the probability that a point gets wet during one time step given that it was wetted in the previous time step. • Possible experiments: • .Perform a control simulation. • .Perform simulations that assume a fractional storm coverage of ranging from .1 to .9 (i.e., set FRAC. WET = x, where x ranges from .1 to .9). • .Perform simulations that assume a fractional storm coverage of .1 and a time step to time step autocorrelation that ranges from .1 to .9. (i.e., set FRAC. WET=0.5 and PRCP CORRELATION=x, where x ranges from .1 to .9). • Questions to answer (Choose 1): • 1. How does runoff ratio (runoff / precipitation) change with the assumed fractional coverage? • 2. How do runoff ratios change when temporal autocorrelations are included?
Storm area = 0.01 Effect of fractional storm area Effect of correlation in storm position. (Storm area = .1) Correlation = 0.9 Correlation = 0.