**1. Introduction**

Desertification is an environmental issue of global concern, especially in arid and semiarid regions. Artificial afforestation has been considered an effective ecological means for combating desertification in many arid desert regions worldwide [1,2]. Great challenges, however, have appeared when the afforestation is conducted in arid desert regions due to the lack of freshwater and extreme environmental conditions, including severe drought [1,3,4]. Due to less rainfall in arid desert regions, water scarcity has become a worldwide issue of increasing severity [1,5]. The lower-quality saline–alkaline groundwater is widely applied [6–8]. Unfortunately, saline water irrigation normally leads to greater salinity hazards to plant growth and survival in groundwater extraction [9,10]. Therefore,

**Citation:** Liu, J.; Zhao, Y.; Zhang, J.; Hu, Q.; Xue, J. Effects of Irrigation Regimes on Soil Water Dynamics of Two Typical Woody Halophyte Species in Taklimakan Desert Highway Shelterbelt. *Water* **2022**, *14*, 1908. https://doi.org/10.3390/ w14121908

Academic Editor: Guido D'Urso

Received: 13 April 2022 Accepted: 10 June 2022 Published: 14 June 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

establishing a suitable irrigation regime for artificial vegetation growth and survival is crucial to saving groundwater utilization and reducing the salinity hazards.

The Taklimakan Desert, called the "Dead Sea", is the second-largest mobile desert in the world. To improve transportation for the exploitation of petroleum resources, the Taklimakan Desert Highway, the longest highway across a shifting desert in the world, was completed in 1995 [11]. To overcome the frequent sand burial of the highway, the Taklimakan Desert Highway shelterbelt was constructed through a biological engineering project in 2003 [12]. The mobile dunes on both sides of the highway were effectively stabilized by introducing drought- and salt-tolerant plants [13], such as *Haloxylon Bunge*, *Calligonum Linn*, and *Tamarix Linn.*

Drip irrigation from saline groundwater is one of the most efficient ways to support artificial shelterbelts [14]. Owing to the differences in the adaptability or adaptive strategies of drought- and salt-tolerant plants, irregular or insufficient irrigation will lead to different responses of plants to water stresses [14]. Therefore, it is crucial to determine suitable irrigation regimes to ensure plant survival in the drip irrigation process. *Haloxylon ammodendron* and *Calligonum mongolicunl* are the two main species in the Taklimakan Desert Highway shelterbelt. It is reported that moderate irrigation intervals are beneficial to the growth of the two species in the Taklimakan Desert Highway shelterbelt since they can save water and support the plants' water demands [14].

In recent years, many researchers have studied the soil water dynamics and irrigation regimes of desert plants. For example, Ding et al. reported that the soil moisture at 0–120 cm depth presents an apparent single-peak curve. The salt accumulation phenomenon is evident at 45–60 cm, and the salt content reaches 10–20 g kg−<sup>1</sup> in the Taklimakan Desert [15]. Li et al. found that saltwater irrigation did not produce salt stress on the plant roots in the Taklimakan Desert [1]. *Haloxylon ammodendron*'s roots are mainly distributed between 20 and 80 cm, while the salt is mainly concentrated in the 0–20 cm surface layer. On the contrary, saline water irrigation is beneficial for increasing soil nutrients. Fu et al. pointed out that *Haloxylon ammodendron* mainly utilizes shallow soil moisture (20–40 cm) and deep soil moisture (100–350 cm) and underground water in May, but deep soil moisture (160–350 cm) and underground water in August in the southern edge of Gurbantunggut Desert [16]. Zhang et al. indicated that the structural heterogeneity of the soil layer has a retarding effect on water content [17]. The soil layers with more clay and silt particles are more prone to salt accumulation at the southern edge of the Gurbantunggut Desert [17].

The previous studies focused mainly on saline water irrigation and its influence on soil properties and plant growth [6,18,19]. However, the effects of saline water irrigation regimes on the soil water dynamics of desert plants have been ignored. Meanwhile, the number of irrigation intervals and periods of different desert plants remain unknown. This study aims to (1) examine the effects of irrigation regimes on the soil water dynamics of two typical woody halophyte species and (2) quantify the irrigation intervals and periods of the two species based on a field test of precision irrigation control. This will provide the theoretical basis for developing water-saving irrigation measures in the shelterbelt.

#### **2. Materials and Methods**

#### *2.1. Study Area*

This study was carried out in the Taklimakan Desert Highway shelterbelt, which was built from Xiaotang to Minfeng, being 436 km long and 72–78 m wide (Figure 1a). It is characterized by an extremely high temperature, less rainfall, and strong evaporation. According to the Tazhong meteorological station (83◦48014.16900 E,38◦5402.03800 N), the average annual temperature in this area is 12.4 ◦C. The extreme minimum temperature is −22.2 ◦C, and the extreme maximum temperature reaches 45.6 ◦C. The average annual rainfall is 24.6 mm, while the average annual evaporation is 3639 mm. The average relative humidity is only 29.4% (Figure 1b). The average annual wind speed is 2.5 m s−<sup>1</sup> , and the maximum instantaneous wind speed reaches 20 m s−<sup>1</sup> [15]. The soil type is mobile aeolian

*Water* **2022**, *14*, x FOR PEER REVIEW 3 of 13

*Water* **2022**, *14*, x FOR PEER REVIEW 3 of 13

soil, and the soil salt content is 1.26–1.63 g kg−<sup>1</sup> [19]. The soil physical properties in this area have been described by Zhang et al. [2]. soil, and the soil salt content is 1.26–1.63 g kg−1 [19]. The soil physical properties in this area have been described by Zhang et al. [2]. soil, and the soil salt content is 1.26–1.63 g kg−1 [19]. The soil physical properties in this area have been described by Zhang et al. [2].

is 24.6 mm, while the average annual evaporation is 3639 mm. The average relative humidity is only 29.4% (Figure 1b). The average annual wind speed is 2.5 m s−1, and the maximum instantaneous wind speed reaches 20 m s−1 [15]. The soil type is mobile aeolian

is 24.6 mm, while the average annual evaporation is 3639 mm. The average relative humidity is only 29.4% (Figure 1b). The average annual wind speed is 2.5 m s−1, and the maximum instantaneous wind speed reaches 20 m s−1 [15]. The soil type is mobile aeolian

**Figure 1.** Taklimakan Desert Highway (**a**), and atmospheric temperature (Ta), air humidity (RH), and wind speed (WS) at the height of 2 m in the middle of the Taklimakan Desert in 2016 (**b**). **Figure 1.** Taklimakan Desert Highway (**a**), and atmospheric temperature (Ta), air humidity (RH), and wind speed (WS) at the height of 2 m in the middle of the Taklimakan Desert in 2016 (**b**). **Figure 1.** Taklimakan Desert Highway (**a**), and atmospheric temperature (Ta), air humidity (RH), and wind speed (WS) at the height of 2 m in the middle of the Taklimakan Desert in 2016 (**b**).

The vegetation community is very sparse, and most areas have no vegetation [18]. *Haloxylon ammodendron*, *Calligonum mongolicum*, and *Tamarix L* are the main three species in the Taklimakan Desert Highway shelterbelt (Figure 2). Saline groundwater was used for drip irrigation. The salinity of irrigation water was 4.03 g L−1; the irrigation period was one time every 10 days in July and August, and 15 days in other months, while there was no irrigation in winter (from November to February of the following year). The irrigation amount was 35 mm each time [11]. The vegetation community is very sparse, and most areas have no vegetation [18]. *Haloxylon ammodendron*, *Calligonum mongolicum*, and *Tamarix L* are the main three species in the Taklimakan Desert Highway shelterbelt (Figure 2). Saline groundwater was used for drip irrigation. The salinity of irrigation water was 4.03 g L−<sup>1</sup> ; the irrigation period was one time every 10 days in July and August, and 15 days in other months, while there was no irrigation in winter (from November to February of the following year). The irrigation amount was 35 mm each time [11]. The vegetation community is very sparse, and most areas have no vegetation [18]. *Haloxylon ammodendron*, *Calligonum mongolicum*, and *Tamarix L* are the main three species in the Taklimakan Desert Highway shelterbelt (Figure 2). Saline groundwater was used for drip irrigation. The salinity of irrigation water was 4.03 g L−1; the irrigation period was one time every 10 days in July and August, and 15 days in other months, while there was no irrigation in winter (from November to February of the following year). The irrigation amount was 35 mm each time [11].

**Figure 2.** Three main species in Taklimakan Desert Highway shelterbelt (adapted from Zhang et al. **Figure 2.** Three main species in Taklimakan Desert Highway shelterbelt (adapted from Zhang et al. **Figure 2.** Three main species in Taklimakan Desert Highway shelterbelt (adapted from Zhang et al. [2]).

#### [2]). *2.2. Experiment Design and Data Processing*

[2]).

*2.2. Experiment Design and Data Processing*  This experiment selected six rows of well-growing trees in the Taklimakan Desert Highway shelterbelt for irrigation treatment, with each row of approximately 100 m. *Ha-2.2. Experiment Design and Data Processing*  This experiment selected six rows of well-growing trees in the Taklimakan Desert Highway shelterbelt for irrigation treatment, with each row of approximately 100 m. *Ha-*This experiment selected six rows of well-growing trees in the Taklimakan Desert Highway shelterbelt for irrigation treatment, with each row of approximately 100 m. *Haloxylon ammodendron* and *Calligonum mongolicum* were planted 8 years ago. The same irrigation amount was adopted, and a small switch was installed on the drip irrigation pipe to control the irrigation period. The distance between two rows of plants was at least 5 m. Two rows of trees were planted in each row with a 1 m spacing. Drip irrigation pipes were laid close to the trunk, and the water outlet holes were spaced by 1 m. Three irrigation levels were set, W1 = 17.5 mm, W2 = 25 mm, and W3 = 35 mm, and the three irrigation

periods were F1 = 10 d, F2 = 20 d, F3 = 40 d. The samples were taken from different plants in each treatment plot. Three plants were selected from each treatment plot and one sample was taken from each plant for a total of three repetitions. Moreover, in the experiment design, the W1F1, W1F2, W1F3, W3F1, and W3F2 treatments were selected in order to analyze the effects of different combinations of irrigation amount and irrigation period on soil moisture under the same total amount of irrigation over 40 days. The total irrigation amount of the W1F1 and W3F2 treatments in 40 days was 70 mm. The total irrigation amount of the W1F2 and W3F1 treatments was 35 mm in 40 days, and that of the W1F3 treatment was 17.5 mm in 40 days. The combination of irrigation amount and irrigation period for W1 and W3 only allowed the same combination of total irrigation amount. In contrast, the irrigation level of W2 was different from the total irrigation amount of other treatments within 40 days, so the W2 treatment was not selected in this study. The specific field configuration is shown in Figure 3. periods were F1 = 10 d, F2 = 20 d, F3 = 40 d. The samples were taken from different plants in each treatment plot. Three plants were selected from each treatment plot and one sample was taken from each plant for a total of three repetitions. Moreover, in the experiment design, the W1F1, W1F2, W1F3, W3F1, and W3F2 treatments were selected in order to analyze the effects of different combinations of irrigation amount and irrigation period on soil moisture under the same total amount of irrigation over 40 days. The total irrigation amount of the W1F1 and W3F2 treatments in 40 days was 70 mm. The total irrigation amount of the W1F2 and W3F1 treatments was 35 mm in 40 days, and that of the W1F3 treatment was 17.5 mm in 40 days. The combination of irrigation amount and irrigation period for W1 and W3 only allowed the same combination of total irrigation amount. In contrast, the irrigation level of W2 was different from the total irrigation amount of other treatments within 40 days, so the W2 treatment was not selected in this study. The specific field configuration is shown in Figure 3.

*loxylon ammodendron* and *Calligonum mongolicum* were planted 8 years ago. The same irrigation amount was adopted, and a small switch was installed on the drip irrigation pipe to control the irrigation period. The distance between two rows of plants was at least 5 m. Two rows of trees were planted in each row with a 1 m spacing. Drip irrigation pipes were laid close to the trunk, and the water outlet holes were spaced by 1 m. Three irrigation levels were set, W1 = 17.5 mm, W2 = 25 mm, and W3 = 35 mm, and the three irrigation

*Water* **2022**, *14*, x FOR PEER REVIEW 4 of 13

**Figure 3.** Schematic diagram of field experiment design. (C. means *Calligonum mongolicum*; H. means *Haloxylon ammodendron*; hereinafter the same*.*) **Figure 3.** Schematic diagram of field experiment design. (C. means *Calligonum mongolicum*; H. means *Haloxylon ammodendron*; hereinafter the same).

Before the test, three water outlet holes were selected under the drip irrigation pipe in each plot, and a measuring cylinder was placed under each outlet hole to measure the average water output per unit time of each drip irrigation pipe. From 10 August to 20 September 2015, soil samples drilled at 30 cm from the root in each cell were taken at a depth of 0–200 cm (0–10 cm, 10–20 cm, 20–40 cm, 40–60 cm, 60–80 cm, 80–100 cm, 100–120 cm, 120–140 cm, 140–160 cm, 160–180 cm, 180–200 cm) on the 1st, 4th, 9th, 14th, 19th, 29th, and 39th days after irrigation, to evaluate the differences in soil moisture changes under different irrigation strategies with the same irrigation amount of 40 days. Before the test, three water outlet holes were selected under the drip irrigation pipe in each plot, and a measuring cylinder was placed under each outlet hole to measure the average water output per unit time of each drip irrigation pipe. From 10 August to 20 September 2015, soil samples drilled at 30 cm from the root in each cell were taken at a depth of 0–200 cm (0–10 cm, 10–20 cm, 20–40 cm, 40–60 cm, 60–80 cm, 80–100 cm, 100–120 cm, 120–140 cm, 140–160 cm, 160–180 cm, 180–200 cm) on the 1st, 4th, 9th, 14th, 19th, 29th, and 39th days after irrigation, to evaluate the differences in soil moisture changes under different irrigation strategies with the same irrigation amount of 40 days.

To analyze the change in soil moisture in the whole plot in different months, a 1-m long aluminum tube instrument with LNW-50A neutron probes (CAS, Nanjing, Jiangsu, CHN, 1986) was used for measurement. The neutron tubes were buried in the middle of each plot in 2016. The horizontal distance between the neutron tube and the dropper was 70 cm, and the buried depth was 320 cm. Before the measurement, the neutron instrument was calibrated in layers (0–40 cm and 40–300 cm). It was drilled at a horizontal distance of 30 cm from the dropper to a depth of 200 cm with each 20 cm layer. Each day before irrigation was selected to measure the soil water moisture by drilling and neutron meter in May, July, August, and September of 2016. To analyze the change in soil moisture in the whole plot in different months, a 1-m long aluminum tube instrument with LNW-50A neutron probes (CAS, Nanjing, Jiangsu, CHN, 1986) was used for measurement. The neutron tubes were buried in the middle of each plot in 2016. The horizontal distance between the neutron tube and the dropper was70 cm, and the buried depth was 320 cm. Before the measurement, the neutron instrumentwas calibrated in layers (0–40 cm and 40–300 cm). It was drilled at a horizontal distanceof 30 cm from the dropper to a depth of 200 cm with each 20 cm layer. Each day beforeirrigation was selected to measure the soil water moisture by drilling and neutron meter inMay, July, August, and September of 2016.

After the drilling samples were taken, they were packed into a numbered aluminum box and weighed immediately. Then, the samples were returned to the laboratory and dried in an oven until the constant weight was determined again. The soil mass moisture content was calculated by the calibration curve.

Data statistics and analysis were performed with Excel and SPSS, and plotted with Origin software. Inter-treatment comparisons were compared by one-way ANOVA.

**Treatment** 

**Average Value** 

**Standard Deviation** 

#### **3. Results**

#### *3.1. Effect of Irrigation Period on the Spatial Distribution of Soil Water Content*

To compare the effects of different irrigation regimes on soil water content, the W1F1, W1F2, W3F1, W3F1, W3F2, and W3F3 treatments in July 2016 (10 d after irrigation) have been selected to analyze the spatial distribution characteristics of soil water content (Table 1).

**Table 1.** Experimental treatments.


Under the same irrigation regime, the soil moisture availability shows a clear difference due to the different root growth distributions of *Haloxylon ammodendron* and *Calligonum mongolicum* (Figure 4). In the vertical depth, the soil water content of *Calligonum mongolicum* at 0–100 cm is greater than 100–200 cm in all treatments. The soil water content of *Haloxylon ammodendron* at 0–100 cm in the W3F1 and W3F2 treatments is less than 100–200 cm, and vice versa under other treatments. At the horizontal distance, the soil water content of *Calligonum mongolicum* at 30 cm under the W3F1 and W3F2 treatments is less than 70 cm, and vice versa under other treatments. The soil water content of *Haloxylon ammodendron* at 30 cm under the W3F1, W3F2, and W1F1 treatments is less than 70 cm, and vice versa under other treatments. *Water* **2022**, *14*, x FOR PEER REVIEW 6 of 13

**Figure 4.** Spatial variation in soil moisture content under different irrigation periods. The soil moisture content in the low (W1), medium (W2), and high (W3) irrigation amounts. F means irrigation frequency with F1 = 10 d, F2 = 20 d, F3 = 40 d. C and H represent *C. mongolicum* and *H. ammodendron*, respectively. **Figure 4.** Spatial variation in soil moisture content under different irrigation periods. The soil moisture content in the low (W1), medium (W2), and high (W3) irrigation amounts. F means irrigation frequency with F1 = 10 d, F2 = 20 d, F3 = 40 d. C and H represent *C. mongolicum* and *H. ammodendron*, respectively.

One-way ANOVA is used to compare the differences in soil water content between

mal distribution according to the Kolmogorov–Smirnov test. In addition, the F-test of equality of variances with the Least Significance Difference (LSD) test shows the homogeneity of variance. Therefore, one-way ANOVA is suitable for statistical analysis. One-way ANOVA displays the differences in soil water content between treatments at the 0.05 significance level. The mean soil water content of *Calligonum mongolicum* at 0–200 cm is W3F1 > W1F1 > W3F3 > W1F2 > W3F2, while the mean soil water content of *Haloxylon ammodendron* at 0–200 cm is W3F1 > W3F3 > W3F2 > W1F1 > W1F2. Under a 35 mm single irrigation amount, the soil moisture content varies significantly between F1 and F2, F1 and F3, but not significantly between F2 and F3 under a 17.5 mm (strain grade) single irrigation amount; the soil water content is F1 > F2, but this is not significant during the study periods. For the same total water irrigation amount, the difference in soil water content between W3F2 and W1F1, W3F3 and W1F2 are insignificant in the *Calligonum mongolicum*, but both are significant in the *Haloxylon ammodendron*. Under a 70 mm total water irrigation amount, the difference in soil water content of *Haloxylon ammodendron* is mainly focused at the 100–200 cm layer. Under 35 mm water irrigation, the soil moisture of the two

species is significantly different between 0–100 cm and 100–200 cm.

**Soil Water Content % 0–100 cm 100–200 cm 0–200 cm** 

**Value Standard Deviation Average Value Standard Deviation** 

**Table 2.** Soil moisture content between different treatments (*p* < 0.05).

**Average** 

One-way ANOVA is used to compare the differences in soil water content between the treatments shown in Table 2. The soil samples taken at different soil depths obey normal distribution according to the Kolmogorov–Smirnov test. In addition, the F-test of equality of variances with the Least Significance Difference (LSD) test shows the homogeneity of variance. Therefore, one-way ANOVA is suitable for statistical analysis. One-way ANOVA displays the differences in soil water content between treatments at the 0.05 significance level. The mean soil water content of *Calligonum mongolicum* at 0–200 cm is W3F1 > W1F1 > W3F3 > W1F2 > W3F2, while the mean soil water content of *Haloxylon ammodendron* at 0–200 cm is W3F1 > W3F3 > W3F2 > W1F1 > W1F2. Under a 35 mm single irrigation amount, the soil moisture content varies significantly between F1 and F2, F1 and F3, but not significantly between F2 and F3 under a 17.5 mm (strain grade) single irrigation amount; the soil water content is F1 > F2, but this is not significant during the study periods. For the same total water irrigation amount, the difference in soil water content between W3F2 and W1F1, W3F3 and W1F2 are insignificant in the *Calligonum mongolicum*, but both are significant in the *Haloxylon ammodendron*. Under a 70 mm total water irrigation amount, the difference in soil water content of *Haloxylon ammodendron* is mainly focused at the 100–200 cm layer. Under 35 mm water irrigation, the soil moisture of the two species is significantly different between 0–100 cm and 100–200 cm.


**Table 2.** Soil moisture content between different treatments (*p* < 0.05).

Note: A lowercase letter indicates a significant difference (*p* = 0.05) after irrigation at the different intervals. The markers (a, b, c) with different letters differ significantly (*p* < 0.05) according to Least Significance Difference test.

*3.2. Temporal Variation in Soil Moisture Profile during Irrigation Period under Different Irrigation Regimes*

According to the W1F1, W1F2, W3F1, W3F2, and W3F3 treatments, the dynamic changes in soil water content during the irrigation period are analyzed. The soil surface water content reaches the maximum on the first day after irrigation. The water moisture at the 0–60 cm layer gradually decreases with the temporal variation after irrigation. Subsequently, the variation range of soil moisture decreases with the increase in soil depth. The soil moisture above 60 cm is greatly affected by meteorological factors, while the soil moisture below 60 cm is mainly influenced by water redistribution and water absorption by roots (Figure 5).

Note: A lowercase letter indicates a significant difference (*p* = 0.05) after irrigation at the different intervals. The markers (a, b, c) with different letters differ significantly (*p* < 0.05) according to Least

*3.2. Temporal Variation in Soil Moisture Profile during Irrigation Period under Different Irriga-*

According to the W1F1, W1F2, W3F1, W3F2, and W3F3 treatments, the dynamic changes in soil water content during the irrigation period are analyzed. The soil surface water content reaches the maximum on the first day after irrigation. The water moisture at the 0–60 cm layer gradually decreases with the temporal variation after irrigation. Subsequently, the variation range of soil moisture decreases with the increase in soil depth. The soil moisture above 60 cm is greatly affected by meteorological factors, while the soil moisture below 60 cm is mainly influenced by water redistribution and water absorption

W3F1-C 5.73 a 1.08 2.93 a 1.29 4.33 a 1.84 W3F2-C 2.72 b 1.10 1.15 b 0.46 1.93 b 1.15 W3F3-C 3.37 b 1.26 1.70 b 0.45 2.53 b 1.26 W1F1-C 3.92 b 1.41 1.93 b 1.01 2.93 b 1.58 W1F2-C 2.95 b 2.45 1.57 b 0.81 2.26 b 1.95 W3F1-H 3.97 ab 1.03 4.90 a 0.70 4.44 a 0.99 W3F2-H 2.69 bc 0.91 3.43 b 1.17 3.06 b 1.11 W3F3-H 4.43 a 1.42 2.43 c 0.34 3.43 b 1.44 W1F1-H 3.09 ab 1.32 0.86 d 0.71 1.97 c 1.54 W1F2-H 1.76 c 1.01 1.20 d 0.59 1.48 c 0.87

Significance Difference test.

*tion Regimes* 

by roots (Figure 5).

**Figure 5.** Temporal variation in soil moisture content in different treatments. The soil moisture content in the low (W1), medium (W2), and high (W3) irrigation amounts. F means irrigation frequency with F1 = 10 d, F2 = 20 d, F3 = 40 d. C and H represent *C. mongolicum* and *H. ammodendron*, respectively. **Figure 5.** Temporal variation in soil moisture content in different treatments. The soil moisture content in the low (W1), medium (W2), and high (W3) irrigation amounts. F means irrigation frequency with F1 = 10 d, F2 = 20 d, F3 = 40 d. C and H represent *C. mongolicum* and *H. ammodendron*, respectively.

Therefore, the soil is divided into two layers for analysis: 0–60 cm (shallow layer) and 60–200 cm (deep layer). The soil water content of 0–60 cm is significantly greater than that of 60–200 cm and decreases rapidly after irrigation, with a decrease rate greater than the soil water content for 60–200 cm. The F1 and F2 irrigation periods decrease rapidly at 1–9 d after irrigation, and the decline rate of deep soil water content is slow or unchanged on the ninth day. The water content in the shallow and deep soil is similar. The treated soil water content under the F2 irrigation period is in a slow decline period at 9–19 d. The soil moisture in the shallow layers under the F3 treatment decreases rapidly during 1–4 d, decreases slowly during 4–9 d, and remains relatively stable in both the shallow and deep layers during 9–39 d (Figure 6).

#### *3.3. Response of Soil Moisture to Irrigation Regime in Different Months*

One-way ANOVA is used to compare the soil moisture content at 0–300 cm under different treatments in the same month from the 0–300 cm layer on 22 May, 12 July, 20 August, and 20 September 2016. As shown in Figure 7, the soil water content of *Calligonum mongolicum* and *Haloxylon ammodendron* at 0–300 cm decreases from May to July, and increases from July to September. The lowest water content is observed in July, while the highest is in September.

layers during 9–39 d (Figure 6).

Therefore, the soil is divided into two layers for analysis: 0–60 cm (shallow layer) and 60–200 cm (deep layer). The soil water content of 0–60 cm is significantly greater than that of 60–200 cm and decreases rapidly after irrigation, with a decrease rate greater than the soil water content for 60–200 cm. The F1 and F2 irrigation periods decrease rapidly at 1–9 d after irrigation, and the decline rate of deep soil water content is slow or unchanged on the ninth day. The water content in the shallow and deep soil is similar. The treated soil water content under the F2 irrigation period is in a slow decline period at 9–19 d. The soil moisture in the shallow layers under the F3 treatment decreases rapidly during 1–4 d, decreases slowly during 4–9 d, and remains relatively stable in both the shallow and deep

**Figure 6.** Variation in soil moisture content with the number of days after irrigation. **Figure 6.** Variation in soil moisture content with the number of days after irrigation.

**Figure 7.** Monthly changes in soil moisture content at 0–300 cm under different treatments. (Lowercase letters indicate differences between treatments in the same month. The markers (a, b, c) with different letters differ significantly (*p* < 0.05) according to Least Significance Difference test.) **Figure 7.** Monthly changes in soil moisture content at 0–300 cm under different treatments. (Lowercase letters indicate differences between treatments in the same month. The markers (a, b, c) with different letters differ significantly (*p* < 0.05) according to Least Significance Difference test).

ments, but decreases under other treatments. The difference in water storage at 0–300 cm is F1 > F3 > F2, which are all positive values. Under W2 treatment, the water storage is F2

**Table 3.** Differences in soil water storage between May and September at 0–300 cm soil layer under each treatment (C and H refer to *Calligonum mongolicum* and *Haloxylon ammodendron*, respectively).

**Treatment Difference in Soil Water Storage in September and May (mm)** 

W1F1-C 51.00 −8.03 −6.62 36.35 W1F2-C 26.71 −3.88 −3.75 19.08 W1F3-C 34.13 −2.22 −5.70 26.21 W2F1-C 7.47 −10.02 −9.44 −11.99 W2F2-C 26.59 −4.44 −13.93 8.22 W2F3-C −10.95 −3.58 −8.28 −22.81 W3F1-C 34.62 −2.94 4.70 36.37

**0–100 cm 100–200 cm 200–300 cm 0–300 cm** 

> F1 > F3, and the water storage is reduced under the W2F1 and W2F3 treatments.

Table 3 shows the differences in soil water storage between May and September at the 0–300 cm soil layer under each treatment. The water storage of *Calligonum mongolicum* at 0–100 cm increases, except for the W2F3 treatment. The water storage at 100–200 cm increases by 2.26 mm only in the W3F3 treatment, while values for other treatments de-

Table 3 shows the differences in soil water storage between May and September at the 0–300 cm soil layer under each treatment. The water storage of *Calligonum mongolicum* at 0–100 cm increases, except for the W2F3 treatment. The water storage at 100–200 cm increases by 2.26 mm only in the W3F3 treatment, while values for other treatments decrease. The water storage at 200–300 cm increases slightly under the W3F1 and W3F2 treatments, but decreases under other treatments. The difference in water storage at 0–300 cm is F1 > F3 > F2, which are all positive values. Under W2 treatment, the water storage is F2 > F1 > F3, and the water storage is reduced under the W2F1 and W2F3 treatments.


**Table 3.** Differences in soil water storage between May and September at 0–300 cm soil layer under each treatment (C and H refer to *Calligonum mongolicum* and *Haloxylon ammodendron*, respectively).

The water storage of *Haloxylon ammodendron* at 0–100 cm decreases under the W1F1, W1F2, W1F3, and W2F2 treatments, while values for the other treatments increase. At 100–200 cm, the water storage volume increases under the W1F1, W1F2, W2F2, and W3F3 treatments. The other treatments' values are reduced. The water storage from 200 to 300 cm is reduced under the W1F3, W2F1, and W2F2 treatments, while values for all other treatments increase. The water storage difference of *Haloxylon ammodendron* at 0–300 cm is F3 > F1 > F2 under W2 > W3 irrigation amount, which are positive values. Under the W1 irrigation amount, the water storage difference is F2 > F1 > F3, and the treatment for W2F3 has a higher value than that for W2F1.

#### **4. Discussion**

In arid and semi-arid regions, the water resources directly affect the distribution and growth of plants. In the Taklimakan Desert, the climate is extremely dry, and the annual precipitation (36.6 mm) is far from sufficient to meet the evapotranspiration demand (3638.6 mm). The groundwater depth is more than 10 m, the replenishment effect of groundwater on soil water is negligible, and the main source of soil water is irrigation water [1]. The plants are facing the danger of long-term water shortage. After the Taklimakan Desert Highway shelterbelt was built, the tree species in the shelterbelt were mainly salt-tolerant. The water for plant growth came from underground high-salinity water drip irrigation. At present, although the current drip irrigation system can basically satisfy the growth of tree species in the shelterbelt, the utilization efficiency of plants for irrigation water is low [2]. The study of a reasonable saline irrigation becomes the basis of shelterbelt management and its sustainable existence.

In this study, the water supply rate is specific, and soil water infiltration is determined by soil water infiltration capacity, which is associated with soil wetness and porosity [20,21]. The surface soil moisture content is low from June to early August, resulting in slow transverse water transfer and little change due to intense surface evaporation. The growth of plant roots increases the non-capillary pores, improves the soil water conductivity, and lays a foundation for efficient soil water transport and storage [22]. Due to the downward growth of plant roots and the formation of soil macropores, the irrigation water can quickly reach the deep soil, causing periodic changes in soil water content. The lower limit of soil water evaporation is 40–60 cm, and excessive surface water content increases water evaporation, which is not conducive to water storage. The deeper the soil layer, the lower the soil moisture. Affected by atmospheric evaporation and water absorption by plant roots, the soil water storage variation coefficient is smaller. The results are consistent with the previous studies, such as Li et al. [1] and Zhang et al. [2].

After irrigation, the variation in soil water can be divided into a period of rapid water decline, slow water decline, and a relatively stable water level. Due to the loose soil in the sandy land, the shallow soil water can quickly infiltrate after irrigation. Subsequently, the strong evaporation effect leads to a dry sand layer forming on the soil surface, which significantly inhibits soil moisture evaporation [23]. The water absorption of plants mainly causes a decrease in soil moisture. When the soil moisture is relatively stable, soil moisture at the 0–200 cm layer cannot meet the needs of plant growth, and the water demand of plants mainly comes from the deep soil water supply, while soil moisture maintains a relatively stable state [2].

We noted the largest differences in soil water content between the two plant types in July, followed by August, May, and September. Soil water dissipation mainly includes soil evaporation and plant transpiration. In July, the temperature is the highest, and the water requirement for plant transpiration and soil evaporation is the largest, so the difference between treatments is the largest. In September, the temperature decreases, plant growth slows down, the water requirement decreases, and the difference between treatments is the least.

Under the irrigation regime with a 35 mm irrigation amount, the plants grow well, and part of the soil evaporation is reduced by shading. Irrigation can not only meet the needs of plant growth but also replenish soil water. Under the irrigation regime with a 17.5 mm irrigation amount, the plant growth is weakened due to drought stress in the early stage, although the irrigation amount is small. The leaf transpiration and root water absorption are reduced, and water dissipation is weakened. Moreover, insufficient soil water partially inhibits soil evaporation [24]. Therefore, this study highlights that an irrigation regime with a 35 mm irrigation amount is beneficial to soil water storage. More irrigation water will infiltrate into the deep soil, which will not be absorbed and utilized by plants, resulting in water waste. Increasing the single irrigation amount and prolonging the irrigation period can allow the more effective use of irrigation water. This study highlights that saline groundwater irrigation provides potential advantages for desert plants' survival under reasonable irrigation regimes.

The saving of water and improvement of water use efficiency are undoubtedly fundamental problems associated with such drought regions to avoid lowering the groundwater levels and to prevent ecological degradation. Although desert plants have strong resistance to water and saline stresses and different stress adaptation mechanisms at different growth stages [25], this study evaluated water dynamics and irrigation regimes only under preset irrigation combinations. Further work should conduct additional studies to examine whether or not an appropriate smaller amount of repeated irrigation will increase the water use efficiency of plants and reduce the ineffective evaporation of water. In addition, soil water infiltration and evaporation lead to salinity storage in the soil. The presence and accumulation of salinity affect plants' physiological ecology. The plants' adaptation to

saline water irrigation and their responses to the different irrigation regimes should be considered in a future study.

#### **5. Conclusions**

Based on a field test of precision irrigation control in the Taklimakan Desert Highway shelterbelt, this study examined the effects of irrigation regimes on the soil water dynamics of two typical woody halophyte species (*Haloxylon* and *Calligonum*), and quantified the irrigation intervals and periods. The effects of saline water irrigation regimes on the soil water dynamics of two typical woody halophyte species (i.e., *Calligonum mongolicum* and *Haloxylon ammodendron*) show that:


This study highlights that saline groundwater irrigation is advantageous for supporting desert plants' survival and preventing ecological degradation under reasonable irrigation regimes. Future work should focus on the plants' adaptation to saline water irrigation and their responses to the different irrigation regimes and water-saving irrigation measures in the desert shelterbelt construction.

**Author Contributions:** Conceptualization and methodology, J.L., Y.Z., J.Z. and J.X.; data analysis, J.L. and Q.H.; writing—original draft preparation, J.L.; writing—review and editing, Y.Z., J.Z. and J.X. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was supported by the National Natural Science Foundation of China (41977009, 41877541, 41471222, and 42071259), the original innovation project of the basic frontier scientific research program, Chinese Academy of Sciences (ZDBS-LY-DQC031), the Natural Science Foundation of Xinjiang Uygur Autonomous Region (2021D01E01), the Third Batch of Tianshan Talents Program of Xinjiang Uygur Autonomous Region (2021–2023), the Youth Innovation Promotion Association of the Chinese Academy of Sciences (2019430), and the State Key Laboratory of Desert and Oasis Ecology, Xinjiang Institute of Ecology and Geography, Chinese Academy of Sciences (G2018-02–08).

**Data Availability Statement:** The data are available from the corresponding author upon reasonable request.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


**Chuanyu Ma <sup>1</sup> , Luobin Tang <sup>1</sup> , Wenqian Chang <sup>1</sup> , Muhammad Tauseef Jaffar <sup>1</sup> , Jianguo Zhang 1,2,\*, Xiong Li <sup>1</sup> , Qing Chang <sup>2</sup> and Jinglong Fan <sup>2</sup>**


**Abstract:** To explore the impact of artificial shelterbelt construction with saline irrigation on the soil water characteristic curve (SWCC) of shifting sandy soil in extreme arid desert areas, three treatments including under the shelterbelt (US), bare land in the shelterbelt (BL) and shifting sandy land (CK) in the hinterland of the Taklimakan Desert were selected. The age of the shelterbelt is 16, and the vegetation cover is mainly *Calligonum mongolicum*. The soils from different depths of 0–30 cm were taken keeping in view the objective of the study. The SWCCs were determined by the centrifugal method and fitting was performed using various models such as the Gardner (G) model, Brooks–Corey (BC) model and Van Genuchten (VG) model. Then, the most suitable SWCC model was selected. The results showed that electrical conductivity (EC) and organic matter content of BL and US decreased with the increasing soil depth, while the EC and organic matter content of CK increased with the soil depth. The changes in soil bulk density, EC and organic matter of 0–5 cm soil were mostly significant (*p* < 0.05) for different treatments, and the differences in SWCCs were also significant among different treatments. Moreover, the construction of an artificial shelterbelt improved soil water-holding capacity and had the most significant impacts on the surface soil. The increase in soil water-holding capacity decreased with increasing soil depth, and the available soil water existed in the form of readily available water. The BC model and VG model were found to be better than the G model in fitting results, and the BC model had the best fitting result on CK, while the VG Model had the best fitting result on BL with higher organic matter and salt contents. Comparing the fitting results of the three models, we concluded that although the fitting accuracy of the VG model tended to decrease with increasing organic matter and salinity, the VG model had the highest fitting accuracy when comparing with BC and G models for the BL treatment with high organic matter and salinity. Therefore, the influence of organic matter and salinity should be considered when establishing soil water transfer function.

**Keywords:** aeolian sandy soil; physiochemical properties; soil water-holding capacity; model fitting; Taklimakan Desert

#### **1. Introduction**

Soil moisture is an important factor affecting plant growth and is a major driving force for the sustainability of many terrestrial ecosystems. Moisture changes have significant impacts on vegetation and soil properties [1–3]. Especially in arid and semi-arid regions, soil moisture availability is one of main factors limiting the type and quantity of vegetation, and the water deficiency can lead to severe degradation of vegetation and reduce vegetation cover [4,5]. Therefore, it is important to understand the soil moisture change pattern for the maintenance of vegetation in arid desert areas [6].

**Citation:** Ma, C.; Tang, L.; Chang, W.; Jaffar, M.T.; Zhang, J.; Li, X.; Chang, Q.; Fan, J. Effect of Shelterbelt Construction on Soil Water Characteristic Curves in an Extreme Arid Shifting Desert. *Water* **2022**, *14*, 1803. https://doi.org/10.3390/ w14111803

Academic Editor: Ognjen Bonacci

Received: 27 April 2022 Accepted: 31 May 2022 Published: 2 June 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

The Taklamakan Desert is located in the hinterland of the Eurasian continent in Chinese southern Xinjiang, and is the largest desert in China and the second largest shifting desert in the world. Southern Xinjiang is one of the poorest regions with an extreme drought climate in China. To accelerate the development of regional social economy, the Taklimakan Desert Highway (TDH) was completed in 1995. TDH is 522 km long running across the Taklimakan Desert from north to south and is the longest highway crossing a moving desert in the world. However, serious sand disasters are a great threat to TDH. Thus, in 2003, a shelter–forest belt (the Taklimakan Desert Highway Shelterbelt (TDHS)) having a length of 436 km was constructed on both sides along the highway, which was dominated by shrubs to protect the highway from shifting sand [7]. Low-quality saline groundwater is used for irrigation in order to ensure the survival of shelterbelt plants [8–11]. The regional ecological environment on both sides of the highway has been greatly improved by the fixation of the shifting sand dunes with the artificial shelterbelt [12]. However, saline groundwater irrigation may aggravate soil salinization in future and is harmful for shelterbelt plants [13,14]. Therefore, it is necessary to study the water and salt transport of shelterbelt soils which is helpful for the sustainable utilization and management of TDHS [15].

SWCC, as a component of the unsaturated soil mechanics framework, provides the information needed to characterize the properties of unsaturated soils [16], is an interpretation of the basic constitutive relationships of unsaturated soil phenomena [17], and is an important tool for studying the properties of unsaturated soils [18]. SWCC is one of basic hydraulic properties that simulate water and solute transport under unsaturated conditions [19]; is universally used in agriculture [20], soil physics [21], soil chemistry [22], mineralogy research [23], geotechnical engineering [24]; and is widely used in the soil– plant–atmosphere continuum (SPC) [18,25] and other fields. Due to its importance in soil hydrodynamics and solute transport modeling, many SWCC models, both numerical and theoretical, have been developed [21]. A good SWCC model should have simple and clear parameters and be easy to use. It can satisfy the three characteristics of accuracy, universality and simplicity as much as possible at the same time [26]. The VG model [27], BC model [28] and G model [29] have relatively few parameters, and can accurately describe the SWCC of various soil textures [30–38]. SWCC is influenced by soil texture, bulk density, organic matter, salinity, temperature, etc. [19,30–33,39–42]. Therefore, the fitting results of the models are often different in various study areas. For example, some scholars [43,44] pointed out that the G model can accurately fit SWCC; however, others [45] came to the opposite conclusion, pointing out that the G model cannot accurately fit SWCC. Matlan et al. [26] compared four models and found that the BC model has the most accurate description of the SWCC of sandy soils. Li et al. [46] also proved this, pointing out that the BC model is more suitable than the VG model on soils with high sand content and low clay content. For the SWCC of aeolian sand covered by biocrust, the fitting effect of the VG model is better than the BC model [47].

The literature on SWCC in arid and semi-arid regions is relatively limited as compared to farmland and forest ecosystems. The influence of artificial shelterbelts and long-term saline irrigation on SWCC of sandy soils and the applicability of SWCC models is still unclear. This hinders our understanding towards the soil water-holding capacity and water availability of shelter forests. In this study, we firstly assumed that artificial shelterbelts and long-term saline irrigation have impacts on the water-holding capacity and water availability of different soil layers. In addition, we assumed that the VG model, BC model, and G model have different accuracies in fitting SWCC. Therefore, our study collected soils from 0 to 30 cm layers under the shelterbelt (US), bare land in the shelterbelt (BL) and shifting sandy land (CK) in the hinterland of the Taklimakan Desert, and their SWCCs were determined by the centrifugal method. Combined with its bulk density, organic matter, salinity and other properties, SWCC, pore distribution and soil moisture were analyzed, and the VG model, BC model and G model were used and compared to fit the SWCCs. We aim to reveal the impact of artificial shelterbelt construction on SWCC on shifting sandy

soil in extreme arid deserts under saline drip irrigation, so as to provide a basis for desert shelterbelt construction and sustainable management.

#### **2. Materials and Methods**

#### *2.1. Study Area*

The sampling area is located at the Taklamakan Desert Research Station of the Chinese Academy of Sciences in the hinterland of Taklamakan Desert (39◦010 N, 83◦360 E, 1100 m.a.s.l.), Xinjiang, China (Figure 1). The study area belongs to a warm temperate arid climate. The annual average temperature is 12.4 ◦C. December is the coldest month with an average monthly temperature of −8.1 ◦C and July is the hottest month with an average monthly temperature of 28.2 ◦C. Annual precipitation is 24.6 mm, average relative humidity is 29.4% and annual potential evaporation is up to 3638.6 mm. Annual average wind speed is 2.5 m·s −1 , and the maximum instantaneous wind speed is up to 20 m·s −1 . The soil type is Xeric Quartzipsamments [48] (Soil Survey Staff 2014), derived from shifting eolian sand, and the basic properties are shown in Table 1.

**Figure 1.** Location of the study area.

**Table 1.** Basic physiochemical properties of the collected soil samples.


Note: CK: shifting sandy land; BL: bare land without vegetation cover in the shelterbelt; US: under the shelterbelt. Different lowercase letters represent the significant differences between different treatments (*p* < 0.05). (Soil texture was classified according to USDA standards based on actual measured soil mechanical composition).

#### *2.2. Sample Collection and Determination*

In August 2020, undisturbed soil samples were collected with a cutting ring of 100 cm<sup>3</sup> from the soil layers of 0–5 cm, 5–10 cm, 10–20 cm and 20–30 cm under US, BL and CK, and all samples had three replicates (the pictures of the sampling points are shown in Figure 2). The ring knife (with soil) samples were soaked in deionized water for 24 h until saturation, weighed, and then the SWCCs were determined using a high-speed centrifuge (CR 22G III model, Hitachi, Japan). Suction values were determined in lab using a centrifuge with different speeds and time settings. The selected suction values were 10.2 cm (310 r/min, 10 min), 30.6 cm (540 r/min, 12 min), 51 cm (690 r/min, 17 min), 71.4 cm (820 r/min, 21 min), 102 cm (980 r/min, 26 min), 204 cm (1200 r/min, 28 min), 612 cm (2190 r/min, 49 min), 1020 cm (3100 r/min, 57 min), 4080 cm (6200 r/min, 77 min) and 8160 cm (8770 r/min, 87 min), respectively. The weight of the ring knife sample was weighed after the completion of centrifugation at each suction value. After centrifugation, the ring knife was dried in an oven (105 ◦C) and then weighed, and the water content at the corresponding suction value was obtained from the difference in weights. The moisture under different suction values was plotted according to the SWCC. Meanwhile, soil samples of corresponding layers were also collected with a soil drill and brought to the laboratory and air dried. Part of the samples was passed through a 2 mm sieve to determine soil pH, EC (with an EC500 pH/Conductivity Meter (ExStik, Boston, MA, USA)), and soil mechanical composition was determined using the hydrometer method and soil texture was classified according to USDA standards. Soil particles were divided into sand particles of 2–0.05 mm, silt particles of 0.05–0.002 mm and clay particles of ≤0.002 mm, and the left soils were used to measure soil organic matter content with potassium dichromate external heating method after passing through a 0.25 mm sieve. The physical and chemical properties of the soil samples are shown in Table 1.

**Figure 2.** Photos of sampling points.

#### *2.3. Calculations of Soil Equivalent Pore Sizes and Moisture Constants*

The equivalent pores are divided into six levels based on their diameters, including narrow micropores (≤0.3 µm), micropores (0.3–5 µm), fine pores (5–30 µm), medium pores (30–75 µm), macropores (75–100 µm) and interstices (≥100 µm) [47,49].

Based on the measured and fitted SWCC parameters, saturated water *θ<sup>s</sup>* , field capacity *θf* , wilting coefficient *θ<sup>r</sup>* , available water content *θ<sup>a</sup>* and readily available water content *θra* were obtained. While *θ<sup>s</sup>* is the soil water content when water suction is zero, *θ<sup>f</sup>* is the soil water content of when *pF* = 1.8 (*pF* is expressed as the logarithm of the centimeter height of the water column of the soil water potential), and *θ<sup>r</sup>* is the soil water content when *pF* = 4.2, *pF* = 3.8 represents temporary wilting coefficient. The water content in the range of *pF* 1.8–4.2 is *θa*, and the water content in the range of *pF* 1.8–3.8 is *θra* [50].

*2.4. SWCC Modeling*

VG model, G model and BC model were used to describe the SWCC.

VG model:

$$\theta\_{(h)} = \theta\_r + \frac{\theta\_s - \theta\_r}{(1 + |\alpha h|^{\
\eta})^m} \tag{1}$$

where *θ(h)* is the volume water content, cm<sup>3</sup> cm−<sup>3</sup> ; *θ<sup>s</sup>* is the saturated water content, cm<sup>3</sup> cm−<sup>3</sup> ; *θ<sup>r</sup>* is the residual water content, cm<sup>3</sup> cm−<sup>3</sup> ; *h* is the matric suction, cm; *α* is a scale parameter that is related to the inverse of the air entry suction, cm−<sup>1</sup> ; *m*, *n* are the fitting parameters; *m* = 1 − 1/*n* (*n* > 1).

G model:

$$
\theta = ah^{-b} \tag{2}
$$

where *θ* is the volume water content, cm<sup>3</sup> cm−<sup>3</sup> ; *h* is the matric suction, cm; *a* and *b* are fitting parameters.

BC model:

$$S\_{\varepsilon} = \frac{\theta - \theta\_r}{\theta\_s - \theta\_r} = \begin{cases} \left( ah \right)^{-\lambda} & ah > 1 \\ 1 & ah \le 1 \end{cases} \tag{3}$$

where *S<sup>e</sup>* is the effective saturation; *λ* is the curve shape parameter; the physical meanings of other parameters are the same as above.

#### *2.5. SWCC Model Fitting Accuray Assesments*

The coefficient of determination (*R* 2 ), root mean square error (*RSME*) and relative error (*RE*) were used to quantitatively evaluate the fitting effect of the models. Grey correlation analysis allows ranking the importance between different influencing factors [51,52]. The soil physiochemical properties and model parameters were quantitatively analyzed by grey correlational method, and the correlational degree was sorted. The calculation formulas are as follows:

$$\mathcal{R}^2 = \frac{\sum\_{i=1}^{N} \left(\theta\_i - \overline{\theta\_i}\right)^2}{\sum\_{i=1}^{N} \left(\beta\_i - \overline{\theta\_i}\right)^2} \tag{4}$$

$$RMSE = \sqrt{\frac{\sum\_{i=1}^{N} \left(\theta\_i - \beta\_i\right)^2}{N}} \tag{5}$$

$$RE = \frac{|\theta\_{\mathbf{i}} - \theta\_{\mathbf{i}}|}{\theta\_{\mathbf{i}}} \times 100\tag{6}$$

$$\gamma\_{0i} = \frac{1}{n} \sum\_{k=1}^{n} \frac{m + \rho M}{\Delta\_{\rm i}(k) + \rho M} \tag{7}$$

In Equations (4)–(6), *N* is the total number of samples of matric suction; *θ<sup>i</sup>* represents the measured value of soil moisture corresponding to the *i*th pressure value; *θ<sup>i</sup>* represents the average value of measured soil moisture; *β<sup>i</sup>* is the fitted value of soil moisture corresponding to the *i*th pressure value. In Equation (7), *γ*0*<sup>i</sup>* is the correlation degree, ∆*<sup>i</sup> (k)* is the difference sequence, *M* is the maximum difference sequence, *m* is the minimum difference sequence, and *ρ* is the resolution coefficient, which is generally 0.5 in the models [51,52].

RETC software was used to solve and fit the parameters, Excel2010 and Origin2018 were used for data processing and mapping and SPSS18.0 was used for one-way ANOVA and multiple comparisons (*a* = 0.05, LSD).

#### **3. Results**

#### *3.1. Effects of Artificial Shelterbelt Construction on Soil Physiochemical Properties*

Soil physiochemical properties were improved after artificial shelterbelt construction. From Table 1, according to USDA system, it is clearly stated that all soil layers of CK are

sandy soil, all soil layers of BL treatment are loamy sandy soil, 0–5 cm soil of US treatment was sandy loam soil and 5–30 cm soil was loamy sandy soil. There were significant differences in bulk density and EC among CK, BL and US in the same soil layer (*p* < 0.05). Except for CK, the soil bulk density under BL and US treatments increased with the increasing soil depth, and the soil bulk density of 0–5 cm under BL and US treatments was significantly lower than CK; the bulk density of the 5–30 cm soil layer was the highest in BL, except for 10–20 cm; the bulk density of US treatment was lower than that of CK and BL. The bulk density of the 0–5 cm soil layer in US treatment was the smallest (1.25 g/cm<sup>3</sup> ), and the highest bulk density at all was 1.6 g/cm<sup>3</sup> . With the increase in soil depth, the EC of CK gradually increased, while EC of BL and US treatments decreased with the increasing depth. The difference in EC between BL, US and CK in the 0–5 cm soil layer was the largest, which was 18,764 µS/cm and 9724 µS/cm, respectively. The EC of CK at 10–30 cm was higher than that of BL and US.

#### *3.2. Screening of Soil Water Characteristic Curve Models*

As listed in Table 2, *R* <sup>2</sup> of the fitting values of VG, BC, and G models was ranged between 0.884 and 0.998, which showed the correlations of three models with the measured data were high. As shown in Figure 3, the fitted results of the G model were all higher than the measured points, and the fitted values of the VG model and BC model for BL in the 5–20 cm soil layer were smaller than the measured values.

**Figure 3.** Fitting curves of soil water characteristics of VG model, BC model and G model. (**a**–**l**) is the serial number of the figure, (**a**–**d**) is the model fitting result of CK processing 0–5 cm, 5–10 cm, 10–20 cm, 20–30 cm soil layer; (**e**–**h**) is BL processing 0–5 cm, 5–10 cm, 10–20 cm and the model fitting results of the 20–30 cm soil layer; (**i**–**l**) is the model fitting result of the US treatment of 0–5 cm, 5–10 cm, 10–20 cm and 20–30 cm soil layers.

The fitting effect of the VG model and BC model for different soil layers of various treatments were always better than that of the G model. The fitting effect of the VG model and BC model was similar for the 10–30 cm soil layers of each treatment. The fitting effect of the BC model for the 0–10 cm soil layer of CK was better than the VG model, while the results of BL were opposite. For US, the BC model had a good fitting effect for the 0–5 cm soil layer, but the fitting effect was contradictory for the 5–10 cm soil layer. For all treatments, the fitting errors of three models increased with the increasing water potential and tended to be stable (Figure 4). The relative errors of the VG model and BC model for different treatments and soil layers were always lower as compared with relative error of the G model.

**Figure 4.** Relative errors of SWCC models under CK, BL and US treatments. (The smaller the circle, the smaller the relative error).

**Table 2.** *R* <sup>2</sup> and *RMSE* of different fitting models.


In addition, *R* <sup>2</sup> and *RMSE* of the VG and BC models showed that the simulating results of CK were better than those of BL and US (Table 2). The results indicated that the increase in salinity and organic matter content may affect the fitness of each model. Table 3 lists the parameters of the three models in SWCC modeling. From grey correlation calculation (Table 4), it could be concluded that VG model parameters *a* and *n* have a higher degree of correlation with the soil physiochemical parameters (BC and G models were not included, because the VG model had the best fitting effect). The grey relational degrees between pH, EC, organic matter, bulk density, sand content, silt content, clay content and model parameters (*a* and *n*) were all larger than 0.6, and the order of correlation showed as: bulk density > sand content > pH > clay content > organic matter > silt content > EC.


**Table 3.** Parameters in the modeling of soil water characteristic curves.

Note: *θ<sup>s</sup>* : saturated water content; *θ<sup>r</sup>* : wilting coefficient; α is a scale parameter that is related to the inverse of the air entry suction; *n* is the fitting parameter; *λ* is the curve shape parameter; *a* and *b* are fitting parameters; *a* × *b* is the specific water capacity when the soil water suction is 1 Bar.

**Table 4.** Correlation degree analysis between VG model parameters and soil basic physiochemical parameters.


*3.3. Effects of Artificial Shelterbelt Construction on Soil Water Retention Performance*

Soil porosity was significantly changed after artificial shelterbelt construction. As listed in Figure 5, few micropores were found under all treatments. Compared with CK, the contents of interstices and macropores in the 0–10 cm soil layer under BL and US were much lower, and the contents of fine pores and micropores were higher. The medium pores of the BL increased the most, with 0–5 cm increased by 3.74% and 5–10 cm increased by 4.94%. Fine pores of US increased the most: 0–5 cm increased by 4.54% and 5–10 cm increased by 2.32%.

**Figure 5.** Distribution of soil equivalent pores in different soil layers under different treatments. (**a**–**d**) for 0–5 cm, 5–10 cm, 10–20 cm and 20–30 cm soil layer, respectively. CK: shifting sandy land; BL: bare land without vegetation cover in the shelterbelt; US: under the shelterbelt.

Soil water parameters of each treatment are shown in Table 5. In the 0–5 cm soil layer, *θs* , *θ<sup>f</sup>* , *θ<sup>r</sup>* , *θ<sup>a</sup>* and *θra* under BL and US were higher than those under CK. Compared with CK, *θ<sup>s</sup>* of BL and US increased by 4.42% and 12.67%, *θ<sup>f</sup>* increased by 68.9% and 70.41%, *θ<sup>r</sup>* increased by 32.84% and 69.47%, *θ<sup>a</sup>* increased by 87.84% and 70.97% and *θra* increased by 87.73% and 70.43%, respectively. In the 5–10 cm soil layer, *θ<sup>s</sup>* of CK was the largest, but *θ<sup>f</sup>* , *θr* , *θ<sup>a</sup>* and *θra* of US and BL were higher than those of CK. In the 10–20 cm soil layer, higher *θf* , *θ<sup>r</sup>* , *θ<sup>a</sup>* and *θra* were observed under US than in BL and CK. In the 20–30 cm soil layer, higher *θ<sup>s</sup>* , *θ<sup>f</sup>* , *θ<sup>r</sup>* , *θ<sup>a</sup>* and *θra* were observed under BL than in US and CK. Meanwhile, the contents of available water in each treatment were almost equal to that of readily available water content.


**Table 5.** Soil water parameters calculated from the VG models.

Note: *θ<sup>s</sup>* : saturated water content; *θ<sup>f</sup>* : field capacity; *θ<sup>r</sup>* : wilting coefficient; *θa*: available water content; *θra*: readily available water content. (BC and G model were not included because the VG model had the best fitting effect).

As shown in Figure 6, the shapes of SWCC of each soil layer were similar, but the variations in soil moisture under per unit suction were clearly different. Under 0–10 cm and 1000–10,000 cm suctions, soil water was lost slowly with the increase in suction. However, under the suction of 10–1000 cm, the curve trended to be steep, and the soil water decreased rapidly with the increase in suction. Figure 6 clearly showed the differences in the course of the water retention curves of the three sampled areas. Soil water-holding capacity of 0–5 cm soil layer was highest under US, followed by BL and CK (Figure 6a). The difference

in water-holding capacity among treatments in the remaining soil layers gradually became smaller. At the suction value corresponding to the field water-holding capacity (*p<sup>F</sup>* = 1.8, i.e., 63 cm water column), the moisture of the 5–10 cm soil layer was as follows: BL ≥ US≥ CK; 10–20 cm soil layer was as follows: US ≥ CK ≥ BL; and 20–30 cm soil layer was as follows: BL ≥ CK ≥ US. When reaching the suction value corresponding to the temporary wilting coefficient (*p<sup>F</sup>* = 3.8, i.e., 6309 cm water column), the water content of the 5–10 cm soil layer behaved as follows: BL ≥ CK ≥ US; the water content of the 10–20 cm soil layer behaved as follows: US ≥ BL ≥ CK; and the water content of the 20–30 cm soil layer behaved as follows: BL ≥ CK ≥ US.

**Figure 6.** Water characteristic curves of each soil layer under different treatments. The fitted values are based on VG model; (**a**) 0–5 cm soil layer; (**b**) 5–10 cm soil layer; (**c**) 10–20 cm soil layer; (**d**) 20–30 cm soil layer.

#### **4. Discussion**

#### *4.1. Artificial Shelterbelt Construction Greatly Changed the Soil Physiochemical Properties*

SWCC is affected by various soil properties such as texture, bulk density, porosity, organic matter and salinity [42,53–55]. In our study, compared with shifting sandy land, the soil properties of BL and US changed obviously (Table 1), which caused the transformation of SWCC (Figure 6). Soil physiochemical properties of the 0–5 cm soil layer were changed most significantly. Shelterbelt construction under saline irrigation significantly decreased the soil bulk density, increased EC and organic matter content, and the soil texture changed from sandy soil to loamy sand and sandy loam (Table 1). After shelterbelt construction, the soil bulk density of 0–5 cm decreased significantly as compared with CK, which is primarily due to the continuous input of plant litters in the surface soil that lead to loosening the soil particles [56,57]. The bulk density of BL was the highest in the 5–30 cm soil layer, indicating that long-term saline irrigation would increase the soil bulk density [58,59]. The change in soil texture was mainly reflected by the increase in silt and clay particles. The main reason is that the shelterbelts reduced the wind speed, which promotes the precipitation of sand and the accumulation of dust fall [60,61]. The second reason is the accumulation of litter and the role of microorganisms. The volume and average radius of soil macropores increased with the increase in volumetric rock fragment content [62,63]. This directly led to the reduction in large pores and increase in small pores under BL and US (Figure 5a), which improved soil water-holding capacity and changed the SWCC. In addition to the influence of rock fragments, soil macropores are also controlled by biological factors [64].

Therefore, when the content of rock fragments is similar, the distribution of soil macropores will also be different (Figure 5b–d).

Under long-term saline irrigation, the salts were added into the soil and resulted in surface accumulation due to strong evaporation rate [13,65]. Therefore, the soil ECs of BL and US in the 0–10 cm soil layer were significantly higher than shifting sandy land, and the EC of BL and US decreased with the increasing soil depth. Vegetation cover resulted in less evaporation under US than BL, and more intense surface salt accumulation under BL. Therefore, soil EC of BL was higher than US in the 0–5 cm soil layer, and EC of BL was lower than US in the 0–10 cm soil layer. Wind erosion has an important impact on the cycle of soil organic carbon. The fine particles in the sand and dust adsorb soil organic carbon. Under the action of wind erosion, soil organic carbon is redistributed along with the movement of sand and dust [65,66]. Therefore, the increase in soil organic matter was mainly due to the accumulation of litters, and atmospheric dustfall also played a certain role in promoting it [67].

#### *4.2. Artificial Shelterbelt Construction Increased the Soil Water-Holding Capacity*

Soil water moves in pores and its transfer rate is directly determined by the size and distribution of pores. Soil bulk density is negatively correlated with soil porosity. Changes in soil primary particles such as sand, silt and clay also affect the distribution of pores. Organic matter and salinity have a direct impact on soil structure and adsorption. These soil physical properties may directly or indirectly affect soil water conductivity. Therefore, SWCC is affected by bulk density, texture, organic matter, porosity, aggregate stability, salinity and other properties [19,30–33,40–42].

Compared with CK for the 0–5 cm soil layer, vegetation coverage resulted in the decline in bulk density and increase in the salinity, organic matter and silt content under BL and US. On the one hand, salt contents in soil will occupy the pore space, and it will cause some soil particles to flocculate together, increase soil pores and enhance soil water-holding capacity [68]. Plant litters will reduce soil bulk density and increase soil organic matter content, soil saturated water content and water conductivity [69,70]. Therefore, the water contents of BL and US were higher than those of CK under the same suction, and with the increase in suction, the water contents of BL and US decreased, which means the SWCC integrally moved upwards and the trend slowed down. The contents of soil interstices and macropores in the 5–10 cm of CK were more than that in BL and US (Figure 5b). SWCC reflects the dehumidification process of soils, during which water is stored in pores and water retained in macropores is preferentially expelled as suction increases. Therefore, the change in water content per unit suction of CK was higher than that of BL and US. Therefore, the curve was the maximum under 0–10 cm suction, and then decreased rapidly to the minimum. The soil physiochemical properties of each treatment in the 10–20 cm soil layer had little difference, and their water-holding capacities were similar. For the 20–30 cm soil layer, the bulk density of BL was significantly higher than that of CK and US, and the water-holding curve was at its peak. Our results are consistent with Lipiec et al. [54]. Comparing the soil water changes over the entire suction range, the water-holding capacity of the surface 0–5 cm soil increased the most significantly, and the water-holding capacity of BL and US increased over the entire suction range compared to CK. The increase in water retention capacity of shelterbelt soils under long-term saline water irrigation is mainly due to the following reasons. Soil texture controls the physical, hydrological and chemical properties of the soil and has a strong influence on water flow paths, residence times and the magnitude and location of salt accumulation. Soil texture governs the water and solute transport. Under irrigation conditions, finer soils limit water infiltration, and coarse-grained soils retain significantly less water than fine-grained soils [71,72]. The accumulation of soil salts is dominated by sodium salts originated from the irrigation water. Excessive concentration of sodium ions in the soil solution will disperse and swell the soil structure, leading to the reduction and blockage of connected pores, thus reducing the permeability and hydraulic conductivity of the soil [73,74]. The salinity of pore water

also affects the development of the diffusion double layer (DDL) around soil particles, which controls the microstructural changes in soil particles during hydration. As the salt content in irrigation water increases, the interlayer space between DDLs expands and soil aggregates are disrupted by swelling and clay dispersion [41,75]. The more significant salinity damages the soil structure and the salt stress results in higher absorption of soil water [76].

#### *4.3. Screening of SWCC Models for Artificial Shelterbelt*

The fitting *R* 2 results of VG model and BC model were significantly higher than those of the G model, and the *RE* was generally lower than that of G model. Therefore, the fitting results of the VG and BC models were better for SWCC, while the BC model had better fitting effect for CK. However, the fitting effects of the VG model were better than the BC model for bare BL and US, especially for the surface soil. This is because the VG model considers more influencing factors when predicting. Therefore, it has a higher accuracy on soils with more complex physiochemical properties. The soil salt content and organic matter content of BL and US were maximum, so the fitting effect of the VG model is better than the BC model. Therefore, comparing the fitting results of the three models, the results concluded that the VG model is the best choice in this regard.

When using RETC software to predict the parameters of the VG model, only the influence of bulk density and texture can be considered. In our study, it can be found that the increase in organic matter and salinity reduce the fitting accuracy of the VG model, and through grey correlation analysis, we found that the soil bulk density, sand content, pH, clay content, organic matter, silt content, EC and other physiochemical indicators of the study area had a larger correlation degree with the parameters *a* and *n* of the VG model (Table 4). Meanwhile, the basic parameters such as organic matter, pH and EC were used as input variables to predict the parameters of the VG model, and our results found that the fitting effect of the VG model could be improved [77,78]. Therefore, pH, EC, organic matter and other indicators should also be used as input variables for calculation when predicting the parameters of the VG model.

#### **5. Conclusions**

The construction of an artificial shelterbelt with long-term saline water irrigation increased the water retention capacity and the content of fine pores in the soils, and reduced the content of macropores. Artificial shelterbelt construction had the greatest impact on the surface soil, and the impact gradually decreased with the increasing soil depth. Compared with the shifting sandy land (CK) in the 0–5 cm soil layer, the saturated water content *θ<sup>s</sup>* of the soil under BL and US increased by 4.42% and 12.67%, the field capacity *θ<sup>f</sup>* increased by 68.9% and 70.41%, and the available water content *θ<sup>a</sup>* increased by 87.84% and 70.97%. In the 5–10 cm soil layer, the *θ<sup>f</sup>* , *θ<sup>r</sup>* , *θ<sup>a</sup>* and *θra* of US and BL were higher than those of CK. In the 10–20 cm soil layer, the soil of US had the best water-holding performance. In the 20–30 cm soil layer, the soil water-holding performance of BL was the best. The available water content of each treatment was in the form of readily available water content. Comparing the coefficient of determination (*R* 2 ), root mean square error (*RSME*) and relative error (*RE*) of the three models, it was found that the G model always overestimated the soil water content and had a lower prediction accuracy, while the BC and VG models had higher prediction accuracies. Although both BC and VG models are suitable for fitting SWCC of the shelterbelt, the VG model is more effective. The parameters *a* and *n* of the VG model had a higher degree of correlation with the soil EC and organic matter. In summary, the increase in organic matter and salinity reduced the fitting accuracy of the model. When predicting the parameters of the VG model and establishing the soil transfer function, soil EC, organic matter and other indicators should be calculated as input variables.

**Author Contributions:** Conceptualization, C.M. and J.Z.; Data curation, L.T., M.T.J., Q.C. and J.F.; Formal analysis, C.M. and W.C.; Investigation, C.M., X.L., Q.C. and J.F.; Methodology, C.M., J.Z. and X.L.; Writing—original draft, C.M.; Writing—review and editing, J.Z. All authors have read and agreed to the published version of the manuscript.

**Funding:** This study was supported by the National Natural Science Foundation of China (No. 41877541, 41471222), Key Scientific and Technological Project of Shaanxi Province (2022ZDLNY02-03).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** All data reported here are available from the authors upon request.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


**Zhiwei Zhang 1,2, Huiyan Yin 1,2, Ying Zhao <sup>3</sup> , Shaoping Wang <sup>4</sup> , Jiahua Han 1,2, Bo Yu 1,2 and Jie Xue 5,6,7,\***


**Abstract:** Soil moisture is a vital factor affecting the hydrological cycle and the evolution of soil and geomorphology, determining the formation and development of the vegetation ecosystem. The previous studies mainly focused on the effects of different land use patterns and vegetation types on soil hydrological changes worldwide. However, the spatial heterogeneity and driving factors of soil gravimetric water content in alpine regions are seldom studied. On the basis of soil sample collection, combined with geostatistical analysis and the geographical detector method, this study examines the spatial heterogeneity and driving factors of soil gravimetric water content in the typical alpine valley desert of the Qinghai–Tibet Plateau. Results show that the average value of soil gravimetric water content at different depths ranges from 3.68% to 7.84%. The optimal theoretical models of soil gravimetric water content in 0–50 cm layers of the dune are different. The nugget coefficient shows that the soil gravimetric water content in the dune has a strong spatial correlation at different depths, and the range of the optimal theoretical model of semi-variance function is 31.23–63.38 m, which is much larger than the 15 m spacing used for sampling. The ranking of the influence of each evaluation factor on the alpine dune is elevation > slope > location > vegetation > aspect. The interaction detection of factors indicates that an interaction exists among evaluation factors, and no factors are independent of one another. In each soil layer of 0–50 cm, the interaction among evaluation factors has a two-factor enhancement and a nonlinear enhancement effect on soil gravimetric water content. This study contributes to the understanding of spatial heterogeneity and driving factors of soil moisture in alpine deserts, and guidance of artificial vegetation restoration and soil structure analysis of different desert types in alpine cold desert regions.

**Keywords:** geographical detector; alpine dunes; spatial heterogeneity; soil moisture; Qinghai— Tibet Plateau
