Land surface spinup for episodic modeling

Soil moisture strongly controls the surface fluxes in mesoscale numerical models, and thereby influences the boundary layer structure. Proper initialization of soil moisture is therefore critical for faithful simulations. In many applications, such as air quality or process studies, the model is run for short, discrete periods (a day to a month). This paper describes one method for soil initialization in these cases – self-spinup. In self-spinup, the model is initialized with a coarse-resolution operational model or reanalysis output, and run for a month, cycling its own soil variables. This allows the soil variables to develop appropriate spatial variability, and may improve the actual values. The month (or other period) can be run more than once if needed. The case shown is for the Boundary Layer Late Afternoon and Sunset Turbulence experiment, conducted in France in 2011. Self-spinup adds spatial variability, which improves the representation of soil moisture patterns around the experiment location, which is quite near the Pyrenees Mountains. The self-spinup also corrects a wet bias in the large-scale analysis. The overall result is a much-improved simulation of boundary layer structure, evaluated by comparison with soundings from the field site. Self-spinup is not recommended as a substitute for multiyear spinup with an offline land data assimilation system in circumstances where the data sets required for such spinup are available at the required resolution. Self-spinup may fail if the modeled precipitation is poorly simulated. It is an expedient for cases when resources are not available to allow a better method to be used.


Introduction
Episodic runs of atmospheric mesoscale numerical models are commonly used for air quality and research process studies, among other applications.Modeling for episodes of one to several days often uses coarse-resolution analysis as initial conditions and interpolates them to fine grids (1-6 km).This type of modeling presents special challenges not faced by continuous operational systems.One of these is the proper initialization of soil moisture, which is critical to the simulation of surface fluxes, and consequently of the atmospheric boundary layer, a first-order control on simulated pollutant concentrations.Here, we describe a case where the usual approach fails, and demonstrate a self-spinup method to improve the results.
The preferred method for soil initialization is to "spin up" the soil for several years with an offline land data assimilation system (Koster et al., 2010;LeMone et al., 2008;Chen et al., 2007;Kumar et al., 2006).This must in general be done at the resolution to be used by the target model (Santanello et al., 2011), and therefore requires highly resolved atmospheric fields and precipitation analyses.In case such a system and data are not available, spinning up the soil with the atmospheric model itself can improve results.This "open loop" approach was shown by Di Giuseppe et al. (2011) to improve the results of 3-month simulations.
The purpose of soil spinup is to provide soil moisture and temperature that are appropriate for the target model.Because the mesoscale target model likely has different soil levels, soil properties, vegetation types, and tuning constants than the larger-scale model, the appropriate soil moisture and temperature are not necessarily the same in the two models.The criterion for appropriateness is that the target model produces correct surface fluxes.Even if the physics of the two models is identical, the greater spatial variability in the finer-grid target model requires some spinup.Our hypotheses, then, are that (1) self-spinup increases spatial variation of soil moisture and temperature (downscaling), and (2) selfspinup removes biases overall and/or in specific locations, land uses, or soil types.
The example used here is from the Boundary Layer Late Afternoon and Sunset Turbulence (BLLAST) field campaign conducted in June and July 2011 in France (Lothon et al., 2014).Extensive measurements were taken from 14 June to 8 July 2011 in and around Lannemezan, southern France, near to the Pyrenees Mountains.The campaign site extended over an area of approximately 100 km 2 covered with heterogeneous vegetation including grass, corn, moor or forest.The site is quite near the Pyrenees and thus poses a challenge for modeling, requiring reasonably fine resolution to avoid confusing the local terrain elevation and land uses with the mountains.The project objective is to understand the main physical processes controlling the afternoon transition of the boundary layer from convective instability to nocturnal stability.Mesoscale models will be used, evaluated, and improved as part of the project.However, before the modeled afternoon transition can be compared with observations, the preceding conditions must match the observations to a reasonable degree.This study therefore emphasizes getting the model into a state where its boundary layer structure looks like the observations in the afternoon.Soil initialization is the key factor in achieving that agreement.

Modeling strategy
The Advanced Research core of the Weather Research and Forecasting model (WRF-ARW) was used in this study.Nested 9 and 3 km horizontal grids were used, each 100×100 points, centered roughly at the experiment site.Two-way nesting was used so that the 3 km solution influenced the outer grid.Fifty vertical levels strongly compressed near the surface were used (lowest eta levels 0.999, 0.998, etc. incrementing by 0.001 up to 0.990).Physics options were WRF single-moment three-class microphysics, RRTM-G shortand longwave radiation, and Kain-Fritsch convection on the 9 km grid only.References and details can be found in Skamarock et al. (2008).The outer 9 km domain covers roughly −5 to 6 degrees longitude and 38.5 to 47 degrees latitude, including much of France, Spain, the Bay of Biscay, and the northeastern Mediterranean Sea.All results shown are from the inner 3 km domain only, and its extent can be seen in the figures.
Vertical mixing was handled by the MYJ (Janjic, 2002) planetary boundary layer (PBL) scheme and its matching sur-face layer.Mixing in MYJ depends only on the local gradient of buoyancy, and is parameterized by prognostic turbulent kinetic energy and a length scale.A number of studies, for example (Shin and Hong, 2011;Angevine et al., 2012;Garcia-Diez et al., 2013), have found that the scheme mixes too little and entrains too little, resulting in boundary layers that are shallower, cooler, and moister than observations.However, MYJ is widely used, and other schemes also have biases.Moreover, it is difficult to distinguish biases in the PBL schemes from other compensating errors even in carefully designed studies (Garcia-Diez et al., 2013).
This study primarily addresses how the land surface in the model should best be initialized.The Noah land surface model (LSM) (Chen and Dudhia, 2001) is used.The Noah LSM has four soil layers (thicknesses 0.1, 0.3, 0.6, 1 m).The soil levels are the same for moisture and temperature and are constant throughout the domain.The atmosphere is initialized each day at 00:00 UTC with the ECMWF Interim Reanalysis (ERA-Interim) (Dee et al., 2011;Albergel et al., 2012), and the lateral boundary conditions for the outer (9 km) domain come from the same analysis.The ERA-Interim fields are on a 0.75 × 0.75 degree grid.
We chose 1 month as the basic spinup period.This is admittedly somewhat arbitrary.The optimal period and length of spinup cannot be specified generally a priori (Yang et al., 2011).In this case, the month of June corresponded approximately to the BLLAST experimental period.During June, the land state including soil moisture and snow in the higher elevations changes quickly, so including earlier periods would probably be less beneficial.The deep soil is weakly coupled to the surface, so its state does not contribute very much to short-term simulations, and spinning up its state requires a longer spinup period.Ultimately the wisdom of the choices we made here must be judged by the result, which for our objectives was satisfactory.
Three model runs were made: one with initialization only from the ERA-Interim analysis, one with soil moisture and temperature cycled once through the month of June, and one with two cycles.For "uncycled" or "ERA" runs, the soil moisture and temperature were initialized from ERA-Interim at 00:00 UTC each day.For "cycle 1", the soil moisture and temperature were initialized from ERA-Interim at 00:00 UTC on 1 June, then self-cycled for 30 days.That is, each day a run was started at 00:00 UTC with the soil temperature and moisture from the 24 h forecast of the previous day's run."Cycle 2" started at 00:00 UTC on 1 June with WRF soil moisture and temperature from cycle 1 at 00:00 UTC on 1 July.

Results
Runs with simple interpolation of the ERA-Interim soil variables into WRF ("uncycled" or "ERA" runs) had a drastic cool, moist bias compared to soundings at the experiment site (Fig. 1, red lines).The figure shows the mean difference between each model run and soundings taken at approximately 17:00 UTC on 7 days.The days included were intensive observation periods (IOPs) of the campaign, relatively undisturbed by synoptic systems, and were chosen for various of the other analyses of the BLLAST project (Lothon et al., 2014).Soundings at 17:00 UTC were used on the premise that analysis of the afternoon transition depends on correctly simulating midday conditions.The bias was largest below ∼ 500 m, indicating that the land surface is a likely source of this bias, most likely due to a soil that is too moist.This results in too much of the incoming solar radiation being allocated to evaporation and too little to sensible heating, reducing PBL warming, growth and entrainment.After cycling the soil for most of a month (cycle 1) or almost 2 months (cycle 2), the results were much improved.The Noah LSM required two cycles for optimum performance.This solved the problem for further studies with mesoscale models in BLLAST, but we wanted to understand the mechanism more thoroughly.Figure 2 shows soil moisture in the second layer of the inner domain for the different simulations.Throughout this paper "soil moisture" means soil moisture content as a volume ratio.We show the second layer for clarity, because the thin first layer has very strong daily cycles that obscure the trends we want to highlight.Due to the coarse ERA-Interim data resolution (0.75 × 0.75 degree grid), the Pyrenees are severely under-resolved (Fig. 2a), and after interpolation to the WRF grid, the Lannemezan site is in a broad area of moist soil representative of the mountains.After only 1 day of simulation, the spatial detail of soil moisture is visibly improved.After 30 days of free running (one cycle), the detail is fully resolved to the WRF grid.Moisture at and near Lannemezan is somewhat reduced.After two cycles (60 days) the overall moisture is somewhat further reduced, but the level of detail is similar.Figure 3 shows soil temperature in the same format as Fig. 2. The initial soil temperature has more detail than moisture, because the WRF initialization modifies the temperature according to the terrain height.Soil temperature in the domain as a whole increases throughout the cycling runs, as would be expected for June, but the temperature in the immediate neighborhood of Lannemezan changes rather little.
The temporal and spatial variation of soil moisture are shown in Fig. 4 (note that each day of June is shown twice on the horizontal axis to illustrate the evolution through two cycles).The second soil layer (0.1-0.4 m) is shown in order to avoid the strong diurnal cycle in the top layer.The ERA-Interim moisture dries somewhat after approximately 8 June, and its spatial variation (standard deviation) has no trend.In the cycled WRF runs, moisture decreases except in early June when it rains.Rain causes the spatial variation to increase in WRF.
Soil moisture itself is difficult to evaluate because observations are sparse.In the BLLAST domain, observations from the SMOSMANIA network can be used (Albergel et al., 2008;Calvet et al., 2007).Figure 5 shows the locations of eight sites within the inner WRF domain that had full data for June 2011.Figure 6 shows the average soil moisture from those sites.The distribution of sites is not uniform, and there is a great deal of site-to-site variation, which is not shown.Observations are a weighted average over all soil levels, the deepest of which is at 0.3 m.The model soil moisture is taken from the nearest grid point to each site.The observations show an abrupt increase in moisture on 6-7 June due to precipitation, which is also present to differing degrees in each of the models.ERA-Interim is more moist and dries more slowly than the observations.This bias was also found by Greve et al. (2013) and Albergel et al. (2012), and discussed by De Rosnay et al. (2013).Cycling WRF runs respond less to the precipitation episodes on 6-7 and 12 June.After 15 June, the drying rate with Noah is similar to the observations, but cycle 2 starts lower, so the two Noah cycles bracket the observations.
For purposes of atmospheric modeling as opposed to hydrology, the important output of the LSM is not the soil mois- ture itself but the partitioning between sensible and latent heat fluxes (Santanello et al., 2013).Figure 7 shows the sensible heat flux from the Noah LSM at 14:00 UTC on the same days as Fig. 2. The heat flux from the Noah LSM has detail immediately despite the smooth soil moisture because the flux varies according to the soil and vegetation types, which have grid-scale variations.At the end of the first simulated month of June (day 30), the pattern is quite different, particularly in the northeast quadrant of the domain where the flux is larger.After another cycle, day 60 has a similar spatial pattern but magnitude has increased further.At and in the immediate neighborhood of Lannemezan, the flux is nearly unchanged.
Precipitation in the model determines whether self-spinup succeeds or fails.After initialization, modeled precipitation is the only source of water to the soil.It must be at least approximately correct in amount and spatial distribution.Figure 8 shows analyzed precipitation for the part of the domain that falls within France, and modeled precipitation, as totals  for the month of June.The modeled precipitation is for the first cycle.The analysis is from the SAFRAN system (Habets et al., 2008;Quintana-Segui et al., 2008).It is a dense but not regularly gridded analysis, but it has been plotted with large enough pixels to fill the space and enable visual comparison with the model output.Correspondence between the analysis and the simulations is not perfect, but the overall amount and distribution is similar.WRF has too little rain overall and especially in the western part of the domain.Note that the amounts over the Pyrenees have been allowed to saturate on the color scale so that the smaller amounts away from the mountains can be seen.The modeled precipitation in cycle 2 is not shown, but differs little from cycle 1.Noah produces 2-3 % less precipitation in cycle 2 than in cycle 1.We also  note that Lannemezan (red X) is in a region of strong gradients in precipitation as well as soil moisture and temperature.
How well did self-spinup work to achieve our objective of simulating reasonable conditions for further research on the BLLAST project? Figure 9 shows potential temperature profiles from frequent soundings (Legain et al., 2013) at experiment site 2 on 25 June, a day identified for further analysis.Noah cycle 2 has nearly no temperature bias in the boundary layer, although its boundary layer height is somewhat low (about 20 %).This confirms the composite result in Fig. 1.The boundary layer is deeper and warmer in cycle 2 even though the heat flux at the site on that particular day is about the same.The same flux leads to a warmer and deeper boundary layer because the air upwind and the local soil are both warmer.On this basis we have chosen the Noah cycle 2 soil variables as the basis for further work on mesoscale modeling for BLLAST.

Conclusions
We have shown that spinning up the soil in a mesoscale model by running the model itself with cycling soil variables can improve simulations over the default method of simply interpolating soil data from a coarser analysis.We call the method "self-spinup"; it is equivalent to the "open loop" method of Di Giuseppe et al. (2011).The somewhat arbitrary choice of a 1-month spinup period (repeated twice) is vindicated by the results.
Self-spinup increases spatial variation of soil variables; that is, it downscales the soil data to the finer grid.It removes biases overall, in this case correcting a cool and moist bias.Evaluation against soil moisture data throughout the WRF domain indicates that the bias existed in the soil moisture itself, not primarily in other aspects of the model.In general, however, even if simulated soil moisture perfectly matches observations the atmospheric simulation may still be degraded by model errors, as implied by de Rosnay et al. (2013) and by Hacker and Angevine (2012).
A 1-day simulation initialized with ERA-Interim soil moisture reduced by 33 % at each point, based on Fig. 6, produced similar soundings to the cycled result in the early afternoon, but not as good in the late afternoon (not shown).Given a similar comparison of simulated to observed soil moisture, users could consider this simpler technique.
The surface heat and moisture fluxes are the primary control on boundary layer structure.However, that control is not necessarily exercised in the local column or at the time of evaluation.It is interesting to note that the heat fluxes in the BLLAST experiment area vary much less with spinup than they do for the domain as a whole, and yet the resulting boundary layer profiles change (improve) substantially.This indicates that the improvement is coming from parts of the domain outside the immediate area, and is probably generally true for any complex terrain situation.Experiment designs should take into account the need for measurements of key quantities over mesoscale domains (at least).
Self-spinup cannot be recommended if spinup with a land data assimilation system can be done with the skills, resources, and (above all) data that are available.The bias reduction found here could be fortuitous.We have no way of knowing whether self-spinup would have corrected a dry bias, for example.The optimal period and length of spinup cannot be specified generally (Yang et al., 2011).However, most episodic modeling experiments are conducted after the fact, so improved performance can be verified with data.Precipitation from the mesoscale model is the most likely source of errors in self-spinup.Here we found that the precipitation the model produced was broadly reasonable and apparently sufficient to improve the results.

Figure 2 .
Figure 2. Soil moisture (dimensionless) in the second layer on the inner domain at different times in the cycling experiments.Initial time is 00:00 UTC on 1 June 2011.The red X is the Lannemezan measurement site.Elevation contours at 200 and 1000 m a.s.l. are plotted.

Figure 3 .
Figure 3. Soil temperature (K) in the second layer on the inner domain at different times in the cycling experiments.Initial time is 00:00 UTC on 1 June 2011.The red X is the Lannemezan measurement site.Elevation contours are plotted at 200 and 1000 m a.s.l.

Figure 4 .
Figure 4. Soil moisture mean (top) and spatial standard deviation (bottom) during June 2011.Each day of June appears twice on the horizontal axis to show the evolution through two cycles.Color scheme: red is ERA-Interim, blue Noah cycle 1, and green Noah cycle 2. Solid lines are averages over whole inner domain; dashed lines are averages over 10×10 grid points extending northeast from the experiment site.

Figure 5 .
Figure 5. Terrain height on inner WRF domain with the eight SMOSMANIA sites used in Fig. 6 plotted as red X marks.Color bar shows the terrain height scale in m a.s.l.

Figure 6 .
Figure 6.Average soil moisture for eight sites in the SMOSMANIA network.Each day of June appears twice on the horizontal axis to show the evolution through two cycles.Color scheme: black is observations, red ERA-Interim, blue Noah cycle 1, and green Noah cycle 2.

Figure 7 .
Figure 7. Sensible heat flux (W m −2 ) on the inner domain for the cycling experiments with Noah LSM.Plots are shown at 14:00 UTC on the same days as in Fig. 2.

Figure 8 .
Figure 8. Analyzed (left) and modeled precipitation total (mm) for the month of June in the first cycle.

Figure 9 .
Figure 9. Potential temperature profiles observed by four frequent soundings at BLLAST site 2 on 25 June 2011 (black) and simulations as shown in the legend.Model results are from the nearest grid point.