Numerical study of surface energy partitioning on the Tibetan plateau: comparative analysis of two biosphere models

The Tibetan Plateau is a critical region in the research of biosphere-atmosphere interactions on both regional and global scales due to its relation to Asian summer monsoon and El Nĩ no. The unique environment on the Plateau provides valuable information for the evaluation of the models’ surface energy partitioning associated with the summer monsoon. In this study, we investigated the surface energy partitioning on this important area through comparative analysis of two biosphere models constrained by the in-situ observation data. Indeed, the characteristics of the Plateau provide a unique opportunity to clarify the structural deficiencies of biosphere models as well as new insight into the surface energy partitioning on the Plateau. Our analysis showed that the observed inconsistency between the two biosphere models was mainly related to: 1) the parameterization for soil evaporation; 2) the way to deal with roughness lengths of momentum and scalars; and 3) the parameterization of subgrid velocity scale for aerodynamic conductance. Our study demonstrates that one should carefully interpret the modeling results on the Plateau especially during the pre-monsoon period.


Introduction
Soil-vegetation-atmosphere interaction on the Tibetan plateau is influential on energy and water cycles on both regional and global scales and thus the Plateau has been the main topics in various research areas.The Plateau's surface energy balance (SEB) particularly plays an important role in the Asian summer monsoon and global climate change and, in turn, its unique environment is vulnerable to climate Correspondence to: Jinkyu Hong (jkhong@yonsei.kr)change (e.g., Denman et al., 2007).Due to rapid changes in land cover and large population pressure for economic growth, the Tibetan Plateau has been undergoing significant environmental changes over the past several decades.Yet the lack of our understanding of underlying feedback mechanism among soil, vegetation and the atmosphere on the Plateau hinders the correct prediction of future climate and its impact on sustainable economic growth.To enhance our understanding of water and energy cycles in this critical area, numerical modeling studies have been conducted for decades (e.g., Peylin et al., 1997;Takayabu et al., 2001;Gao et al., 2004;Hong and Kim, 2008;Yang et al., 2008;van der Velde et al., 2009).Despite many pioneering studies based on recently developed biosphere models, large uncertainties still remain in the simulation of surface energy partitioning (Takayabu et al., 2001;Hong and Kim, 2008).We do not have lucid explanations on how the environmental conditions on the Plateau influence the land-atmosphere interactions compared with other low altitude areas and the impact of its unique environment on the physical parameterization in biosphere models.
With its high elevation, the Tibetan plateau has unique characteristics: 1) downward shortwave radiation is very large because of its location in a high elevation region (∼4000 m); 2) vegetation cover fraction is clustered with large areas of exposed soil; 3) daytime upward longwave radiation is much larger than downward longwave radiation and thus the Plateau is a heat source to the atmosphere (e.g., Chen et al., 1985;Li and Yanai, 1996); and 4) radiative coupling, which is expressed by difference between net radiation and isothermal net radiation (Raupach, 2001), is not negligible (Hong and Kim, 2008).Under these conditions, it would be a challenging task to simulate SEB using biosphere models based on our current understanding.Conversely, such a distinctive environment on the Plateau is unique enough to merit further investigation on models' performance and Published by Copernicus Publications on behalf of the European Geosciences Union.
parameterizations.Indeed, the characteristics of the Plateau provide a unique opportunity to clarify the structural deficiencies of biosphere models as well as new insight into the surface energy partitioning on the Plateau (e.g., Hong and Kim, 2008;Yang et al., 2009).
The objectives of this study are to better characterize the performance of two biosphere models to simulate surface energy and water fluxes on the Tibetan Plateau; and then to elucidate the characteristics of surface energy balance on the Plateau through comparative analysis of two biosphere models that are constrained by in-situ measurements.All acronyms and notations used in this study are explained in the Appendix A.

Field observation
The intensive field observation of the GEWEX Asian Monsoon Experiment (GAME)-Tibet was conducted to monitor and to understand energy and water cycles on regional scale on the Tibetan Plateau in 1998.For this study, we used the data collected at the BJ station in Naqu (31.37 • N; 91.90 • E, 4580 m above m.s.l.) from 30 May to 14 September, 1998.The site was flat with a fetch over 1 km depending on wind direction.The soil surface was sparsely covered with short grasses with an average canopy height of <0.05 m and leaf area index (LAI) of <0.5.The atmospheric pressure was on average about 580 hPa during the study period.
The maximum of 30-min averaged downward shortwave radiation was 1235 Wm −2 , and 30-min averaged downward longwave radiation varied from 188 to 388 Wm −2 during the simulation period.Average downward shortwave and longwave radiations during total simulation period were 250 Wm −2 and 310 Wm −2 respectively.During the simulation period, the total precipitation at the site was 388 mm over the study period (Fig. 1).Eddy-covariance system consisted of a three dimensional sonic anemometer (CSAT3, Campbell Scientific Inc., Logan, Utah, USA), krypton hygrometer (KH20, Campbell Scientific Inc., Logan, Utah, USA) and a fine-wire thermocouple and their measurement height was 2.9 m above the ground surface (Table 1).The sampling rate was 20 Hz and the data were processed every 30 min for calculating turbulence statistics.Other supporting measurements included radiative fluxes, soil heat fluxes, soil temperature, soil water content, and the profiles of other meteorological variables (i.e., wind speed, humidity, and temperature) (Fig. 1).Detailed descriptions for the site are available in Choi et al. (2004), Yang et al. (2008) and the GAME-Tibet website (http://monsoon.t.u-tokyo.ac.jp/tibet/).

Model description
To simulate water and energy cycles on the Plateau, we have used two representative land surface models (LSMs): Simple Biosphere Model 2 (SiB2) and Noah LSM (version 2.7) (Table 2).The SiB2 incorporates a photosynthesis-conductance model and calculates soil temperature through the forcerestore method (two layer model).The soil moisture is computed using the diffusion equation at the soil surface layer, root zone, and deep soil layer.For calculating radiative transfer within the canopy, SiB2 uses the two-stream approximation, and albedo of soil surface and vegetation depends on soil and vegetation properties as well as on wavelength of incident radiation.
The Noah LSM adopts the Zilitinkevich parameter to consider the difference of roughness length between momentum and heat (Zilitinkevich, 1995).It calculates both soil temperature and soil moisture content at the same depth, and the number of soil layers is changeable from 2 to 20.Albedo is prescribed by users on a monthly basis before running Noah LSM and this prescribed monthly albedo is linearly interpolated into daily values.
For the calculation of turbulence fluxes, both models use the Monin-Obukhov (MO) similarity.The main contrast between the two models is a way to calculate canopy conductance.Noah LSM adopts the empirical functions to limit maximum canopy conductance (g c ) proposed by Jarvis (1976).In SiB2, canopy conductance is computed from canopy photosynthesis-conductance model (Sellers et al., 1996b).SiB2 also requires aerodynamic parameters such as the zero-plane displacement height and roughness length for describing turbulent exchanges above and within canopy.These aerodynamic parameters were estimated using MOMPT program in this study (Sellers et al., 1996b).The numerical integration was done every 30 min using the observed driving meteorological variables and initial conditions for soil temperature and soil moisture were defined based on the field observations (Table 1).See Fig. 2 in Sellers et al. (1996b) and Fig. 1 in Chen and Dudhia (2001) for the schematic diagram of the model structures.

Input parameters and initial conditions
Initial values and phenological parameters for running the two models were taken from the data observed at the site.In this study, the types of vegetation and soil at the site were assigned to "Agricultural/C3 grassland" and "sandy loam" in the models and default vegetation/soil-type dependent parameters were adopted.No parameter calibration was al-lowed so that we scrutinized structural deficiencies of the models.Detailed values of soil/vegetation-type dependent parameters are available at Zobler (1986) and Sellers et al. (1996a) for SiB2 and "http://www.emc.ncep.noaa.gov/mmb/gcp/noahlsm/Noah LSM USERGUIDE 2.7.1.doc",Cosby et al. (1984), andvan der Velde et al. (2009) for Noah model.Based on the field observations, canopy top and bottom heights for SiB2 were set to 0.05 m and 0.005 m, respectively.Concurrently, rooting depth was set to 0.3 m and the depth of each soil layer in SiB2 (soil surface layer, root zone, and deep soil layer) was 0.07, 0.23, and 1.2 m, respectively.The depth of the first soil layer in SiB2 was obtained from the observed soil temperature profiles for the force-restore method (Arya, 2001).In Noah LSM the number of soil layers was set to four, and the depth of soil layer was 0.07, 0.23, 0.5, and 0.7 m, respectively.The observed monthly averaged surface albedo was pre-defined to match the simulated reflected shortwave radiation to the observed values for Noah LSM.However, for SiB2, the default albedo was used with vegetation and soil type as well as the different wavelength.

Roughness length
Roughness length was allowed to have temporal variability in SiB2 with seasonal march of summer monsoon on the Plateau and had an order of 0.001 m which is approximately 2% of the canopy height.In Noah LSM, by contrast, www.biogeosciences.net/7/557/2010/Biogeosciences, 7, 557-568, 2010 observed monthly averaged surface albedo was pre-defined to match the simulated reflected short-110 wave radiation to the observed values for Noah LSM.However, for SiB2, the default albedo was used with vegetation and soil type as well as the different wavelength.

Roughness length
Roughness length was allowed to have temporal variability in SiB2 with seasonal march of summer 115 monsoon on the Plateau and had an order of 0.001 m which is approximately 10 % of the canopy height.In Noah LSM, by contrast, roughness length should be fixed throughout the simulation period and was assigned to mean roughness length in this study.To examine the sensitivity of the simulated SEB to roughness length, the models were driven by two different roughness lengths (0.001 and 0.01 6 roughness length should be fixed throughout the simulation period and was assigned to mean roughness length in this study.To examine the sensitivity of the simulated SEB to roughness length, the models were driven by two different roughness lengths (0.001 and 0.01 m) (Fig. 2).
With an increase of the roughness length (i.e., EXPNR compared to EXPN), both sensible (H) and latent heat fluxes (LE) slightly increased on a monthly scale, whereas surface temperature decreased in Noah LSM (i.e., reduction of OLR).Consequently, outgoing longwave radiation (OLR) reduced and hence net radiation (R n ) increased in Noah LSM.Our interpretation is that the increase of roughness length enhanced turbulent fluxes and such enhanced turbulent fluxes reduced surface temperature.By contrast, SiB2 showed different responses to the increase of roughness length probably due to the lack of coherence among the aerodynamic parameters in the subroutine EXPSR.In the EXPSR of SiB2 simulation, roughness length was forced to be a constant irrespective of the variation of vegetation cover fraction and LAI.In general, our interpretations below do not depend on such variation of roughness lengths and therefore we discussed the results from EXPS and EXPN with the observation data in the following sections.

Friction velocity
Figure 3 demonstrates the simulated friction velocity (u * ) against the observed values.In general, the two models reproduced the observed friction velocity reasonably well.Nevetherless, two features are worth noting: (1) Noah LSM simulated larger u * than SiB2 particularly when u * was small, thus reproducing more comparable u * to the observation; and (2) Noah LSM was not successful to replicate u * < 0.07 ms −1 .
We noted that the u * values simulated by SiB2 were smaller than those by Noah LSM particularly during the premonsoon season and then tended to become comparable as the monsoon progressed (not shown here).Our analysis showed that such a difference was related to the parameterization of convection velocity (w * ) which is used in calculating u * .Typically, u * is computed as: where C M is the drag coefficient.Equation ( 1) is adopted in SiB2 with a minimum wind speed of 0.25 ms −1 .In Noah LSM, however, w * is added in computing u * with a minimum wind speed of 0.01 ms −1 and a minimum friction velocity of 0.07 ms −1 :  We noted that the u * values simulated by SiB2 were smaller than those by Noah LSM particularly during the pre-monsoon season and then tended to become comparable as the monsoon progressed (not shown here).Our analysis showed that such a difference was related to the parameterization of

Friction velocity
where w * is given as where β is an empirical parameter (=1.2) and C H is the turbulent exchange coefficient for heat.The convection velocity w * was proposed for stable numerical integration because mesoscale horizontal wind causes a finite effective friction velocity even in the limit of zero mean wind when local free convection prevails in the surface layer (Businger, 1973;Högstrom, 1996).To remedy these conditions in the model, there is the prescribed minimum wind speed (u ≥ 0.25 ms −1 ) in SiB2.However, Noah LSM assigns the minimum friction velocity (u * ≥ 0.07 ms −1 ) and also convection velocity for the numerically stable simulation in calm conditions, which explains the minimum u * manifested in Fig. 3.
In the nighttime conditions, the convection velocity turns off.In the daytime, however, the magnitude of w * was sometimes compatible to the mean wind speed on the Plateau when the wind was calm and vertical gradient of temperature was large, i.e., a typical condition in the pre-monsoon season on the Plateau (Fig. 4).Indeed, the convection velocity is proportional to H, and therefore u * differences between the two models were larger during the pre-monsoon period when H was dominant surface fluxes.Without the convection velocity parameterization, the drag coefficient for momentum in Noah LSM became similar to that in SiB2 (not shown).
In Eq. ( 3), we also note that by default, surface temperature was set to 270 K, and the atmospheric boundary layer (ABL) height 1000 m in Noah LSM.Therefore, the convection velocity was exaggerated in Eq. (3) when the surface temperature and ABL height on the Plateau were higher than those fixed values, especially during the daytime in the summer season.Over the ocean surface where the diurnal variation of surface temperature is weak, the constant ABL height of 1000 m may not induce substantial bias (e.g.Zeng et al., 1998).Also, such an overestimation of w * helps numerical integration to be stable in the atmosphere-vegetation coupled model without the bias of model results due to negative feedback (McNaughton and Jarvis, 1991;Kaimal and Finnigan, 1994).There is, however, no regulation process in the off-line simulation and the bias in w * clearly makes an impact when the biosphere models were driven using the observed data particularly on the Tibetan Plateau during the pre-monsoon season.
The w * parameterization was derived for the case when the time average of wind speed (u) was replaced by the speed of the time-averaged wind vector ( √ u 2 + v 2 ) for the stable time integration of numerical models (Beljaars, 1995;Mahrt and Sun, 1995).In the present study, the time averaged wind vector from a 3-dimensional sonic anemometer was used for driving the two models.When a biosphere model is driven by the time average of wind speed from a cup anemometer, however, it may be inappropriate to include w * in the aerodynamic bulk formula.Such differences of mean wind speed have not been clearly pointed out in most model studies but become important when the wind meandering is substantial.

Sensible heat fluxes
The two models overestimated H against the observation.In particular, H in Noah LSM was larger than in SiB2 (Fig. 2).Plausible explanations about the differences in H between the models are 1) the difference in simulated T s which is expressed by OLR and 2) the difference in aerodynamic conductance.
The sensible heat flux is proportional to surface temperature in the models and thus larger T s in Noah LSM had a direct influence on H. Our analysis revealed that this larger www.biogeosciences.net/7/557/2010/Biogeosciences, 7, 557-568, 2010 inute averaged mean wind speed and convection velocity calculated in Noah LSM before and oon period.
arameterization was derived for the case when the time average of wind speed (u) was the speed of the time-averaged wind vector ( √ u 2 + v 2 ) for the stable time integration l models (Beljaars, 1995;Mahrt and Sun, 1995).In the present study, the time averaged from a 3-dimensional sonic anemometer was used for driving the two models.When model is driven by the time average of wind speed from a cup anemometer, however, appropriate to include w * in the aerodynamic bulk formula.Such differences of mean have not been clearly pointed out in most model studies but become important when the ering is substantial.

le heat fluxes
dels overestimated H against the observation.In particular, H in Noah LSM was larger 2 (Fig. 2).Plausible explanations about the differences in H between the models are ence in simulated T s which is expressed by OLR and 2) the difference in aerodynamic e. T s in Noah LSM mainly resulted from the underestimation of soil evaporation.As will be shown later, the Noah LSM underestimated soil evaporation and the T s difference largely disappeared when the soil evaporation in the models was adjusted.Also, using the same models, Hong and Kim (2008) showed that inaccurate assignments of surface albedo and emissivity can make substantial bias in T s simulation between the two models on the Plateau.In addition, the overestimated H during the nighttime in Noah LSM can be attributed to the overestimation of aerodynamic conductance necessitating an adjustment of the minimum wind and friction velocities before testing the model.The inhomogeneity of transfer mechanisms between momentum and heat is another factor contributing to the difference between the drag coefficients for heat (C H ) and momentum (C M ).This difference results from the surface heterogeneity as well as the pressure terms in the momentum conservation equation, which are handled differently in the two models.The different transfer mechanism between momentum and heat has been one of the critical issues in bound-ary layer meteorology (e.g., Garratt, 1992;Brutsaert and Mawdsley, 1996;Chen et al., 1997;Lhomme et al., 2000;Schüttemeyer et al., 2008;Yang et al., 2008) and the difference is commonly expressed by the Zilitinkevich parameter, C (Zilitinkevich, 1995) or Stanton number, B: Equation ( 4) suggests that the reduction of C increases C H by increasing z ho .Surface temperature becomes lower, but R n and H increase with the reduction of C.
The parameter C is related to soil status, land cover type, and vegetation coverage.Its value is within the range from 0.1 to 0.4 (Chen et al., 1997) and closely related to the model performance.For example, Beljaars and Viterbo (1994) found that an overestimation of evaporation in the winter is reduced when z m0 /z h0 = 10.Hopwood (1995) concluded that z m0 /z h0 would be ∼80 for an inhomogeneous land surface.Chen et al. (1997) also pointed out that the effect of large C becomes more apparent when the soil gets drier and the skin temperature higher (i.e., similar to the environment of the Plateau).Furthermore, z m0 /z h0 is likely to become very small for sparse canopy because the roughness elements in sparse canopy absorb momentum efficiently.It is also reported that, the ratio of z m0 /z h0 shows diurnal variation with atmospheric conditions (e.g., Su et al., 2001;Yang et al., 2003Yang et al., , 2008)).
In Noah LSM, C was prescribed as 0.4.The parameterization of C is, however, not considered explicitly in SiB2 because aerodynamic resistance for heat is calculated from aerodynamic conductance for momentum, causing the differences in g h (≡ H ρC p 1 θ−θ s ) between the two models.Figure 5 shows the sensitivity of the simulated surface fluxes on C in Noah LSM.When C diminished from 0.4 (EXPN) to 0.1 (EXPN01), R n , H, and LE slightly increased, but T s (so OLR) decreased.The g h of Noah LSM was larger than that of SiB2 for a given roughness length and the difference in g h diminished with increasing C. In contrast, the reduction of g h increases the bias of R n between the models.The simulation of surface temperature, sensible heat flux, and soil heat flux on the Plateau has improved by considering the diurnal variation of C in the model (e.g., Yang et al., 2009).In our study, however, the alteration of C did not produce concurrent reductions in inter-model biases of all the surface flux components including R n and LE on a monthly scale.
The two biosphere models assume that the source/sink distributions between temperature and water vapor are the same and the turbulent transfer satisfies the MO similarity (i.e., uniform distribution of scalars and no impact of entrainment process).On the Plateau, however, the source/sink distribution of water vapor is not the same as that of heat due to seasonal variations associated with the Asian monsoon (Choi et al., 2004).Moreover, scalar fluxes can be enhanced from the prediction of the MO similarity because of the role of outer turbulence, that do not contribute to momentum fluxes but scalar fluxes in the surface layer (Mc-Naughton and Brunet, 2002;Hong et al., 2004;McNaughton et al., 2007).Hence, the models, which are based on the MO similarity, can underestimate scalar fluxes.Recently, Mc-Naughton et al. (2007) proposed the outer dissipation rate to remedy the MO similarity and further study should be done to resolve this issue in the future.

Latent heat fluxes
The two models overestimated soil water content φ but reproduced the observed diurnal and seasonal patterns fairly well (Fig. 6).With the bias in φ, the two models showed different partitioning of rainfall into evapotranspiration (Fig. 7).
To quantify the effect of soil moisture bias in the models and the observation, we conducted control simulations by forcing the observed φ into the models at each time step (i.e., reducing φ in the models).In this control simulation, R n and LE decreased but H increased as expected.Desipite the similar φ and smaller aerodynamic conductance, SiB2 produced consistently larger LE, for which different explanation was necessary.
LE was divided into soil evaporation (LE g ) and transpiration (LE c ). Figures 7 and 8 show that the two models showed similar LE c and the difference was originated mainly from LE g .It is also noted that the ratio of LE g to total LE was biased against the observation data (Fig. 8).We found that the difference in LE g resulted from the use of a different formula for computing LE g .Soil evaporation in Noah LSM is formulated as: The bare soil evaporation exponent, F is one of the sitedependent parameters and was set to a default value of 2.0 in this study.In SiB2, different formula is used: ) Sellers et al. (1992) included r soil to prevent the excessive soil evaporation rate and this resistance was estimated from the time series of station data from the FIFE (First ISLSCP Field Experiment) project using a series of inverse-mode runs with SiB2: Again, r soil is a site-specific parameter.In particular, R n and G were forced to match with observations (for energy budget closure) (Sellers et al., 1992) The δ 18 O isotope measurement at the same site suggested that approximately 27% of the precipitation be lost by evaporation from soil surface during monsoon (Tsujimura et al., 2001).For estimating the impact of soil evaporation submodels in the models, the two models have been re-simulated after adapting soil evaporation in such a way that the simulated LE g was same with the observed LE g .Such an adjustment produced an increased θ in SiB2 but a decrease in Noah LSM.This simulation result leads us to conclude that soil evaporation is a main factor for consistent simulation of SEB over the Plateau (Fig. 9).It is important to note that we adapted the soil evaporation using the soil water contents overestimated by the two models, suggesting the need for the improvement of soil hydrology in the biosphere models.Despite such adaptation, the simulated soil water contents were still larger than the observed, and therefore other components of water budget (e.g., runoff) should be carefully assessed.
It is also noteworthy that the two models reproduced the similar surface energy partitioning after adjusting soil evaporation and gave us the similar Bowen ratio (H/LE), compared to the surface energy partitioning observed by Choi et al. (2004).Choi et al. (2004) reported the observed seasonal variation of the energy partitioning at the same site in terms of the Bowen ratio.They showed that the Bowen ratio decreased rapidly after the onset of the monsoon and remained low until mid September.Then, the Bowen ratio began to increase with the withdrawal of the monsoon.In the simulations without any model calibration, two models produced inconsistent surface energy partitioning (Figs. 2 and  8a).After the adjustment of soil evaporation, the two models showed better agreement with the observation (Fig. 8b).Indeed, similarly to the observation, in the two models, sensible heat flux from the surface was the major source of heating over the plateau prior to the onset of the monsoon, whereas the latent heat released by condensation became the primary source with the onset of the monsoon.

Soil heat fluxes
The soil heat flux (G), exceeding >200 Wm −2 during the daytime in the pre-monsoon season, plays a substantial role in surface energy partitioning on the Plateau (e.g.Choi et al., 2004;Gao, 2005).In general, both models reproduced the observed diurnal variations as well as monthly means of G (Figs. 2 and 10).The monthly averaged energy residual in Noah LSM was about −2∼−4 Wm −2 , and the sum of residual and G in Noah LSM was nearly the same as G in SiB2.We note, however, that the two models simulated soil temperature and soil water contents differently and G in SiB2 was more sensitive to the radiative forcing than Noah LSM.Part of the reason may be explained by smaller soil heat capacity in SiB2.
Different formulations for soil heat capacity are used in the two models.In SiB2, soil heat capacity is calculated following Camilo and Schmugge (1981): In Noah LSM, soil heat capacity is calculated as: www.biogeosciences.net/7/557/2010/Biogeosciences, 7, 557-568, 2010 The two biosphere models assume that the source/sink distributions between temperature and water vapor are the same and the turbulent transfer satisfies the MO similarity (i.e., uniform distribution 230 of scalars and no impact of entrainment process).On the Plateau, however, the source/sink distribution of water vapor is not the same as that of heat due to seasonal variations associated with the Asian monsoon (Choi et al., 2004).Moreover, scalar fluxes can be enhanced from the prediction It is also noteworthy that the two models reproduced the similar surface energy partitioning after djusting soil evaporation and gave us the similar Bowen ratio (H/LE), compared to the surface nergy partitioning observed by Choi et al. (2004).Choi et al. (2004) reported the observed seasonal ariation of the energy partitioning at the same site in terms of the Bowen ratio.They showed that the owen ratio decreased rapidly after the onset of the monsoon and remained low until mid September.
hen, the Bowen ratio began to increase with the withdrawal of the monsoon.In the simulations ithout any model calibration, two models reproduced inconsistent surface energy partitioning (Fig. and 8a).After the adjustment of soil evaporation, the two models showed better agreement to the bservation (Fig. 8b).Indeed, similarly to the observation, in the two models, sensible heat flux rom the surface was the major source of heating over the plateau prior to the onset of the monsoon, Note that the observed G is based on Eq. ( 9).The C soil in SiB2 is larger (within 5%) than that in Noah LSM for a given soil moisture content.The two models, however, simulated θ differently.Consequently, soil heat capacity in SiB2 was 20% smaller than that in Noah LSM, thereby increasing temporal variations of G in SiB2.The SiB2 calculates G using the force-restore method whereas the Noah LSM calculates G using the diffusion equation independently of other energy budget components.If we neglect the energy transfer due to phase changes in water and snow, the force-restore method is formulated as: ntercomparison of the ratio of soil evaporation (LEg) to canopy transpiration (LEc) during the whole on period.
artitioning of evapotranspiration into soil evaporation and transpiration in the models.soil heat capacity in SiB2.
rent formulations for soil heat capacity are used in the two models.In SiB2, soil heat capacity lated following Camilo and Schmugge (1981): LSM, soil heat capacity is calculated as: at the observed G is based on Eqn. 9.The C soil in SiB2 is larger (within 5 %) than that LSM for a given soil moisture content.The two models, however, simulated θ differently.
14 ntly, soil heat capacity in SiB2 was 20 % smaller than that in Noah LSM, thereby increasral variations of G in SiB2.
2 calculates G using the force-restore method whereas the Noah LSM calculates G using ion equation independently of other energy budget components.If we neglect the energy ue to phase changes in water and snow, the force-Restore method is formulated as: mulation, soil is divided into two layers and G is calculated at each time step after com-, H and LE.
tically, the depth of soil surface layer in the force-restore method is half of the damping In this formulation, soil is divided into two layers and G is calculated at each time step after computing R n , H and LE.Theoretically, the depth of soil surface layer in the forcerestore method is half of the damping depth for the diurnal soil thermal wave (Garratt, 1992;Arya, 2001).The damping depth depends not only on soil physical properties, but also  on the soil moisture content.In reality, the damping depth varied from 0.05 to 0.14 m on the Plateau as the monsoon progressed (not shown).The depth of the first soil layer in SiB2 was, however, constant during the simulation period.Indeed, an increase of the depth of the soil surface layer results in an increase of the Bowen ratio.Based on the observed damping depth on the Plateau, the depth of the soil surface layer in SiB2 tended to partition more energy into LE during the monsoon period.Consequently, one should consider the variations in the depth of soil surface layer with time for long-term simulation.

Summary and conclusions
Two representative biosphere models (SiB2 and Noah LSM) were examined using the in-situ observation data on the Tibetan Plateau, and the models' performance to simulate surface energy balance was scrutinized in this study.Compared to the observation data, the two biosphere models did not provide convergent estimates of surface fluxes throughout the simulation period including the seasonal march of the summer monsoon.The biases between the models were related mainly to 1) different aerodynamic conductance due to the convection velocity parameterization; 2) inhomogeneity of turbulent transfer between momentum and scalars; and 3) different formulations of soil evaporation.The structural deficiencies of the model seen in this study can be manifested in simulating SEB in other altitudes and latitudes as well.
The 18 O stable isotope data provided critical information on validating and analyzing the simulated evapotranspiration by the biosphere models.The two models produced dissimilar surface energy partitioning without any calibration due to the propagation of uncertainties of several input parameters.After adjusting key parameters, the two model better reproduced the observed surface energy partitioning consistently.Our findings reaffirm the importance of careful interpretations of the modeling outputs due to uncertainties inherent in the models and the significance of in-situ field observations in modeling SEB on the Tibetan Plateau.For an accurate modeling of SEB on the Plateau, more attention should be given to retrieving information about the soil properties, which also emphasizes the difficulty in interpreting the modeling outputs on the Plateau without constraints that are based on field observations.

Fig. 2 .
Fig. 2. Monthly averaged surface flux at the site.EXPS and EXPN denote the simulation of SiB2 and Noah LSM with 0.001 m of roughness length, which is calculated from MOMOPT program.EXPSR and EXPNR denote the simulation of SiB2 and Noah LSM with 0.01 m of roughness length.

Fig. 2 .
Fig. 2. Monthly averaged surface flux at the site.EXPS and EXPN denote the simulation of SiB2 and Noah LSM with 0.001 m of roughness length, which is calculated from MOMOPT program.EXPSR and EXPNR denote the simulation of SiB2 and Noah LSM with 0.01 m of roughness length.

Fig. 3 .
Fig. 3. Comparison of simulated and observed 30-minute averaged friction velocities.Dash line is a one-to-one line.

Figure 3
Figure 3 demonstrates the simulated friction velocity (u * ) against the observed values.In general, the two models reproduced the observed friction velocity reasonably well.Nevetherless, two features are worth noting: (1) Noah LSM simulated larger u * than SiB2 particularly when u * was small, 135

140 7Fig. 3 .
Fig. 3. Comparison of simulated and observed 30-min averaged friction velocities.Dash line is a one-to-one line.

Fig. 4 .
Fig. 4. 30-min averaged mean wind speed and convection velocity calculated in Noah LSM before and during monsoon period.

Fig. 7 .
Fig. 7. Intercomparison of the ratio of soil evaporation (LEg) to canopy transpiration (LEc) during the who simulation period.

Fig. 7 .
Fig. 7. Intercomparison of the ratio of soil evaporation (LE g ) to canopy transpiration (LE c ) during the whole simulation period.

Fig. 8 .
Fig. 8. Partitioning of evapotranspiration into soil evaporation and transpiration in the models.

Fig. 9 .
Fig. 9.The simulated surface energy partitioning: the simulation results before (A) and after (B) adjusting soil evaporation during the whole simulation period.
heat of air e a atmospheric vapor pressure e s saturated atmospheric vapor pressure g gravity constant g c canopy conductance g h aerodynamic conductance for heat h soil relative humidity at soil surface k von Karman constant = 0.4 r d aerodynamic resistance between ground and canopy air surface r soil bare soil surface resistance t time (second) u streamwise wind speed u * friction velocity z m0 roughness length for momentum z h0 roughness length for heat w * convective velocity scale B Stanton number C Zilitinkevich parameter C air heat capacity of air C d heat capacity of deep soil layer in SiB2 C g heat capacity of surface soil layer in SiB2 C H turbulent exchange coefficient for heat C M drag coefficient for momentum C s heat capacity of mineral soil C soil heat capacity of soil C w heat capacity of water E p potential evaporation G soil heat flux H sensible heat flux LE latent heat flux LE c latent heat flux from vegetation LE g latent heat flux from soil Re Reynolds number R n net radiation T d deep soil temperature in SiB2 T g surface soil temperature in SiB2 T s surface temperature W g vegetation cover fraction W 1 soil wetness fraction ≡ φ/φ s θ mean potential temperature of air θ s mean potential temperature of ground surface φ soil water content φ s porosity φ th soil water content threshold in Noah LSM ρ air density τ d day length ≡ 86400 seconds

Table 1 .
The observed meteorological variables.

Table 2 .
Summary of the model structure and soil/canopy parameters of Noah and SiB2.