Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

The distribution of Van Genuchten model parameters on soil-water characteristic curves in Chinese Loess Plateau and new predicting method on unsaturated permeability coefficient of loess

  • Shiyue Fang,

    Roles Conceptualization, Formal analysis, Funding acquisition, Methodology, Project administration, Resources, Supervision, Writing – review & editing

    Affiliation College of Geology and Environment, Xi’an University of Science and Technology, Xi’an, Shaanxi, China

  • Pengfei Shen ,

    Roles Conceptualization, Data curation, Formal analysis, Methodology, Project administration, Resources, Software, Supervision, Visualization, Writing – original draft

    408517463@qq.com

    Affiliation College of Geology and Environment, Xi’an University of Science and Technology, Xi’an, Shaanxi, China

  • Xinhai Qi,

    Roles Investigation, Resources

    Affiliation Gansu Shaanxi Branch of West to East Gas Transmission Company, PipeChina, Xi’an, Shaanxi, China

  • Fan Zhao,

    Roles Investigation, Validation

    Affiliation College of Geology and Environment, Xi’an University of Science and Technology, Xi’an, Shaanxi, China

  • Yue Gu,

    Roles Validation

    Affiliation College of Geology and Environment, Xi’an University of Science and Technology, Xi’an, Shaanxi, China

  • Jiaxin Huang,

    Roles Validation

    Affiliation College of Geology and Environment, Xi’an University of Science and Technology, Xi’an, Shaanxi, China

  • Yan Li

    Roles Data curation

    Affiliation North Shaanxi Mining Hong Liu Lin Company Limited, Yulin, Shaanxi, China

Abstract

The unsaturated permeability coefficients are often used to solve geotechnical problems associated with unsaturated soils. But it is very difficult to measure. However, the unsaturated permeability coefficients can be predicted by the Soil-water Characteristic Curves (SWCCs). The Van Genuchten Model (VG model) is very rife as it’s smooth and good fitting, thus, it has the most research data. Therefore, the research data on VG model parameters (α, n, θs and θr) of Malan loess in Chinese Loess Plateau are collected in the past two decades to obtain the spatial distribution characteristics of parameters. The trend surface analysis method is employed to clarify the regional scale distribution and the variation regular pattern on ArcGIS. Then the linear regression method is utilized to fit the relationship between suction and water content in three different regions of Chinese Loess Plateau, which is divided according to the properties and particle gradation. By using this relationship and the trend surface analysis contour map, the unsaturated permeability coefficient of the sample can be predicted after measuring the saturated permeability coefficient. The example verification shows that the difference between the prediction results and the experimental results is very small when the sample has the lower saturation, and the deviation is slightly larger if it has the higher saturation, but they are all within the acceptable range. This method not only saves the test cost, but also considers the physical properties of the loess in the three different regions of the Loess Plateau. With the improvement of data and the gradual improvement of sampling density, the prediction accuracy will gradually improve. It can provide convenience for solving the engineering problems of loess and water and other engineering applications.

Introduction

The permeability coefficient, also known as the hydraulic conductivity, represents the seepage velocity of water in different media under the unit hydraulic gradient. The unsaturated permeability coefficient is related to the water content, and the value of the unsaturated permeability coefficient changes continuously with the change of the matrix suction, and the variation range can span several orders of magnitude, so it is difficult to accurately measure the unsaturated permeability coefficient. Unsaturated permeability coefficients are often used to solve geotechnical problems associated with unsaturated soils. In view of the difficulty of its measurement, many scholars have done extensive research on the prediction of unsaturated permeability coefficient. Mualem [1] made a systematic analysis of the calculation models of unsaturated soil permeability coefficient, and divided these calculation models into the following three categories: (1) empirical models; (2) statistical models; (3) macroscopic models.

Empirical models can be roughly divided into two categories. Expressing the permeability coefficient of unsaturated soils with moisture content, such as the research of Gardner [2]、Davidson et al. [3]and Campbell [4]. Or expressing the permeability coefficient of unsaturated soils with matrix suction, Richards [5], Christensen [6], Wind [7], Gardner [2], Rijtema [8], Philip [9]and Leong & Rahardjo [10]have carried out more in-depth research in this category.

The statistical models is derived from the Hagen-Poiseuille equation [11] and the Kelvin capillary model. This type of model established by the corresponding hydraulic parameters which is based on the probability of connected pores in soil under seepage [1, 12, 13]. Childs&Collis-Geroge [14] established the statistical model of unsaturated permeability coefficient in 1950, after that, Marshall [15]&Kunze [16] made improvements to this statistical model. In the process of calculating the unsaturated permeability coefficient, the soil-water characteristic curve and the volumetric moisture content test data are firstly divided into equal parts, then the water moisture is calculated based on the median value of each division, and finally apply the permeability coefficient equation to calculate the permeability coefficient.

The macroscopic models relates the permeability coefficient of unsaturated soils to the effective saturation. Such as, Gardner [2] Brooks & Corey [17] Mualem [18] have done related research, of which the Mualem model is the most representative model, which is describes the relationship between saturated permeability coefficient and soil-water characteristic curves. Because the saturated permeability coefficient is easy to measured, so the soil-water characteristic curve plays an important role in predicting unsaturated permeability coefficient. Nowadays, it has been widely used to the estimation of unsaturated soil properties.

The Soil–water Characteristic Curves (i.e. SWCCs) describe the relationship between water content and matrix suction. The water content can be volumetric moisture content (θ), gravimetric moisture content (ωW) or degree of saturation (Sr). Matrix suction can also be called capillary pressure. It can be represented by (uauw), where ua represents pore air pressure and uw represents pore water pressure [13, 19, 20]. A typical SWCC is shown in Fig 1 [21]. The SWCCs (both drying and wetting) are sigmoidal in shape, the two transition points on the drying curve are the air-entry value and residual suction value [22]. The air-entry value(G1) is defined as the matric suction at which air enters into the largest pores by drying. The residual water content (θr) is the water content at which water phase is discontinuous. The suction value corresponding to θr is called the residual suction value(G2). The air-entry value and residual suction value divide the drying SWCC into three identifiable stages, namely, boundary effect zone, transition zone and residual zone [19, 23].

Generally, there are two ways to get the SWCC for a soil sample [24], direct measure methods and indirect predict methods. However, direct measure methods have some defects. Some test instrumentations are restricted by their design parameters, so, some could not measure all suction range [25], some have shortcomings of the time-consuming, special demands on the testing environment, and some need expensive testing cost [26]. So far, there are three different approaches to predict the SWCCs .1) database mining for the SWCC from similar soils, 2) prediction from the grain-size distribution(i.e. GSD) curve, 3) correlations relating parameters in the SWCC equations to soil properties [22]. Numerous studies have been performed on the last two approaches. ARYA and PARIS [27] were the first to propose a physico-empirical model for predicting the SWCC from the GSD. The GSD is translated into the pore-size distribution (i.e. PSD), which is then related to the SWCC through the capillary theory. The pore radius prediction was based on the assumption of spherical particles and cylindrical pores. The ARYA and PARIS [27] model was later modified and improved by considering the random packing of spherical particles and the influence of soil structure [28]. Recently, the PSD was directly measured and used to predict the SWCC and changes in the PSD during drying and wetting were also considered. SIMMS and YANFUL [29] developed a model for predicting the drying SWCC based on the evolution of measured PSDs for a compacted clayey soil. HU et al. [30] proposed a hysteretic SWCC model based on changes in the PSD taking account of the influence of volumetric deformation on the SWCC. There are other physico-empirical models that do not have to translate the GSD into the PSD. And The methods which relate the fitting parameters in mathematical equations to soil properties. The widely-known SWCC models, such as the Grander [2] model, BROOKS and COREY [17] model, Van Genuchten [13] model and FREDLUND & XING [31] model. Among them, the Van Genuchten [13] model is widely used, it can better fit the soil-water characteristic curves of various soils. which is: (1) Where θω is the volumetric moisture content, cm3/cm3; θs is the saturated water content, cm3/cm3; θr is the residual water content, cm3/cm3; ψ is the matrix suction, kPa; α is the parameter related to the air entry value, kPa; n is the parameter controls the slope of the soil-water characteristic curve; .

The mathematical form of the Van Genuchten Model (i.e. VG model), which accounts for an inflection point, thus it allows greater flexibility than the other model over a wider range of suction and captures the sigmoidal shape of typical curves better. Smooth transitions at the air-entry value and approaching the residual suction value condition are more effectively captured. Parameter m (i.e. ) reduces the flexibility of the VG model but significantly simplifies it. Thus resulting in greater stability during parameter optimization and permitting closed-form solution of the hydraulic conductivity function [21].

Loess soils are widely distributed in arid and semi-arid regions, covering about 10% land area of the world. Many countries in Asia, Europe, North America, South America and Africa have loess soil deposits. In China, the renowned Loess Plateau extends over 6.4×105 km2, accounting for over 7% land area of the country [22, 3236]. The Chinese Loess Plateau originates from aeolian dust deposit since 2.6 million years ago, which spans over the Pleistocene period (i.e. 2.4–0.01 Ma) and Holocene period (i.e. recent 0.01 Ma). Since the structure of Holocene loess and Malan loess is relatively loose, the engineering geological problems of loess are common in these two layers [34]. However, Holocene loess is usually cultivated layer, and engineering construction is mostly in Malan loess, such as retaining wall foundation, oil and gas pipeline trench foundation, etc. Therefore, the study of loess engineering properties is generally concentrated in Malan loess.

Malan loess is typically in a state of unsaturated condition. The unsaturated permeability coefficient is usually used to address geotechnical problems associated with Malan loess [37]. In order to solve engineering geotechnical problems in Malan loess better, research on the regional distribution of SWCCs, especially, VG model and the unsaturated permeability coefficient become more and more important.

In view of this, the VG model four parameters (α, n, θs and θr) of Malan loess SWCCs tested by experts and scholars will be collected in Loess Plateau for recent years. The regional distribution of the VG model parameters will be elucidated by employing the trend surface analysis and the linear regression methods. And through the connection between the parameters of the VG model and the unsaturated permeability coefficient, a new method for calculating the unsaturated permeability coefficient of loess is proposed. In this way, the engineering problem related to water in loess will be solved easily and provide convenience for other engineering applications.

Materials and methods

Materials

SWCCs of Malan loess is the key to solve engineering problems related to loess and water. Hence, the parameters of VG model on SWCCs of Malan loess in Chinese Loess Plateau are gathered and summarized as shown in S1 Table & Table 1.

The sampling ranges are all located in the Malan loess stratum at a distance of 2-5m from the upper strata boundary, according to the soil properties and particle gradation, Fan et al. [38] divided the loess of the Loess Plateau into three zones, which are shown in Fig 2. By sorting out the sampling results, get the difference in particle size composition of each zone, which are shown in Table 2. The sandy loess (Zone I), has the highest sand content, clay content is 6%~11%, powder content is 55%~65% and sand content is 25%~33%. The typical loess (Zone II), contains a large amount of powder and sand and a small amount of clay, clay content is 7%~12%, powder content is 60%~70% and sand content is 20%~27%. The clay loess (Zone III), has the highest clay content and a large amount of powder, clay content is 10%~20%, powder content is 70%~75%, and sand content is 5%~14%. In this way, the collected data on VG model parameters of each zone can be statistically analyzed, which are shown in Table 3.

thumbnail
Fig 2. The zoning of the Chinese Loess Plateau(Only the boundary of the Loess Plateau comes from Huang [32], illustrative purposes only).

https://doi.org/10.1371/journal.pone.0278307.g002

thumbnail
Table 2. Particle gradation components of selected the sampling results.

https://doi.org/10.1371/journal.pone.0278307.t002

According to the data in the Table 3, the parameters show obvious regularity. The average value of Parameter θs and θr in Zone II are both larger than the other two zones. The average value of Parameter α in Zone I is less than that in Zone III and Zone II. The average value of Parameter n is biggest in Zone I, second in Zone II and smallest in Zone III.

Methods

Trend surface analysis method.

Trend surface analysis is a multivariate statistical analysis method. It takes the geological variable Z as the dependent variable and the two-dimensional plane coordinates X and Y as the independent variables. A multi-form (high-order) regression equation can be constructed to form a plane or more complex spatial surface in geometry to reflect the geological variable in space. Therefore, it separates the regional variation characteristics (trends) of environmental variables in space from the local variation characteristics (abnormalities). It can also remove the regional variation component and highlight the local variation component, to study the abnormal points and analyze their causes [71, 72].

Radial basis function method.

Radial basis function (RBF) is a real-valued function whose value depends only on the distance from the origin [73], that is (2)

Or it can be the distance to any point c, where c is called the center point, that is (3)

Any function Φ(x) that satisfies Eq (3) is called a radial basis function, and the distance is calculated by Euclidean distance [74, 75].

The radial basis function is an effective method to deal with the multivariate function approximation problem. The main work of its studies is the expansion of space surface and its properties, and applies this space surface to handle the description problem of the data. It has the advantages of high efficiency, simple operation, easy programming, convenient storage in the computer, isotropy, independent of grids, and high solution accuracy. It has been extensively researched and applied in many fields, such as, function approximation, neural networks, multi-scale analysis, mining, geophysics, surveying and mapping, remote sensing and signal processing and machine learning [76].

Universal Kriging.

The variability of regional variables is represented by three parts: the deterministic part, the relevant part and the random part. The ordinary Kriging method requires variables to satisfy the second-order stationary assumption or inherent assumption, that is, assuming that the deterministic part is constant in space. For non-stationary variables, the deterministic part is not constant in space. It must be assumed that the deterministic part is distributed with space, also known as drift or tendency. The best linear estimation process corresponding to this method is called universal Kriging method [7780].

Because the mean value is no longer a constant in space but a spatial variable, considering a regional variable z(x) defined in the study area A at position x, it can be represented by a deterministic drift m(x) and a residual part r(x): (4)

Matrix representation: (5)

Through the definition of drift, the expectation of z(x) at x is m(x), that is (6)

If the drift function is known, the drift can be minus in the original data. Then, if the remained residual part r(x) satisfies the inherent assumption, the residual part can be estimated by Kriging method. Hence, the estimated value can be obtained by plus the drift to the estimated residual part at relative position [79].

Assuming that m(x) can be expressed by a linear combination of functions Pk (x)(k = 1,2,…,K): (7)

Matrix representation: (8) Where ak is an unknown coefficient; pk (x) is a known function of x, which is often expressed by a polynomial of xk−1. (9) Where P is (n × k) polynomial matrix; n is the number of data point; p is the number of drifts, which can be determined by multiple stepwise regression; A is the polynomial coefficient of drift; R is the theoretical residue.

Assuming that the residual part r(x) satisfies the second-order equilibrium condition, σ(h) is used to represent its covariance coefficient, the covariance of r(x′) and r(x) are only a function of distance xx′: (10)

Matrix representation: (11)

Assume that σ(h) is known. The main problem of universal kriging is to determine the best coefficient ak of the drift part [79].

Operation on ArcGIS.

Importing the sorted data into ArcGIS for analysis, the radial basis function method is used for the measured value contour graph, because it has high efficiency when dealing with multiple problems, can handle large-scale scattered data, and has good approximation capabilities. The universal Kriging method is employed for the prediction value contour graph and the prediction standard error value graph, because it is a geo-statistical method that takes into account the unbiased linear estimator with drift. It can allow the existence of measurement errors.

By using ArcGIS for analysis, we can establish a more realistic way of predicting SWCCs and unsaturated permeability coefficients that can be updated at any time in the Chinese Loess Plateau. Because each piece of data is obtained from actual experiments, it is real and effective, and it is closer to the actual soil properties of the Chinese Loess Plateau. As the sampling density increases, its accuracy will gradually become more accurate.

Through the method below, we can analyze each VG model parameter of the Malan loess SWCCs in Chinese Loess Plateau.

  1. Step I, input the map. Import study area boundaries(the boundary map of the Loess Plateau) and administrative division boundaries into ArcGIS as vector files. Superimposed these two vector files to obtain the administrative division map of the Loess Plateau in Arcmap.(Only the boundary of the Loess Plateau in this article comes from Huang [32] as a basemap for vectorization, illustrative purposes only).
  2. Step II, import the Data. Add the geographic location points on the administrative division map of the Loess Plateau, and add the data of α, n, θs and θr to its attribute table.
  3. Step III, perform the spatial deterministic interpolation. Utilize Radial Basis Function Interpolation Method function in Interpolation Analysis under the Geostatistical Analyst Tool in Arctoolbox. After inputting the data of α, n, θs and θr respectively, employ the boundary layer of the Loess Plateau in environment as the processing range, so as to obtain the measured value contour graph of each parameter of the VG model.
  4. Step IV, perform the geostatistical interpolation. Utilize the Geostatistical Wizard command under Geostatistical Analyst in the toolbar. Select the Kriging/CoKriging interpolation method in the pop-up dialog box, and enter the α, n, θs and θr data respectively. Then, choose universal kriging method in the kriging method type list box, employ prediction and prediction standard error respectively in the output type list box, and then config the required factor to generate the prediction value contour graph and the prediction standard error value graph for each parameter value of the VG model.
  5. Step V, perform the trend surface analysis. Utilize Trend analysis in Explore data under Geostatistical Analyst in the toolbar, and input α, n, θs and θr data to obtain the trend analysis graph of each parameter.

Results

VG model parameter α

From the trend analysis graph (Fig 3A) and the prediction value contour graph (Fig 3C), the trend regular pattern of the VG model parameter α in the Chinese Loess Plateau can be obtained. The measured value contour graph (Fig 3B) can illustrate the spatial distribution of the VG model parameter α. The prediction standard error value graph (Fig 3D) reflects the accuracy of the predicted value of the VG model parameter α and the abnormal situation of the data. Combining these four graphs, the regional distribution of VG model parameter α in the Chinese Loess Plateau is obtained as follows:

  1. The trend value of the parameter α (A) is higher in the southeast and lower in the northwest. It decreases from the west to the east and then increases gradually. Overall, the value in the west is higher than in the east.
  2. The measured value of parameter α (B) is basically higher in Zone II, followed by Zone III, and lower in Zone I. Among them, the values are the higher (0.135~0.2) at Zichang of Shaanxi province and Xining of Qinghai province in Zone II, it is the lowest (0.012~0.02) at Ordos in the Chinese Loess Plateau.
  3. The predicted value of parameter α (C) is higher in Zone II and Zone III, and lower in Zone I. Among them, Ansai, Yan’an and Zichang in north-central Shaanxi province, Sanmenxia in northern Henan province and Xining in Qinghai province have higher values (0.07~0.09). Southern Ningxia province, eastern Gansu province and Ordos have the lowest values (0.01~0.02).
  4. The prediction standard error value of the parameter α (D) has little deviation between the measured value and the predicted value in most areas, especially in the central and southern parts of Shaanxi province in the Chinese Loess Plateau, where the error values are generally less than 0.017. However, in the Luoyang area, the error value is larger (0.03). But most error values are small, which certify the predicted value are accurate in most areas.
thumbnail
Fig 3. Analysis result of parameter α.

Trend analysis graph(A). measured value contour graph(B). prediction value contour graph(C). prediction standard error value graph(D) (Only the boundary of the Loess Plateau comes from Huang [32], illustrative purposes only).

https://doi.org/10.1371/journal.pone.0278307.g003

VG model parameter n

By the same analysis method, we can acquire the regional distribution of VG model parameter n in the Chinese Loess Plateau, which is shown as follows:

  1. From the trend analysis graph (Fig 4A) and the prediction value contour graph (Fig 4C), we can see that the trend value of the parameter n is higher in the north and lower in the south. It decreases from the east to the west and then increases gradually. Overall, it is higher in the northwest than in the southeast.
  2. As shown in Fig 4B, the VG model parameter n is higher obviously in Zone I, followed by Zone II, and lower in Zone III. Among them, the values are higher (2.03~2.11) at Heifangtai and Dongxiang in Gansu province than other areas in Zone II. And Qishan is the lowest (1.02~1.14), which is in Zone III.
  3. From Fig 4C, we can see that the predicted value of parameter n is highest in Zone I, second in Zone II, and lowest in Zone III, and its value varies significantly with the three divided zones. The highest value is in the eastern part of Ordos (1.74~1.85), and the lowest values (1.29~1.46) are in the eastern part of Xi’an and the southern part of Sanmenxia.
  4. As described in Fig 4D, the prediction standard error value of the parameter n reflects a little deviation between the measured value and the predicted value in the majority of areas, especially in the central Shaanxi province and eastern Gansu province, where the error values are generally less than 0.09. However, in Xining of Qinghai province, the error value is the biggest (0.15~0.17). The error values are small and attest predicted value of parameter n are accurate in most areas.
thumbnail
Fig 4. Analysis result of parameter n.

Trend analysis graph(A). measured value contour graph(B). prediction value contour graph(C). prediction standard error value graph(D) (Only the boundary of the Loess Plateau comes from Huang [32], illustrative purposes only).

https://doi.org/10.1371/journal.pone.0278307.g004

VG model parameter θs

From these four graphs, the regional distribution of VG model parameter θs in the Chinese Loess Plateau can be obtained as follows:

  1. As shown in Fig 5A and 5C, the trend of the VG model parameter θs is higher in the southeast and lower in the northwest. It decreases from the west to the east and then increases gradually. Overall, it is higher in the east than in the west. And it increases slowly after gradually decreasing from the north to the south. Overall, it is higher in the south than in the north.
  2. The measured value contour graph (Fig 5B) can illustrate that the VG model parameter θs is higher in Zone II obviously, which the high value is relatively centralized. Then, Zone III is the second, and lower in Zone I. Among them, they are the higher value (0.52~0.56) in the triangle area of Tongchuan, Heyang and Luochuan in Shaanxi province, which is in Zone II. Liulin in Shanxi province and Ordos in Inner Mongolia province are followed (0.47~0.5). In the entire Chinese Loess Plateau, it is the lowest (0.23~0.31) in Yangling, Jingbian in Shaanxi province and Hequ in Shanxi province. In northern Ningxia province and western Inner Mongolia province (0.31~0.4) were lower too.
  3. Fig 5C shows that the predicted value of parameter θs is higher in Zone II, and lower in Zone I and III. Among them, Luochuan, Tongchuan, and Heyang in Central Shaanxi, Liulin in Western Shanxi, and Ordos in Inner Mongolia this band shape area have higher values. Among them, Luochuan, Tongchuan, and Heyang in central Shaanxi have the highest values in the triangle area (0.52~0.56). Other regions, such as Northwestern Shanxi, Southern part of Yangling, and Northern part of Jingbian to Western Inner Mongolia have the lowest values (0.23~0.31).
  4. The prediction standard error value graph (Fig 5D) reflects that the accuracy of the predicted value of the VG model parameter θs has a little deviation between the measured value and the predicted value in most areas, especially in most areas of Shaanxi province, western Shanxi province, and the belt shape region from Xining to Lanzhou, the error values are generally less than 0.06. However, in Yinchuan of Ningxia province and Hohhot of Inner Mongolia province, the error values are larger (0.124~0.129).
thumbnail
Fig 5. Analysis result of parameter θs.

Trend analysis graph(A). measured value contour graph(B). prediction value contour graph(C). prediction standard error value graph(D) (Only the boundary of the Loess Plateau comes from Huang [32], illustrative purposes only).

https://doi.org/10.1371/journal.pone.0278307.g005

VG model parameter θr

From Fig 6, we can get the regional distribution of VG model parameter θr in the Chinese Loess Plateau, which is shown as follows:

  1. From Fig 6A and 6C, the trend regular pattern of the VG model parameter θr in the Chinese Loess Plateau can be obtained. It is higher in the southeast and lower in the northwest. It rises gradually from west to east and then falls rapidly. Overall, it is higher in the east than in the west. In the north-south direction, it descends rapidly from north to south and then rises slowly. Overall, the trend value in the south is higher than in the north.
  2. Fig 6B can illustrate the spatial distribution of the VG model parameter θr. It is higher in Zone II obviously, followed by Zone III, and lower in Zone I. Among them, Baishui of Shaanxi province and Qingyang of Gansu province have the highest values (0.16~0.18), which is in Zone II. While the values in northern Shaanxi province and northwest Shanxi province are lower. The Shenmu in the north of Shaanxi province is the lowest (0.03~0.04).
  3. From Fig 6C, we can see that the predicted value of parameter θr is higher in Zone II and Zone III, and lower in Zone I. Among them, Jingbian, Yulin, and Shenmu in north-central Shaanxi province, Hequ and Xingxian in northwest Shanxi province have lower values (0.03~0.04). The highest values are in eastern Gansu province and southern-central Shaanxi province (0.17~0.18).
  4. Fig 6D reflects that the accuracy of the predicted value of the VG model parameter θr has little deviation in southern-central of Shaanxi province and Lanzhou and it is generally less than 0.019. In Xining of Qinghai province, Yinchuan of Ningxia province and Hohhot of Inner Mongolia province, the error values are larger (0.05). The error value in Shaanxi Province is relatively small, less than 0.03.
thumbnail
Fig 6. Analysis result of parameter θr.

Trend analysis graph(A). measured value contour graph(B). prediction value contour graph(C). prediction standard error value graph(D) (Only the boundary of the Loess Plateau comes from Huang [32], illustrative purposes only).

https://doi.org/10.1371/journal.pone.0278307.g006

Discussion

Analysis on the reasons of distribution

Combined with scholars’ research on SWCC and the prediction map, the property of SWCCs on the Loess Plateau can be analyzed as follows.

The previous work on unsaturated soil has indicated that the particle size composition of loess has an important influence on the SWCCs of loess. If the content of gross particles is high, the pores of the loess are larger and it is easier to form a drainage channel [20]. In addition, the specific surface area of loess is small, so the adsorption capacity of water will be reduced. Therefore, the air intake value of gross particle soil is relatively small, the dehumidification rate is relatively large, and the residual moisture content is relatively low [38, 81]. If the content of fine grains of loess are high and the size of soil particles are small, the pores of the loess are small and the drainage is more difficult [20]. Because of the large specific surface area of the loess, the adsorption capacity for water is strong. Therefore, the adsorption capacity of clay particles is the strongest, powder particles is second, and the sand particles is the weakest. At same time, the higher the fine particle composition is, the greater the air intake value is, the smaller the dehumidification rate is, and the greater the residual moisture content is [38, 81].

For the parameter α, studies have shown that it has a higher correlation with the air entry value. The larger the air entry value is, the smaller the parameter α is [8183]. It also shows that the higher the fine particle content of loess is, the greater the value of the parameter α is. As the content of fine particles is relatively high in the central and southern parts of the Chinese Loess Plateau, while the northern area is mainly sandy loess, hence, the α values of Zone II and Zone III are higher than those in Zone I. The reason that some areas have the abnormal value is that the air entry value is not only affected by the particle composition, but also controlled by the loess porosity [38], and its value is also susceptible to external disturbances during the measurement process.

For the parameter n, studies have shown that it is a parameter related to the moisture reduction rate. If the moisture reduction rate is high, the value of parameter n is large too. [8183]. That is, the higher the fine particle content of the loess is, the smaller the value of the parameter n should be. Because the clay loess area and the typical loess area have more fine particles, and the sand loess area has more gross particles, the value of parameter n is low in Zone II and III and high in Zone I. And, therefore, the prediction map of the parameter n shows that it’s low in the southeast and high in the northwest consistently.

For the parameter θs, its value is mainly related to the porosity ratio of loess. The higher the porosity is, the greater the saturated water content of the loess has. Because of the fine particle content of loess, the void ratio is great. It can be inferred that the tinier the loess particles is, the larger the parameter θs should be. According to the analysis of the particle size division of the Chinese Loess Plateau and the spatial trend of the parameter θs, it can be concluded that the parameter θs is higher in the southeast and lower in the northwest, and the high values are basically distributed on both banks of the Yellow River. It is presumed to be a large amount of fine particle matter brought by the alluvium of river. Therefore, the high value of the parameter θs gradually increases along the direction of the Yellow River.

The higher the content of fine particles is, the stronger the adsorption capacity to water is, and the greater the value of the parameter θr is [20, 81]. Combining the particle size division of the Chinese Loess Plateau and the spatial trend analysis of the parameter θr, it can be concluded that the loess in the southern contains more clay and powder particles, and the northern part contains more sand particles. Therefore, the parameter θr has a larger value in the southern, and it has smaller value in the northern of Chinese Loess Plateau. It is basically consistent with the predicted map of the parameter θr.

Calculation method of unsaturated permeability coefficient

By applying the hydraulic conductivity Eq (12) [18] and combining with the statistical trend data to deduce the unsaturated permeability coefficient of a certain place. (12) Where k is the unsaturated permeability coefficient; ks is the saturated permeability coefficient; Se is the effective saturation, it can be calculated by Eq (13). (13) Where, ψ is the matrix suction, kPa; α and n are VG model fitting parameters.

It can be seen from the formula that only α, n, ψ and ks are required to derive the unsaturated permeability coefficient k.

Import the matrix suction and saturation data that have been tested in the Loess Plateau region into Origin for data transformation, draw a scatter plot (Figs 79), and then perform linear fitting to obtain the fitting formulas between the matrix suction and saturation in three different regions of Chinese Loess Plateau (1416). (14) (15) (16) Where, ψ is the matrix suction, kPa; θ is the volumetric moisture content, cm3/cm3.

thumbnail
Fig 7. Sandy loess (Zone I) degree of saturation and matric suction fitting result.

https://doi.org/10.1371/journal.pone.0278307.g007

thumbnail
Fig 8. Typical loess (Zone II) degree of saturation and matric suction fitting result.

https://doi.org/10.1371/journal.pone.0278307.g008

thumbnail
Fig 9. Clay loess(Zone III) degree of saturation and matric suction fitting result.

https://doi.org/10.1371/journal.pone.0278307.g009

If you want to obtain the unsaturated permeability coefficient of a certain sample in the Loess Plateau, read the value of α, n, θs, θr at the certain place on the prediction map firstly. Then, select one formula (14–16) according to the certain place in which zone in Chinese Loess Plateau and calculate the matrix suction ψ. Thirdly, calculate m through n and substitute it into the formula (13) to obtain the effective saturation Se. After measuring the saturated permeability coefficient ks of the sample and using the formula (12), we can gain the unsaturated permeability coefficient k.

This method is simple and economical, compared with fractal method and artificial intelligence algorithm, this method has a wider application range and can reflect the changes of parameters in the region. And it saves a lot of computational troubles (for example, the fractal method needs to analyze the more troublesome pore and particle characteristics), and the data can be obtained quickly. It can provide a preliminary guidance for engineering construction on a large scale. But due to the various measurement methods of the sampled data and different experimental standards, the loess samples will be disturbed to different degrees, so some values will have a certain degree of error. For the parameter α that is greatly affected by the experiment and the void ratio, it is recommended to use the formula for calculating the dry density of each area to assist in the calculation of the parameter α in the actual application. For example, Ji [84] gave the formula for estimating the parameter α of dry density in the Heifangtai area, Luo [85] gave the formula for estimating the parameter α of dry density in Binxian, Shaanxi, and Zhu [86], Ren [87], similar research has also been done. In the actual situation, we can refer to more literature and choose the best. If the testing data of VG model parameters in more and more areas we have, the prediction accuracy of the unsaturated permeability coefficient will be further improved. In order to solve the engineering problems related to loess and water better, and give the better guidance on engineering construction, future research should further focus on the prediction accuracy.

Calculation example verification

The measured results of other scholars’ unsaturated permeability coefficients in different zones are used to verify the prediction accuracy of this method, which is shown in Table 4.

thumbnail
Table 4. Unsaturated permeability coefficient data source.

https://doi.org/10.1371/journal.pone.0278307.t004

thumbnail
Table 5. The prediction data of VG model parameter by contour graph.

https://doi.org/10.1371/journal.pone.0278307.t005

By using the proposed method, the unsaturated permeability coefficient has been predicted in these three research areas. And compared the predict results with the experimental test data, the results are been shown in Figs 1012 and Table 6.

thumbnail
Fig 10. The predicted unsaturated permeability coefficient graph in Zone III(Yangling).

https://doi.org/10.1371/journal.pone.0278307.g010

thumbnail
Fig 11. The predicted unsaturated permeability coefficient graph in Zone II(Yan’an).

https://doi.org/10.1371/journal.pone.0278307.g011

thumbnail
Fig 12. The predicted unsaturated permeability coefficient graph in Zone I(Jingbian).

https://doi.org/10.1371/journal.pone.0278307.g012

thumbnail
Table 6. The prediction deviation of the prediction results.

https://doi.org/10.1371/journal.pone.0278307.t006

According to Figs 1012 and Table 6, it can be seen that the prediction results are basically accurate. The three test groups in the Zone III have the same saturated permeability coefficient, so they only have one predicted result.

It can be found in Fig 10 and Table 6 that the prediction results are basically in line with the experimental results, and the prediction deviation is basically controlled within ±4.77e- 7 cm/s. When the saturation of soil sample is low, the deviation is small, and it is basically controlled at ±2.22e-8cm/s. When the saturation of soil sample is high, the unsaturated permeability coefficient is large, and the prediction deviation is also large at lower saturation stage, but it is also controlled at ±2.29e-6cm/s. And the transition between the predicted result and the saturated permeability coefficient is relatively smooth, which is more in line with the actual situation.

In Fig 11 and Table 6, comparing the prediction results with the experimental test results, the prediction results are basically accurate, and the prediction deviation is basically controlled within ±2.17e-8cm/s. The predicted results in the lower saturation section of soil sample are in good agreement with the experimental results, with a deviation of ±1.06e-8cm/s. When the saturation of soil sample is higher, the unsaturated permeability coefficient becomes larger, and the prediction deviation is larger than soil sample in lower saturation stage, which is basically controlled at ±3.08e-7 cm/s, and it can be seen that in the experimental test results, some unsaturated permeability coefficient rebounded. The results of this situation are generally caused by experimental deviations, and the prediction results can make up the test deviation.

In Fig 12 and Table 6, since the unsaturated permeability coefficient of sandy loess in Zone I is the largest, the difference between the predicted results and the experimental results is basically ±3.16e-6. However, the prediction results are extremely accurate at low saturation of soil sample, and the deviation is basically ±1.94e-10cm/s. In the higher saturation section of soil sample, the prediction deviation is also larger due to the higher saturation and the influence of the higher permeability coefficient of sandy loess, with a deviation of ±1.44e-5cm/s. For the larger unsaturated permeability coefficient in the sandy loess area, the deviation is within the acceptable range.

By predicting the unsaturated permeability coefficients of the three regions, it can be concluded that: for approximately estimating the unsaturated permeability coefficient of a certain place, the unsaturated permeability coefficient of a certain place in the Loess Plateau can be obtained quickly and economically only by a relatively simple experiment proposed in this paper. The prediction result of the method of permeability coefficient is relatively accurate, the deviation is very small in the lower saturation section, and the deviation is slightly larger in the higher saturation section, but it is also within the acceptable range.

Conclusion

The SWCCs reflects the relationship between water content and matrix suction. It has been widely used to estimate the loess properties. The unsaturated permeability coefficient of loess is the key to solving loess engineering problems. Therefore, research on the SWCCs fitting model is very important. In view of this, the VG model parameters of the SWCCs data in the Loess Plateau region studied by experts and scholars in recent years are collected in this paper, and trend analysis and one-variable linear regression methods are used to study the regional distribution law of the VG model parameter and the calculated method on unsaturated permeability coefficient is developed. The conclusions are shown as the following.

  1. The parameters of the VG model in the Loess Plateau have obvious regularities:
    1. ①For the parameter α, the value is high in the area with high fine particle content, such as the southeast of the Chinese Loess Plateau, and it is low in the area with high gross particle content, like, in the northwest. There are abnormal values in some areas such as Zichang of Shaanxi province and Xining of Qinghai province in Zone II.
    2. ②For the parameter n, the value is low in the areas with high fine particle content, for example, in the southeast, and it is high in the areas with high gross particle content, like in the northwest.
    3. ③For the parameter θs, it is high in the southeast and low in the northwest, and the high value is basically distributed on both banks of the Yellow River. It is speculated that a large amount of fine particle matter is brought about by the alluvium of the river.
    4. ④For the parameter θr, the value is large in the south and small in the north. And this trend is mainly related to the average particle size of loess.
  2. Through linear fitting of matrix suction and saturation, the linear formulas of them are obtained in three different zones on the Chinese Loess Plateau. As long as the ks of a loess sample in a certain area are obtained, the unsaturated permeability coefficient of the sample can be derived by reading the α, n, θs, θr of the certain place on the prediction map and using these fitting formulas. The prediction result of the method of permeability coefficient is relatively accurate, the deviation is very small in the lower saturation section, and the deviation is slightly larger in the higher saturation section, but it is also within the acceptable range. This method can greatly reduce the experiment cost. It will provide convenience for solving the engineering problems related to loess and water and other engineering applications.

Supporting information

S1 Table. The data summary of Van Genuchten model parameters on soil-water characteristic curves in Chinese Loess Plateau.

https://doi.org/10.1371/journal.pone.0278307.s001

(XLS)

Acknowledgments

Our sincere thanks to Pro. Huantong Li of Xi’an University of Science and Technology for advice of figures made and Mr. Yang Li of Xi’an University of Science and Technology for their technical support during the preparation of the manuscript. Finally, we would like to thank anonymous reviewers and an academic editor for their valuable comments and suggestions which helped us on improve this paper.

References

  1. 1. Mualem Y. Hydraulic conductivity of unsaturated soils: Prediction and formulas. Methods of Soil Analysis: Part 1 Physical and Mineralogical Methods. 1986;5:799–823.
  2. 2. Gardner W. Mathematics of isothermal water conduction in unsaturated soils. Highway research board special report. 1958;40:78–87.
  3. 3. Davidson J, Stone L, Nielsen D, Larue M. Field measurement and use of soil©\water properties. Water Resources Research. 1969;5(6):1312–21.
  4. 4. Campbell JD. Pore pressures and volume changes in unsaturated soils [PhD diss.]: University of Illinois at Urbana-Champaign; 1973.
  5. 5. Richards LA. Capillary conduction of liquids through porous mediums. Physics. 1931;1(5):318–33.
  6. 6. Christensen HR. Permeability-capillary potential curves for three prairie soils. Soil Science. 1944;57(5):381–90.
  7. 7. Wind G. A field experiment concerning capillary rise of moisture in a heavy clay soil. Netherlands Journal of Agricultural Science. 1955;3(1):60–9.
  8. 8. Rijtema PE. An analysis of actual evapotranspiration: Pudoc; 1965.
  9. 9. Philip J. Linearized unsteady multidimensional infiltration. Water resources research. 1986;22(12):1717–27.
  10. 10. Leong EC, Rahardjo H. Permeability functions for unsaturated soils. Journal of geotechnical and geoenvironmental engineering. 1997;123(12):1118–26.
  11. 11. Xue H, Fan Q, Shu C. Prediction of micro-channel flows using direct simulation Monte Carlo. Probabilistic Engineering Mechanics. 2000;15(2):213–9.
  12. 12. Samingan AS, Leong E-C, Rahardjo H. A flexible wall permeameter for measurements of water and air coefficients of permeability of residual soils. Canadian geotechnical journal. 2003;40(3):559–74.
  13. 13. GENUCHTEN MTV. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Science Society of America Journal. 1980;44(5):892–8.
  14. 14. Childs EC, Collis-George N. The permeability of porous materials. Proceedings of the Royal Society of London Series A Mathematical and Physical Sciences. 1950;201(1066):392–405.
  15. 15. Marshall T. A relation between permeability and size distribution of pores. Journal of soil science. 1958;9(1):1–8.
  16. 16. Kunze R, Uehara G, Graham K. Factors important in the calculation of hydraulic conductivity. Soil Science Society of America Journal. 1968;32(6):760–5.
  17. 17. Brooks RH. Hydraulic properties of porous media. Hydrology Papers. 1964(24):37.
  18. 18. Mualem Y. A new model for predicting the conductivity of unsaturated porous media. Water Resour Res. 1976(12):513–22.
  19. 19. Vanapalli S.K., Fredlund DG, Pufahl DE. Model for the prediction of shear strength with respect to soil suction. Canadian Geotechnical Journal. 1996;33(3):379–92.
  20. 20. Vanapalli S.K., Fredlund DG, Pufahl DE. The influence of soil structure and stress history on the soil-water characteristics of a compacted till. Geotechnique. 1999;49(2):143–59.
  21. 21. Lu N, Likos WJ. UNSATURATED SOIL MECHANICS. New Jersey: John Wiley & Sons, Inc.; 2004.
  22. 22. Li P, Li TL, Vanapalli S.K.. Prediction of soil-water characteristic curve for Malan loess in Loess Plateau of China. Journal of Central South University. 2018;25(2):432–47.
  23. 23. Nie YP, Ni WK, Li XN, Wang HM, Yuan KZ, Tuo WX. The Influence of Drying-Wetting Cycles on the Suction Stress of Compacted Loess and the associated Microscopic Mechanism. Water. 2021;13(13).
  24. 24. Feng CX, Tian B, Lu XC, Beer M, Broggi M, Bi SF, et al. Bayesian Updating of Soil–Water Character Curve Parameters Based on the Monitor Data of a Large-Scale Landslide Model Experiment. Applied Sciences. 2020;10(16).
  25. 25. Birang NE, Fredlund DG. Assessment of the WP4-T Device for Measuring Total Suction. Geotechnical Testing Journal. 2016;39(3):500–6.
  26. 26. Rahardjo H, Dong XF, Lee DTT, Leong EC, Fong Y. Expedited Soil-Water Characteristic Curve Tests Using Combined Centrifuge and Chilled Mirror Techniques. Geotechnical Testing Journal. 2017;41(1):207–17.
  27. 27. Arya LM, Paris JF. A physicoempirical model to predict the soil moisture characteristic from particle©\size distribution and bulk density data. Soil Science Society of America Journal. 1981;45(6):1023–30.
  28. 28. Haverkamp R, Parlange J-Y. Predicting the water-retention curve from particle-size distribution: 1. Sandy soils without organic matter1. Soil Science. 1986;142(6):325–39.
  29. 29. Simms P, Yanful E. A pore-network model for hydromechanical coupling in unsaturated compacted clayey soils. Canadian Geotechnical Journal. 2005;42(2):499–514.
  30. 30. Hu R, Chen YF, Liu HH, Zhou CB. A water retention curve and unsaturated hydraulic conductivity model for deformable soils: consideration of the change in pore-size distribution. Geotechnique. 2013;63(16):1389–405.
  31. 31. Fredlund DG, Xing A. Equations for the soil-water characteristic curve. Canadian geotechnical journal. 1994;31(4):521–32.
  32. 32. Huang LM, Shao MA. Advances and perspectives on soil water research in China’s Loess Plateau. Earth-Science Reviews. 2019;199.
  33. 33. Gu TF, Zhang MS, Wang JD, Wang CX, Xu YJ, Wang X. The effect of irrigation on slope stability in the Heifangtai Platform, Gansu Province, China. Engineering Geology. 2019;248:346–56.
  34. 34. Deng J, Chong YJ, Zhang D, Ren CJ, Zhao FZ, Zhang XX, et al. Temporal Variations in Soil Enzyme Activities and Responses to Land-Use Change in the Loess Plateau, China. Applied Sciences. 2019;9(15).
  35. 35. Zhao QZ, Ma XW, Liang L, Yao WQ. Spatial–Temporal Variation Characteristics of Multiple Meteorological Variables and Vegetation over the Loess Plateau Region. Applied Sciences. 2020;10(3).
  36. 36. Wang HF, Xiao B, Wang MY, Shao MA. Modeling the Soil Water Retention Curves of Soil-Gravel Mixtures with Regression Method on the Loess Plateau of China. PLoS One. 2013;8(3):e59475. pmid:23555040
  37. 37. Zhang FX, Wilson GW, Fredlund DD. Estimation of Permeability Function for Bulyanhulu Tailings. Geotechnical Testing Journal. 2019;42(3):670–86.
  38. 38. Fan S, Zhang XC, Pei XJ, Zhang S, Chang ZL, Cui XT. Experimental study on soil-water characteristic curve of undisturbed loess in different regions and its application to slope stability evaluation. Water Resources and Hydropower Engineering. 2019;50(1):154–61.
  39. 39. Zeng L. Study on the Influence of the prior infiltration channels on loess slope stability-Take the Lu Caotai landslide in Shaanxi as an example [Master’s thesis]: Chang’an University, Xi’an, Shaanxi, China; 2011.
  40. 40. Zhao AH. Hydraulic Parameters Measurements and Estimations for Loess Soils [Master’s thesis]: Northwest A & F University, Yangling, Shaanxi, China; 2008.
  41. 41. Yang JY. The Establishment and Evaluation of Soil Fractal Characteristics and Pedo-transfer Functions in Hetao Irrigation District [Master’s thesis]: Inner Mongolia Agricultural University, Hohhot, Inner Mongolia, China; 2007.
  42. 42. Wang JW. Study on the Soil Water Characteristic Curve and Strength Characteristics of Unsaturated Loess in East Gansu Province [Master’s thesis]: Lanzhou University of Technology, Lanzhou, Gansu, China; 2017.
  43. 43. Nie KY. Analysis on the Steady Seepage of Unsaturated Soil and the Optimization of Slope Stability Numerical Calculation Method [Master’s thesis]: Northwest A&F University, Yangling, Shaanxi, China; 2018.
  44. 44. Chen C. Experimental Study on Mechanical and Hydraulic Characteristics of Compacted and Lime-Modified Loess [Master’s thesis]: Xi’an University of Technology, Xi’an, Shaanxi, China; 2019.
  45. 45. Qiao PT, Tang YM, Zhang CL, Kang HW, Zhou XL. Catastrophic Analysis of Loess Landslide Considering Unsaturated Seepage under Rainfall Condition. Henan Science. 2019;37(2):255–61.
  46. 46. Shi LY. Experimental study and Numerical analysis on Rainfall Infiltration and Drip Irrigation of Fractured Loess Slope [Master’s thesis]: Chengdu University of Technology, Chengdu, Sichuan, China; 2017.
  47. 47. Yang H. Experimental Research on the Permeation Properties of fissured loess [Master’s thesis]: Chang’an University, Xi’an, Shaanxi, China; 2016.
  48. 48. Jiang C-X. The unsaturated characteristic and mine high slope stability research of Malan loess in Shanzi region [PhD diss.]: China University of Geosciences, Wuhan, Hubei, China; 2016.
  49. 49. Wang JE, Xiang W, Zuo X. Situation and Prevention of Loess Water Erosion Problem along the West-to-East Gas Pipeline in China. Journal of Earth Science. 2010;21(6):968–73.
  50. 50. Bi J. Research on moisture field of unsaturated loess foundation by numerical analysis [Master’s thesis]: Lanzhou University, Lanzhou, Gansu, China; 2015.
  51. 51. Zhang PY. Research on Strength Evolution Mechanism of Remolded Unsaturated Loess and Slope Stability under Drying-wetting Cycles [Master’s thesis]: Lanzhou Jiaotong University, Lanzhou, Gansu, China; 2019.
  52. 52. Jiang Y. Study of the characteristic curves for unsaturated loess [Master’s thesis]: Lanzhou University, Lanzhou, Gansu, China; 2012.
  53. 53. Liu P. Laboratory Experimental Study on Unsaturated Permeability of Remolded Loess [Master’s thesis]: Lanzhou University, Lanzhou, Gansu, China; 2019.
  54. 54. Cao CW. Study on unsaturated infiltration characteristics of loess in Heifangtai, Gansu Province [Master’s thesis]: Chengdu University of Technology, Chengdu, Sichuan, China; 2018.
  55. 55. Gao ZA, Li P, Xiao JJ, Li TL. Evaluation of Shear Strength of Unsaturated Loess Using Direct Shear Test. Journal of Engineering Geology. 2019;28(2):344–51.
  56. 56. Liu PF, Yin YP, Li B, Sun YB. The tests of soil-water characteristical curve and permeability coefficient prediction of the wetting-drying cycles for unsaturated loess. The Chinese Journal of Geological Hazard and Control. 2015;26(4):125–9.
  57. 57. Li P. Study on the regularities of water migration in unsaturated loess [Master’s thesis]: Chang’an University, Xi’an, Shaanxi, China; 2013.
  58. 58. Li P, Li TL, Wang H, Liang Y. Soil-water characteristic curve and permeability prediction on Childs & Collis-Geroge model of unsaturated loess. Rock and Soil Mechanics. 2013;34(s2):184–9.
  59. 59. Zhang YQ. Experimental Research on Infiltration in Unsaturated Loess [Master’s thesis]: Chang’an University, Xi’an, Shaanxi, China; 2016.
  60. 60. Sun F. Study on Influence of Soil Water Characterisitics on Stability of Loess Slope under Rainfall Conditions [Master’s thesis]: Xi’an University of Technology, Xi’an, Shaanxi, China; 2019.
  61. 61. Zhang X, Zhao WW, Wang LX, Liu YX, Liu Y, Feng Q. Relationship between soil water content and soil particle size on typical slopes of the Loess Plateau during a drought year. Science of The Total Environment. 2019;648:943–54. pmid:30144762
  62. 62. Wang XF, Shi HL, Liu SK. Hydraulic Conductivity Curve of Unsaturated Loess in Zichang Area of Shaanxi Province. Water Resources and Power. 2014;32(4):41–5.
  63. 63. Guo H, Li J, Wang PC, Guo R. Estimation of soil-water characteristic curve for undisturbed loess soil. China Sciencepaper. 2018;13(1):37–40.
  64. 64. Guo H, Chen Q, Chen DL, Li J. Study on soil-water characteristics curve fitting for undisturbed loess soil. Yangtze River. 2016;47(19):92–5.
  65. 65. Wang TX, Lu J, Yue CK. Soil-water characteristic curve for unsaturated loess considering temperature and density effect. Rock and Soil Mechanics. 2008;29(1):41–5.
  66. 66. Yang T. Study on anisotropic permeability of loess considering density and freezing-thawing effect [Master’s thesis]: Xi’an University of Architecture and Technology, Xi’an, Shaanxi, China; 2016.
  67. 67. Han JM. Research on Loess Slope Failure Mechanism and Mitigation Technology under Rainfall Infiltration [PhD diss.]: Xi’an University of Science and Technology, Xi’an, Shaanxi, China; 2015.
  68. 68. Yue LB. Study on Seepage Law and Collapse Characteeistics of Loess Foundation under Rainfall Infiltration [Master’s thesis]: Xi’an University of Technology, Xi’an, Shaanxi, China; 2019.
  69. 69. Qiu HB. Study of Rain Infiltration on Stability of Unsaturated Pit Supporting Structure [Master’s thesis]: Xi’an University of Architecture and Technology, Xi’an, Shaanxi, China; 2013.
  70. 70. Wang Y. The laboratory model test and the numerical simulation of instability on loess slope induced by rainfall infiltration [Master’s thesis]: Wuhan Institute of Technology, Wuhan, Hubei, China; 2018.
  71. 71. Tang FY, Cheng YJ, Bao CJ, Hu JL, Liang Q, Wu Y, et al. Spatio-Temporal Trends and Risk Factors for Shigella from 2001 to 2011 in Jiangsu Province, People’s Republic of China. PLoS One. 2014;9(1):e83487. pmid:24416167
  72. 72. Fang SY, Liu Z. Characteristics of sulfur content of raw coal and main controlling factors in Hujiahe mining field. Journal of Xi’an University of Science and Technology. 2017;37(3):339–45.
  73. 73. Yu LY, Wang XY, Hou ZY, Du ZY, Zeng YF, Mu ZY. Path Planning Optimization for Driverless Vehicle in Parallel Parking Integrating Radial Basis Function Neural Network. Applied Sciences. 2021;11(17).
  74. 74. Xie CC, Jia SJ, Li YJ, An C, Yang C. Mapping Based Quality Metrics for Mesh Deformation Algorithms Using Radial Basis Functions. Applied Sciences. 2021;11(1).
  75. 75. Chen XQ. Interpolation of Flow Measurement Section based on Radial Basis Function (RBF) Neural Network. Journal of Changjiang Institute of Technology. 2020;37(2):30–2.
  76. 76. Qiao YY, Wu JL, Feng XL. The Theory Method and Application of MQ Radial Basis Function. Journal of Xinjiang University(Natural Science Edition). 2015;32(4):379–87.
  77. 77. Kowalczyk K, Kowalczyk AM, Chojka A. Modeling of the Vertical Movements of the Earth’s Crust in Poland with the Co-Kriging Method Based on Various Sources of Data. Applied Sciences. 2020;10(9).
  78. 78. Arshid MU, Kamal MA. Regional Geotechnical Mapping Employing Kriging on Electronic Geodatabase. Applied Sciences. 2020;10(21).
  79. 79. Yang GL, Zhang GM, Li SX. Application of universal Kriging interpolation in geomagnetic map. Journal of Chinese Inertial Technology. 2008(2):162–6.
  80. 80. Panday D, Maharjan B, Chalise D, Shrestha RK, Twanabasu B. Digital soil mapping in the Bara district of Nepal using kriging tool in ArcGIS. PLoS One. 2018;13(10):e0206350. pmid:30365521
  81. 81. Bai XL, Zhang B, Wang HX. Experiment and SWCC fitting study of unsaturated and remolded loess. Mineral Resources and Geology. 2020;34(5):988–94. https://doi.org/0.19856/j.cnki.issn.1001-5663.2020.05.019
  82. 82. Luo CC. A New method for simultaneous measurement of soil-water characteristic curves and permeability coefficient for unsaturated soils [Master’s thesis]: Hubei University of Technology, Wuhan, Hubei, China; 2020.
  83. 83. Cai GQ, Han BW, Yang Y, Liu W, Zhao CG. Experimental study on soil-water characteristic curves of sandy loess. Chinese Journal of Geotechnical Engineering. 2020;42(s1):11–5.
  84. 84. Ji XP, Hou XK. The Effect of Degree of Density on the Soil Water Characteristic Curve of Q3 Loess. Journal of Gansu Sciences. 2021;33(3):136–41.
  85. 85. Luo Y. Experimental Study on Joint Seepage Characteristics of Loess and Numerical Analysis of Joint Loess [Master’s thesis]: Xi’an University of Architecture and Technology, Xi’an, Shaanxi, China; 2011.
  86. 86. Zhu LJ. Study of Effects of Rain Infiltration on Stability of Unsaturated Loess Slope [Master’s thesis]: Xi’an University of Architecture and Technology, Xi’an, Shaanxi, China; 2008.
  87. 87. Ren HB. Numerical analysis on moisture field of loess foundation under trench seepage [Master’s thesis]: Xi’an University of Architecture and Technology, Xi’an, Shaanxi, China; 2008.
  88. 88. Hu HJ, Li CH, Cui YJ, Shao ZX, Ai YD. Research on the determination of permeability coefficient of unsaturated remolded loess under wetting condition. Journal of Hydraulic Engineering. 2018;49(10):1216–26.
  89. 89. Wang HM, Ni WK. Prediction model of saturated/unsaturated permeability coefficient of compacted loess with different dry density. Rock and Soil Mechanics. 2022;43(03).
  90. 90. Yang Y. Experimental investigation on water permeability of unsaturated soils under hydro-mechanical coupled conditions [Master’s thesis]: Beijing Jiaotong University; 2019.