Stepwise Assessment of Different Saltation Theories in Comparison with Field Observation Data

: Wind ‐ blown dust models use input data, including soil conditions and meteorology, to interpret the multi ‐ step wind erosion process and predict the quantity of dust emission. Therefore, the accuracy of the wind ‐ blown dust models is dependent on the accuracy of each input condition and the robustness of the model schemes for each elemental step of wind erosion. A thorough evaluation of a wind ‐ blown model thus requires validation of the input conditions and the elemental model schemes. However, most model evaluations and intercomparisons have focused on the final output of the models, i.e., the vertical dust emission. Recently, a delicate set of measurement data for saltation flux and friction velocity was reported from the Japan ‐ Australia Dust Experiment (JADE) Project, which enabled the step ‐ by ‐ step evaluation of wind ‐ blown dust models up to the saltation step. When all the input parameters were provided from the observations, both the two widely used saltation schemes showed very good agreement with measurements, with the correlation coefficient and the agreement of index both being larger than 0.9, which demonstrated the strong robustness of the physical schemes for saltation. However, using the meteorology model to estimate the input conditions such as weather and soil conditions, considerably degraded the models’ performance. The critical reason for the model failure was determined to be the inaccuracy in the estimation of the threshold friction velocity (representing soil condition), followed by inaccurate estimation of surface wind speed. It was not possible to determine which of the two saltation schemes was superior, based on the present study results. Such differentiation will require further evaluation studies using more measurements of saltation flux and vertical dust emissions.


Introduction
Asian dust, which refers to the wind-blown soil dust particles emitted from arid regions in North China and Mongolia and transported eastward on the prevailing westerlies [1], influences the air quality of North-East Asian countries [2][3][4]. In particular, the frequency of the outbreak of Asian dust storms in winter significantly increased in the 2000s. Although the frequency and intensity of the winter season Asian dust storms remain lower than those of spring [5], they are often accompanied by other toxic air pollutants emitted from residential heating and hence can be a considerable threat to public health [6]. In order to predict the dust storms and prepare countermeasures, a well-designed numerical wind-blown dust model, which simulates the wind erosion phenomenon to predict the wind-blown dust emission in arid regions, is required [7].
In wind-blown dust models, dust emission is calculated through a series of sequential steps to interpret the wind erosion. Therefore, the accuracy of a wind-blown dust model depends on the robustness of the schemes used for the elemental steps [8]. The first step of wind erosion The next section provides the information of observation data and a detailed explanation of the elemental steps composing wind-blown dust models. In Section 3, different saltation flux schemes are evaluated by comparing their predictions with observations. The conclusions are presented in Section 4.

Observational Data
The observational data provided by the JADE project  were used for the model evaluation in this study. Observations were conducted in Australia (33° 50′ 42.4″ S, 142° 44′ 9.0″ E) from 23 February to 14 March 2006. The 1-minute raw data for the friction velocity and saltation flux are available to the public, as shown in Error! Reference source not found.. The information on the parameters such as soil properties and the roughness length required for the model evaluation were obtained from the literature: Shao et al. [43], Ishizuka et al. [44,45], and Liu et al. [40]. More detailed information is provided in Section 3.
The model-predicted friction velocities and saltation fluxes were compared with the observation in this study to evaluate the wind-blown dust models. For this purpose, 4 distinct periods when sufficient saltation fluxes for model validation were observed, during the JADE project-Julian days 56~57, 61~62, 65, and 71-were designated as Dust Events 1, 2, 3, and 4, respectively (Error! Reference source not found.).

Domain and Configuration of Meteorological Model
Meteorological fields required for the dust emission simulation were obtained by running the meteorology model WRF v3.8.1. The horizontal resolution of the meteorological model was set to 5 km × 5 km because, according to Liu et al. [40], the size of the observation site was about 4 km from north to south and about 1 km from east to west. The NCEP FNL (final analysis) data were used. The simulation was conducted from 15 February to 16 March 2006, thus that the whole observation period was included. Error! Reference source not found. shows the summarized information of model attributes set for the WRF running. Error! Reference source not found. shows the meteorological model domain.
where (= 0.4) is the von Karman constant, the wind speed at height from the surface in the planetary boundary layer at which a logarithmic wind profile can be applied, and the roughness length due to the nonerodible elements present on the surface. The value of varies depending on the vegetation fraction and the distribution of obstacles [46]. Therefore, reliable information on surface conditions is essential for accurate estimation of the friction velocity. In this study, the effects of land-use and vegetation fraction on the roughness length and the resulting friction velocity were investigated by comparing the model-predicted friction velocity with observed values.

Roughness Length
The surface roughness length, which was an important parameter to determine the friction velocity, was widely used to determine the exchange of heat, momentum, and water between the surface and atmosphere and the characteristics of the planetary boundary layer [12]. Marticorena et al. [47] and Hebrard et al. [48] assumed that the roughness length ( ) was proportional to the physical height of obstacles (ℎ) when the total roughness density ( ) was larger than the threshold value of 0.045. Shao and Yang [13] argued that 0.045 was too small for North-East Asia with various land types such as dense tall vegetation and sparse short vegetation, and suggested a new value of 0.2.
Recently, Foroutan et al. [49] suggested Equation Error! Reference source not found. for the roughness length based on the data reported in previous studies of Shao and Yang [13] and King et al. [50]: where the roughness length decreases with increasing when is larger than the threshold value 0.2 suggested by Shao and Yang [13]. This modification was made to take into account the observation that the disturbance of wind by highly dense roughness elements was weaker than that by roughness elements with an intermediate level of density [13,49].
Equation Error! Reference source not found. was used for the determination of the friction velocity in this study. In order to examine whether the model accounts for the surface properties at the observation site properly, the land-use and vegetation fraction used in the model were compared with field observations, and their effects on the roughness length were investigated.

Threshold Friction Velocity of Dry Bare Soil
The threshold friction velocity ( * ) represents the minimum wind speed for the mobilization of soil grains. Its value for dry bare soil can be calculated by considering the balance among the wind drag, gravity, and interparticle cohesion [15]. Then the effects of soil moisture and roughness elements (drag partitioning) that increase the value of * are taken into account.
Most wind-blown dust models adopt Equation Error! Reference source not found. [15] or Equation Error! Reference source not found. [ [49] and the Sunchon National University (SCNU) dust model [52], whereas Equation Error! Reference source not found. is used in Kok's dust model [54,55] and the NMMB/BSC-Dust model [56]. In this study, both Equations Error! Reference source not found. and Error! Reference source not found. combining the effects of drag partitioning and soil moisture are evaluated by comparing with the measured threshold friction velocity during the JADE project [40].

Drag Partitioning Effect
The drag partitioning effect factor ( ) accounts for the reduction in the influence of wind on soil grains due to the physical obstacles existing on surface [41]. When drag partitioning takes place, a stronger wind is required to mobilize soil grains, which leads to a larger value of the threshold friction velocity. The WRF-DuMo and the CMAQ dust model use the drag partitioning effect factor suggested recently by Darmenova et al. (2009), expressed by Equation Error! Reference source not found., which was a modification of the original suggestion of Raupach et al. (1993) by separating the effects of vegetation and solid (nonvegetation) elements: where the subscripts and represent vegetation and nonvegetation, respectively, and A is the vegetation fraction. The "VEG" value simulated by WRF was used as A .
ln 1 is the roughness density of the vegetated surface with vegetation, the coefficient indicating the direction and distribution of roughness elements, and the roughness density of the nonvegetated surface. According to the measurement of Marticorena et al. (2006), varies in the range of 0.01~ 0.16. Xi and Sokolik [41] suggested different values corresponding to the land-use condition. The values of the constant parameters are σ 1.45, 0.16, 202, σ 1.0, 0.5, 90., and 0.35 [11]. According to Raupach (1992), the value of lies between 0 and 0.5 for a flat erodible surface and is close to 1 for a geomorphologically stable surface. The observation-based values of , , and were determined by Wyatt and Nickling [57] using multiple nonlinear regression of observations. Equation Error! Reference source not found. was used in this study to calculate the drag partitioning effect factor.

Soil Moisture Effect
The soil moisture effect factor ( ) was a key point multiplied to the threshold friction velocity for dry soil to take into account the reduced soil erodibility due to soil moisture [41,58]. Most windblown dust models use the formula suggested by Fecan et al. (1999), expressed by where is the soil moisture (kg/kg), the threshold soil moisture calculated as w′ 0.14 0.17 , and the mass fraction of clay. Equation Error! Reference source not found. was used in this study to calculate the soil moisture effect factor.
The soil moisture used in Equation Error! Reference source not found. must be the value for the top soil (0~1 cm deep) because the dust emission takes place solely from this top soil. However, most current land-surface models adopted in the meteorology models are using a top layer depth of 10 cm. Usually, the top 10 cm average soil moisture is larger than that of the top 1 cm, and hence its use in a wind-blown dust model may lead to considerable overestimation of the threshold friction velocity and underestimation of dust emission. Three different solutions have been suggested to circumvent this problem: 1. Use the Noah's land and surface scheme [59][60][61][62] with 10 cm top soil layer to obtain the soil moisture content and multiply it by a modifying factor smaller than 1 [11]; 2. Use the Pleim-Xiu scheme [63], whose top soil layer is 1 cm thick [49]; or 3. Modify the soil layer-related pretreatment part of the WRF code to construct a 1 cm top soil layer in Noah's scheme [58].
In this study, the 2 nd and 3 rd methods were used to determine the soil moisture effect factor. Two different versions of the Noah's scheme, the Noah land surface model (LSM) [59][60][61][62] and Noah-MP [64], were used for the 3 rd method. The Noah LSM simulates the changes in skin temperature, soil moisture, and snow cover state of the surface based on water and energy balance equations for 4 soil layers. This scheme has been reported to suffer a few problems such as impermeable frozen soil with snow cover and overestimation of runoff [65]. Noah-MP (Noah LSM with Multiple Physics options) is an advanced scheme developed to resolve the problems of the Noah LSM and improve various physics processes such as snow physics and soil hydrology [66]. The Noah-MP scheme also provides multiparameterization options to extend the choice for the vegetation fraction, ground surface albedo, and radiative transfer, depending on local characteristics [67]. However, this advantage may also lead to a disadvantage of having to check the combination of all physical processes to find the optimal combination [68]. The soil moisture effect factors obtained in this way were evaluated by comparing with satellite observation data.

Saltation Flux
Saltation flux ( ) refers to the amount of soil grains moving horizontally by frog leaping when the friction velocity exceeds its threshold value. The value of depends on the friction velocity and on the variables that affect the threshold friction velocity: soil texture, soil moisture content, vegetation fraction, soil grain size distribution, etc. [15,69,70].
Most saltation schemes adopt the mathematical forms that can be represented by Equation Error! Reference source not found.: * According to this formula, the saltation flux increases proportionally to the th power of the friction velocity when the friction velocity is larger than its threshold value. The parameters and can be tuned by comparing the model prediction with the measured saltation flux. In this study, we tried to determine the values of and simultaneously by curve-fitting observed saltation data using Equation Error! Reference source not found..
The most widely used formula for the saltation flux is that suggested by Kawamura [71] and supported by White (1979) Equation Error! Reference source not found. is adopted in the CMAQ dust model [49] and in the SCNU model [52]. When the friction velocity is large enough, Equation Error! Reference source not found. can be approximated by Equation Error! Reference source not found. with 3.
On the other hand, Kok et al. [72] suggested the following formula using theoretical reasoning and validated it through elaborate field experiments [73]: * * * * Equation Error! Reference source not found. corresponds to 2 when compared to Equation Error! Reference source not found.. Both Equations Error! Reference source not found. and Error! Reference source not found. are evaluated in this study by comparing their predictions with the observed saltation flux.

Effect of the Roughness Length on the Friction Velocity
As was explained in the previous section, the friction velocity depends on the roughness length that varies with land-use condition and vegetation fraction. This necessitates checking the land-use type and vegetation fraction provided to the wind-blown dust model.
The gridded land-use type and vegetation fraction are determined in the meteorology model using the 20-category IGBP-Modified MODIS data. Shao et al. [43] reported that the surface of the observation site consisted of 84% bare surface, 1% rock, 1% tree litter, 12% litter, and 2% crust, with little vegetation. The roughness length was reportedly 0.48 mm, which was within the range of the roughness length values reported for East Asian dust source regions 0.00382 ~ 0.642 mm [42]. This indicated that the land-use condition and the roughness length of the observation site were similar to those of the Gobi desert.
On the other hand, the land-use type and vegetation fraction determined for the observation site in the model were shrubland and 25.28%, respectively, showing considerable difference from the report of Shao et al. [43]. To assess the effects of the model input for the land-use type and vegetation fraction on the friction velocity, the roughness lengths and friction velocities were calculated under three different conditions and compared with one another. In CASE 1, the model input land-use type and vegetation fraction were used. The land-use type was changed from shrubland to barren land in CASE 2, whereas in CASE 3 the vegetation fraction was also changed to 2% thus that the roughness length became the same as the 0.48 mm reported by Liu et al. [40]. As shown in Table 2 and Figure 3, CASE1 and CASE2 based on the vegetation fraction and/or land-use type that did not correspond with observed surface conditions resulted in the overestimation of the roughness length and the friction velocity. On the other hand, the index of agreement (IOA) was the highest in CASE 3, in which the roughness length was adjusted to the reported value by using more realistic land-use type and vegetation fraction value. Therefore, the accuracy of the vegetation fraction as well as the landuse type was considered to contribute significantly to the prediction accuracy of the wind-blown dust model.

Drag Partitioning Effect
As is indicated in Equation Error! Reference source not found., the drag partitioning effect factor was also dependent on the vegetation fraction and land-use condition. Therefore, the drag partitioning effect factor was calculated under three conditions shown in Error! Reference source not found.: 3.438, 2.074, and 1.110 for CASE 1, CASE 2, and CASE 3, respectively. The validity of these values was discussed in Section 3.4.

Effect of Land Surface Scheme on Soil Moisture
The three land surface schemes introduced in Section 2.2.6 were evaluated by the soil moisture values for 2016 predicted by them with measured values. The SMAP satellite observation data provide the 5-cm top layer soil moisture contents with a temporal resolution of 3 h and a spatial resolution of 9 km × 9 km. The characteristics of the observation dataset and the land surface schemes evaluated are summarized in Error! Reference source not found.. Fifteen points (5 each from the Gobi desert, grassland, and cropland) were chosen arbitrarily from the main wind-blown dust source regions of North-East Asia, as shown in Error! Reference source not found.. The meteorology model WRF was run with the top soil layer thickness of 5 cm for consistency with observation. Table 3. Information of the two soil moisture observation datasets used for the evaluation of the land surface schemes.

Contents Analysis Site East Asia Australia
Observation data      Therefore, the Noah-MP scheme was selected for further validation through the comparison with the soil moisture measurements made in the JADE project. Liu et al. [40] and Shao et al. [43] reported that the soil moisture contents at the soil depth of 0.

Estimation of the Threshold Friction Velocity
The threshold friction velocity of the bare soil was determined, using Equations Error! Reference source not found. and Error! Reference source not found. independently, to be 0.209 m s -1 and 0.227 m s -1 , respectively. These two values differed by less than 10%, which was in agreement with a previous study that reported a small difference in the threshold friction velocities calculated using these two formulas [52]. It was not possible to determine the superiority of one of the two schemes over the other based solely on the predicted threshold friction velocity. Therefore, Equation Error! Reference source not found., which was simpler and delivers the physical meaning of the formula more intuitively [52], was chosen for further analysis.
The drag partitioning effect factors calculated in Section 3.2 were multiplied to the threshold friction velocity of bare soil (0.209 m s -1 ) to determine the threshold friction velocity of dry surface with nonerodible elements: 0.716 m s -1 for CASE 1, 0.434 m s -1 for CASE 2, and 0.231 m s -1 for CASE 3. Liu et al. [40] and Shao et al. [43] reported that the threshold friction velocities determined for Dust Events 2 and 4 were 0.2 m s -1 and 0.28 m s -1 , respectively. Therefore, even without taking the soil moisture effect into account, the threshold friction velocity of dry surface determined for CASE 1 or 2 was much larger than the observed values. This comparison further demonstrated the need to use the land-use and vegetation fraction (CASE 3), which can match the measured roughness length.
The soil moisture effect factors determined in Section 3.3 were multiplied to the threshold friction velocity of dry surface (0.231 m/s) to determine the threshold friction velocity of wet surface: 0.250 m s -1 , 0.247 m s -1 , and 0.245 m s -1 for 10 cm top soil, 1 cm top soil, and 1 cm top soil with compensation factor, respectively. These three values differed by much less than those caused by the difference in the land-use and vegetation fraction. This result indicates that the Australian desert contains such a low soil moisture that the soil moisture effect was very low and the effect of the top soil layer depth may not be important. However, the dust source regions of North-East Asia are known to have relatively high soil moisture contents. Therefore, the effect of the top soil layer depth may be much more important in North-East Asia. The detailed results for the threshold friction velocity of wet soil determined under various conditions are summarized in Error! Reference source not found..  The soil condition during the JADE project was reported to be largely unchanged, implying that the threshold friction velocity must be represented by a single value. Therefore, in Error! Reference source not found.e we curve-fitted all the measured saltation flux values together. The most adequate fitting was obtained when the threshold friction velocity value was set to 0.28 m s -1 . In this case the optimum exponent was determined to be 4.49. Ishizuka et al. [45] reported the optimum exponent of about 4 for the same project data set using a formula similar to Equation Error! Reference source not found., which was very close to the value obtained in this study.

Comparison of Two Saltation Schemes with the Measured Friction Velocity
Wind-blown dust models show different results when calculating the saltation flux, depending on the exponent of the saltation flux scheme. In this subsection, the two representative saltation schemes introduced in Section 2.2.7, Kawamura and White (Equation Error! Reference source not found.) and Kok (Equation Error! Reference source not found.) expressed as u* 3 and u* 2 , respectively, to focus on the exponent value, were evaluated by comparing their predictions for the saltation flux with the measurements. In Error! Reference source not found., the saltation flux values measured during all the dust event periods were curve-fitted as a function of the measured friction velocity using Equations Error! Reference source not found. and Error! Reference source not found. with three different threshold friction velocity values of 0.2 m/s, 0.24 m s -1 , and 0.28 m s -1 . Again, the best fitting was obtained with the threshold friction velocity value of 0.28 m s -1 . The value of the coefficient k determined for this case was 0.828 for the Kawamura and White scheme and 1.910 for the Kok scheme. Both of these values are smaller than the coefficients suggested originally by the scheme developers: 2.6 [36] and 5 [72,73], respectively. CMAQ used the coefficient of 1.0, rather than the original value of 2.6, when it adopted the Kawamura and White scheme in its wind-blown dust emission module. The result of this study seems to support the coefficient taken by CMAQ. The two schemes showed very similar levels of correlation coefficients and IOAs.

Comparison of Two Saltation Schemes with the Model-Predicted Friction Velocity
In this subsection, a similar analysis was performed to that reported in Section 3.5.2 but the model-predicted threshold friction velocities, instead of the measured ones, were used to calculate the saltation flux using the two saltation schemes.
Error! Reference source not found.a,b show the results obtained using the model-predicted average threshold friction velocity for the whole project period of 0.245, whereas Error! Reference source not found.c,d show the results obtained using a measured threshold friction velocity of 0.28 that provided the best fitting in the previous subsection. Again, the threshold friction velocity of 0.28 led to higher IOA values for both schemes and the two saltation schemes did not show a large difference in their performances. However, both the correlation coefficients and IOAs obtained for both schemes were much smaller than those obtained with measured friction velocities when Error! Reference source not found.(8e,8f) andError! Reference source not found.(9c,9d) were compared. This result indicates that inaccurate estimation of wind field by the meteorology model degrades the accuracy of the saltation schemes and that accurate estimation of the input parameters is more important than the functional form (i.e., the value of the exponent n) of the saltation schemes for accurate prediction of the saltation flux.

Conclusions
The saltation flux measurements obtained from a desert region of Australia, whose surface characteristics were similar to those of the Gobi desert, were used to evaluate two widely used physics-based saltation schemes. The accurate prediction of the threshold friction velocity was shown to be the most important element to enhance the accuracy of the saltation schemes, which requires the realistic representation of land use, vegetation fraction, and soil conditions. The Australian desert where the JADE project was conducted has very low soil moisture content, precluding quantitative assessment of the effects of the soil moisture, which may be much more important in North-East Asian dust source regions. Both schemes tested in this study simulated the saltation process reasonably well when using measured input conditions. Further evaluation of the wind-blown dust models will require comparison with more measurements of the saltation flux and vertical dust emissions in other desert regions, including North-East Asia.
Supplementary Materials: The following are available online at www.mdpi.com/xxx/s1. Figure S1: Scatterplots comparing the soil moisture contents predicted by different land surface schemes with observation at different sites. The information of the land surface schemes used and the observation sites is provided in each sub-figure.