Moisture Migration and Recharge Pattern of Low-Permeability Thick Cohesive Soil in Northern Margin of the Jianghan Plain

: An extremely low hydraulic conductivity of cohesive soil causes a low transport rate of water and solute, with a time-consuming result, as we all know. Stable isotopes (δD and δ 18 O) and in situ monitoring systems of the data about soil water, rainfall, and groundwater were used to analyze the soil moisture migration pattern, using a conceptual model in the field test site, simulated by Hydrus 1D. The results show that multiple rainfalls' accumulations can cause the water to recharge from soil moisture to micro-confined groundwater, gradually. The soil moisture dynamic change is composed of a dehydration period and absorption period; the cohesive soil water content below 5.0 m was affected by the micro-confined groundwater level and dehydrated in advance due to the level decline. The thick cohesive soil profile can be divided into a shallow mixing zone (0–2 m), steady zone (2–5 m), and deep mixing zone (5–15 m). The effective precipitation re-charge was 234 mm and the average infiltration recharge coefficient ( Rc ) was 0.1389, but the water exchange between the cohesive soil moisture and groundwater was 349 mm in two hydrological years. This paper reveals the moisture migration and recharge pattern of low-permeability thick cohesive soil in a humid area with a micro-confined groundwater aquifer; this is of great significance for groundwater resources evaluation and environmental protection in humid climate plain areas.

Cohesive soil water and solute transport studies not only consume a great deal of time but are also hard to do [11].Various studies have attempted to use technologies such as hydrological geochemistry, stable isotopes, in situ monitoring, and numerical simulations to reveal the migration of soil moisture and the conversion relationships among rainfall, soil moisture, and groundwater [12][13][14][15][16][17][18][19][20][21][22][23][24][25][26][27][28].But almost all of these studies were located in arid and semi-arid areas with unconfined groundwater aquifers.The moisture migration and recharge pattern of thick cohesive soil with low hydraulic conductivity covering a micro-confined groundwater aquifer in a humid area is still unclear.
The Quaternary stratum is widely distributed in the Jianghan Plain and is covered with thick cohesive soil at a thickness of more than 10 m.Based on stable isotopes (D, 18 O) technology, using in situ monitoring technology, we found the moisture migration and recharge pattern of thick cohesive soil in the field test site with a micro-confined groundwater level change.Finally, we obtained the precipitation infiltration recharge amount and soil water-groundwater exchange amount through soil moisture transport numerical simulation by Hydrus 1D.This research is of great significance for groundwater resources evaluation and development and groundwater environmental protection in humid climate plain areas such as the Jianghan Plain.

Geological Setting
The scientific field test site (31°03′58.46″N, 113°55′55.47″E) is about 600 m from the urban area of Xiaogan City, which is located in the northern region of the Jianghan Plain.Economic crops such as peanuts, shallots, and cotton had been planted before in the field test site, but the land use type of the field test site had changed to bare land since 1 January 2018.The elevation of the field test site is 33.40 m.According to drilling data [10], we could determine the geological structure (Table 1) and distinguish that the sandgravel layer in the bottom part of the Quaternary Upper Pleistocene (Qp3 al ) stratum, covered with 15 m thick cohesive soil, forms the main aquifer in the field test site, and its initial stable groundwater level was 25.72 m with a 7.32 m pressurized water head.Considering the groundwater recharge, flow, and drainage mode, the Quaternary Upper Pleistocene (Qp3 al ) micro-confined aquifer on the eastern side of the Huan River mainly accepts rainfall infiltration recharge from in front of the eastern mountain and the lateral runoff recharge of Wudang Group (QbW2) weathered fissure water.Then, the micro-groundwater flows through the unconfined pore water aquifer in the form of horizontal runoff and finally discharges into the Huan River [29].The field test site is located at the transition zone between the Dabie Mountain Area and the Jianghan Plain.
Soils at different depths in the field test site were classified into different types (Table 1) by the percentage data of clay, silt, and sand from the soil particle analysis results.Normal physical properties data such as bulk density, hydraulic conductivity, and permeability were also tested in the laboratory and are shown in Table 1.

Monitoring Systems
The field test site included three monitoring systems [12]: a meteorological monitoring system, a soil monitoring system, and a groundwater monitoring system.The ground conditions were all bare land during the monitoring period from 1 January 2018 to 4 November 2020.Equipment models, factors, and monitoring frequency of the field test site are shown in Table 2.The soil monitoring system was aimed at monitoring the volume of soil moisture and soil water potential at different layers of the east side profile of WELL 3 (φ = 2 m) (named 3E) in the unsaturated zone (0-6 m) above the micro-confined groundwater level.The depths of the monitored layers in profile 3E were 0.2 m, 0.5 m, 0.9 m, 1.4 m, 2.0 m, 2.5 m, 3.0 m, 3.5 m, 4.0 m, 4.5 m, 5.0 m, and 6.0 m (Figure 1).The surface heat flux of the soil was also monitored at a depth of 0.1 m.The meteorological monitoring system was located at the northeast side of the site, and was mainly aimed at normal meteorological factors (Table 2) monitoring in order to evaluate the amount of moisture and energy received.The groundwater monitoring system was located at the southwest side of the field test site and was aimed at monitoring the dynamic changes in the micro-confined groundwater level.
The electrical energy to support the operation of the meteorological monitoring system and soil monitoring system was provided by a battery, which was charged from solar panels equipped by the meteorological monitoring system (Figure 1).The groundwater monitoring system had its own battery inside the probe.In situ monitoring data (soil volume water content, soil water potential, surface heat flux) of the soil monitoring system were measured from probes, which were buried horizontally in different monitoring layers in profile 3E.The CS257 probe measured the resistance value (R) when the water potential between the sensor and the soil reached an equilibrium state.The CS650 probe measured the soil bulk dielectric permittivity (Ka) by using the time domain reflectometry (TDR) method.Finally, the program self-contained inside the CR1000 data collector transformed them to soil water potential and soil volume water content.The relationship between dielectric permittivity and volume water content in mineral soils has been described by the Topp Equation (1980) in an empirical fashion using a third-degree polynomial, and it works well in most mineral soils, so a specific soil calibration of the CS650 sensor is usually not necessary.To obtain reliable data, we performed specific soil volume water content calibrations in the laboratory to confirm the equation relating Ka to θv.In situ monitoring data of the soil monitoring system and meteorological monitoring system were collected by a CR1000 data collector (Figure 1) every 10 min.The groundwater level was monitored in situ and recorded by probes (Levelogger III 3001) in SYC-01 every 10 min.The data inside were exported by Solinst software when the probe was connected with the computer.

Sampling
The water samples described in this paper include information about the soil moisture of profile 3E, rainfall from the field test, and the groundwater of SYC-01 from the field test.We also took soil samples for water content testing in a vertical direction inside the field test site.

O) testing
We put a plastic bucket on the open-air roof of resident's house near the field test site and collect rainfall samples after each rainfall.The resting rainfall in the plastic bucket was all poured out after sampling every time.The rainfall sample was collected into clean PET bottles (their volume is 50 cm 3 ) without bubbles after the rain stopped.Then, we sealed the PET bottle with elastic film to ensure that the rainfall sample was isolated from air.All the rainfall samples collected for stable isotopes (D, 18 O) testing were saved at a low temperature, about 3-4 °C.

O) testing
We put down a stainless steel sampler with the function of a water inlet from its bottom into SYC-01 at the depth of 16 m and pulled up a full sampler of groundwater.Then, the groundwater sample was collected from the stainless steel sampler into clean PET bottles (their volume is 50 cm 3 ) without bubbles.Finally, we sealed the PET bottles with elastic film to ensure that the rainfall sample was isolated from air.All the groundwater samples collected for stable isotopes (D, 18 O) testing were saved at a low temperature, about 3-4 °C.
(3) Soil moisture samples for stable isotopes (D, 18 O) testing It was not easy to collect soil moisture stable isotopes (D, 18 O) test samples.First of all, we used a drilling machine (QY-100L, Jiangsu Wuxi Mineral Exploration Machinery General Factory Ltd., Wuxi, JiangSu, China) to obtain a soil column in the vertical direction inside the field test site.In order to avoid taking contaminated samples, we needed to cut a 2 cm thickness of the soil on the surface part of the soil column away before taking the soil samples.Then, we used glass bottles (their volume is 8 cm 3 ) to take clean soil samples from the soil column every 0.1 m from 0 to 14.30 m.The glass bottles full of the soil samples were sealed by elastic film isolated from air and transported to the laboratory at a low temperature, about 3-4 °C.
In the laboratory, we put medical-use cotton balls into the sampled glass bottles (their volume is 8 cm 3 ) in preparation for soil water extraction.Then, we put the sampled glass bottles (their volume is 8 cm 3 ) into an LI-2100 automatic water extraction system (LICA Inc., Beijing, China) for soil water extraction.The LI-2100 automatic water extraction system created a 1500 Pa vacuum environment with the temperature set at 105-110 °C for 300 min for every soil water extraction process.At the same time, we also used a condensation vapor collection part of the same vacuum environment.The soil moisture extracted by LI-2100 included strongly bound water, weakly bound water, capillary water, and gravity water.
We extracted soil moisture samples into new glass bottles (their volume is 2 cm 3 ) without bubbles from the soil samples.Then, we sealed each new glass bottle with elastic film to ensure that the water sample was isolated from air.All the soil moisture samples collected for stable isotopes (D, 18 O) testing were saved at a low temperature, about 3-4 °C.
(4) Soil samples for soil volume water content testing Before taking soil samples, we weighed, recorded, and numbered every standard steel ring knife and every porous aluminum box.We used standard steel ring knives (φ61.8 × 20 mm, their volume is 60 cm 3 ) to take standard soil samples at the same depths of each soil column when collecting soil samples.After cutting a 2 cm thickness of soil on the surface part of the soil column away, we cut down into each vertical soil column (about φ70 mm × 50 mm), and then put a standard steel ring knife fully inside it.Then, we used a soil cutter to flatten both sides of the standard steel ring knife and put it inside a porous aluminum box.Following this, we weighed and recorded the weight of every standard steel ring knife full of a soil sample with a balance (±0.01 g).Finally, we put every standard steel ring knife full of a soil sample inside a numbered porous aluminum box in preparation for volume water content testing.

Testing (1) Stable isotopes value testing
The stable hydrogen and oxygen isotopes were analyzed via LGR water isotope analysis (LWA-45EP, LGR Inc., Scottsdale, Arizona, USA) at the Wuhan Center of Geological Survey.The analysis method of this machine is high-resolution laser spectroscopy absorption measurement analysis.The test items were δD and δ 18 O, which were expressed as δ(‰) = (Rsample/Rstandard − 1) × 1000‰, where Rsample and Rstandard reflect the sample and standard's heavy-to-light isotope ratios, respectively, at the Wuhan Center of Geological Survey, Wuhan, China.The standard calibrated to worldwide reference mate-rials, V-SMOW (Vienna Standard Mean Ocean Water, Vienna, Austria), was utilized for computation, and the analytical accuracies for δ 18 O-H2O and δD-H2O were ±0.5‰ and ±0.1‰, respectively.
We tested the standard water sample supplied by LGR Company before starting the test of the water samples.The respective analytical accuracy for δ 18 O-H2O and δD-H2O of this standard water sample test should be less than ±0.5‰ and ±0.1‰.We ranked standard water samples supplied by LGR Company after every three samples, and then checked the respective analytical accuracy via software self-contained inside the LGR water isotope analysis equipment in order to ensure the reliability of testing data.
(2) Soil volume water content testing First of all, we put the numbered porous aluminum box and standard steel ring knife full of the soil sample inside together into a baking oven to dry at a temperature of about 105-110 °C for about 12 h.We put the soil sample inside a dryer for about 20 min to cool to room temperature and then weighed it.We repeated this process until the differences among the last three weights were less than 0.05 g.Next, we put the numbered porous aluminum box and standard steel ring knife without the soil sample inside together into a baking oven to dry at a temperature of about 105-110 °C for about 4 h.We repeated this process until the differences among the last three weights were less than 0.05 g.
We were able to obtain the weights of the dry soil sample (W2) and natural soil sample (W1) after the processes described above.The volume of the soil sample (V) was 60 cm 3 .The soil volume water content (SVWC) could be obtained by the following formula: where W1 is the weight of the natural soil sample (M); W2 is the weight of the dry soil sample (M); ρ is the water density at 20 °C (ML −3 ); V is the volume of the soil sample (standard steel ring knife) (L 3 ).

Soil Moisture Migration Stratification
The data (soil volume water content and soil water potential) of the 3E profile, rainfall data, and groundwater level data were monitored from 1 January 2018 to 14 November 2020, and were used to show the dynamic changes (Figure 2) of the soil volume water content and soil water potential as they were affected by rainfall and the groundwater level over time.According to the soil volume water content dynamic change diagram (Figure 2), we were able to determine that abundant rainfall can make the soil volume water content remain at a high level all the time, except for the shallow layers during a drought period.We can see that there was an obvious stratification phenomenon within the same period, and the soil profile of the field test site can be mainly divided into zones according to this phenomenon by soil volume water content dynamic changes at different depths in profile 3E (Figure 2).The vadose zone of the field test site included the following zones: the sensitive zone of rainfall infiltration (0-1.4 m), the buffer zone of rainfall infiltration (1.4-3.5 m), and the migration zone of rainfall infiltration (3.5-5.0 m).The changing trend of the soil moisture of soil layers at 5.0 m and 6.0 m was similar to the 4.0 m and 4.5 m soil layers, with an obvious difference, but the same could be observed with the groundwater level's rising and falling.The obvious difference was mainly caused by the groundwater level, which is consistent with the previous research results [15].Finally, we regarded the soil layer at the depths from 5.0 m to 15.0 m as the rainfall infiltration and groundwater level co-influenced zone.Cohesive soil exhibits the characteristic of threshold gradient I0 in many studies [30][31][32], which can cause a pressure head loss.The capillary rise in pores increases the influence range of groundwater to the upper cohesive soil layer.In order to analyze the influence of micro-confined groundwater level changes on cohesive soil volume water content at depths of 5.0 m and 6.0 m, we transformed the micro-confined groundwater level to the elevation of the capillary zone roof (H) to take threshold gradient I0 and the capillary rise in the pores into consideration.The transformation formula is expressed as H = H0/(1 + I0) + H1 + H2; where H0, H1, H2, and H0/(1 + I0) reflect the micro-confined groundwater level, the roof elevation of the micro-confined aquifer (18.40 m), the capillary rising distance, and the thickness of the water-bearing zone of cohesive soil.Not only the geological age, lithology, and geologic structure, but also the physical properties, clay mineral components, and main chemical constituents of cohesive soil in the field test site are almost the same as in the Nierji area [33].Because of these similarities, we decided that the parameter values of cohesive soil (I0 = 0.08, H2 = 1.65 m) in our field test site would be referred to the field test data in the Nierji area.
In order to analyze whether or not the soil moisture in the depth range of 5.0 m to 6.0 m was affected by the micro-confined groundwater level, when threshold gradient I0 and soil moisture in pores were taken into consideration, the soil moisture data (5.0 m layer and 6.0 m layer), rainfall data, and capillary zone roof elevation data (H) during the monitoring period were used to plot a dynamic change diagram (Figure 3) for analysis.

Moisture Migration Rules of Thick Cohesive Soil
Two infiltration pathways occur in cohesive soil layers: piston (or uniform) flow and preferential flow.Piston flow leads to stable wetting fronts parallel to the soil surface, whereas preferential flow results in irregular wetting.Preferential flow usually occurs in the shallow cohesive layer, which includes macropores, cracks, and roots.After reaching a certain depth, the preferential flow disappears, and only piston flow remains below this certain depth [34][35][36].
According to the dynamic change in soil volume water content (Figure 3), we could see the response of the soil volume water content in the shallow layers above 0.9 m to the rain almost immediately, which, by means of preferential flow, occurred within the depth of 0.9 m.Then, the cohesive soil water content change in soil layers below 0.9 m occurred more slowly, with a lag [37], to the rain than did shallow surface layers above 0.9 m.There was only piston flow remaining in cohesive soil layers below the depth of 0.9 m.Dynamic changes in soil volume water content in the profile of 3E showed an "increase-decrease-increase-decrease" wave-cycle trend (Figure 2).We divided the wave-cycle trend of soil volume water content into a "dehydration period" and "water absorption period" during the monitoring period.This indicates that water recharge was absorbed in cohesive soil pores during the water absorption period, and water in cohesive soil pores was lost during the dehydration period.Two transition times between the dehydration period and the water absorption period in the cohesive soil layers in the field test site were included from June 2019 to April 2020.Information about transition time, lag time, response time of soil moisture, and groundwater level (GWL) is shown in Table 3 and Figure 4.The initial time when the soil moisture changed from the water absorption period to the dehydration period was the moment when the free infiltration process started without rain falling.The soil moisture of the cohesive soil layer near the surface began to decrease and continued to migrate downward when the transition started.According to the transition time of each soil layer (Table 3), we found that soil layers above the depth of 4.5 m (Figure 4) continued infiltrating, resulting in water loss from the surface to the bottom in the unsaturated zone during the dehydration period without rainfall.The groundwater level started to fall 1 day and 15 h after the rain stopped, followed by a moisture decrease in soil layers below 4.5 m, then the moisture decreased in soil layers in the depth range of 1.4-4.5 m.Water release from the deeper soil layer to groundwater in the field test site was caused by a pressure head decrease at the bottom of the soil layer [38,39] due to the groundwater level falling.The capillary saturation zone roof moved down, enhancing the soil water's decrease in the unsaturated zone.It can be seen from the dynamic changes in soil water content and groundwater level on 1 July 2019 that the micro-confined groundwater level's falling can affect the soil volume water content of the soil layer at the depth of 5.0 m and cause an earlier transition time for soil layers below 4.5 m than that of the soil layers in the depth range of 1.4-4.5 m.
The initial time when the soil moisture changed from the dehydration period to the water absorption period was the moment when rainfall infiltration was greater than evaporation.In the early stage of infiltration, water was mainly absorbed by soil particles to increase soil moisture content and did not move down until a positive water potential gradient was reached.The soil moisture began to infiltrate and migrate downward after the early stage of infiltration to the lower soil layers.The lower layers repeated the soil moisture change process of the upper soil layer; the soil water moved downward gradually layer by layer.Heavy rain after a long drought period brought layers of thick cohesive soil into the water absorption period, and the soil moisture response showed a lag from surface to bottom in thick cohesive soil.The lag time increased from surface to bottom, indicating that the groundwater level rise on 6 January 2020 contributed little to the increase in soil water content above the depth of 6.0 m.The increase in soil water content above the depth of 6.0 m was mainly caused by the accumulation of rainfall infiltration.
The monthly rainfall increased gradually in the field test site from January 2020 to July 2020.The monitoring layers of thick cohesive soil were all in the water absorption period, a continued increase in soil volume content (Figure 5).The total water potential gradient value (Iusa) of the profile was positive and close to 1 (Figure 5), which means that soil water continued moving downward throughout the whole profile.The climate of the field test site is a subtropical monsoon climate with abundant rain; soil water in the unsaturated zone can infiltrate gradually into the saturated zone and finally recharge the groundwater under multiple rainfalls' accumulation.The "water absorption period" and "dehydration period" cycle of the thick cohesive soil had a significant buffering effect on groundwater recharge.The infiltration recharge was stored in thick cohesive soil in the "water absorption period" and then released into the groundwater aquifer in the "dehydration period".This can help reduce the influence degree on groundwater of drought via the release of soil moisture to maintain the water recharge to the groundwater aquifer, which increases the recharge amount of groundwater in drought periods.On the other hand, it can also help in obtaining much more groundwater recharge from the weathered fissure water aquifer (QbW2) and storing it in thick cohesive soil in the "water absorption period" due to the buffering effect.Many more available water resources can be provided by groundwater aquifers for water resource management in the Jianghan Plain area because of this "water absorption period" and "dehydration period" cycle.This finding will be helpful for lowpermeability medium water recharge evaluation.

Vertical Water Movement Conceptual Model of Thick Cohesive Soil
In order to reveal the variation characteristic details of stable isotopes (D, 18 O) and soil volume water content in the thick cohesive soil profile to a large extent, the vertical interval of the soil samples was set to 10 cm.It was obvious that the thick cohesive soil profile in the test site can be divided into three zones (0-160 cm, 160-670 cm, and 670-1430 cm), with different variation characteristic details (Figure 6).The depth range of the first zone was 0-160 cm (Figure 6).Evaporation effects frequently can enrich δ 18 O and δD values, and frequent precipitation effects can impoverish δ 18 O and δD values [40].There was abundant rain and strong evaporation because of the high temperature at the field test site.The δD and δ 18 O values of soil water changed greatly with no obvious trend because of the effects of the frequent alternation of precipitation and evaporation (Figure 6).The maximum depth of the evaporation effect reached 160 cm.
The depth range of the second zone was 160-670 cm (Figure 6).Soil water content was not affected by evaporation but by precipitation infiltration, which caused stable effective recharge to soil layers below by the means of piston flow [41].The second zone could be divided into two parts because of the δ 18 O and δD values' range caused by the difference in soil particle composition.The δD and δ 18 O values of the first part (160-440 cm) ranged from −49.9‰ to −43.3‰ and from −6.66‰ to −5.99‰, and the second part (440-670 cm) ranged from −51.7‰ to −46.1‰ and from −7.26‰ to −6.60‰.The higher the soil clay content, the more light water molecules with relatively rich δD and δ 18 O values the soil absorbed, replacing the heavy water molecules of bound water with relatively depleted δD and δ 18 O values, and the more obvious was the isotope fractionation [26,42].The higher the soil clay content, the higher the bound water content percentage was and the lower the content percentages of capillary water and gravity water were.Because the soil moisture extracted by the LI-2100 automatic water extraction system included bound water, capillary water, and gravity water [43], soil water δD and δ 18 O values of the first part (160-440 cm) were much richer than those of the second part (440-670 cm).The soil water δD and δ 18 O values inside the same part changed a little.The δD and δ 18 O values represented the mixing results of multiple rainfalls and evaporation effects during the wet season or dry season in a year.
The depth range of the third zone was 670-1430 cm, and the δD and δ 18 O values increased, generally in the range of −50.28‰ to −42.76‰ and −7.31‰ to −5.16‰, with the increase in depth (Figure 6).The elevation range of the groundwater level was 22.10-27.40m with a maximum amplitude of 5.30 m.Groundwater and the soil water above recharged together to the soil water when the groundwater level increased, whereas the soil water recharged to the groundwater.The δD and δ 18 O values of groundwater were richer compared to the soil water of cohesive soil in the test site.Both groundwater infiltration upward and soil water infiltration downward are slow-changing processes.The recharge from groundwater to soil water enriched the δD and δ 18 O values of the soil water after mixing.The closer to the roof of the groundwater aquifer it was, the richer were the δD and δ 18 O values of the soil water that was affected by groundwater in general (Figure 6).There was an obvious depleted isotope value part of the third zone from 1060 cm to 1230 cm that was affected by a higher soil clay content, which was the same as in the second part (440-670 cm) of the second zone in the field test site.The δD and δ 18 O values of soil water from 1300 cm to 1430 cm were relatively stable in the ranges of −45.47‰ to −42.76‰ and −5.94‰ to −5.29‰, close to the groundwater, because of the almost direct influence of the groundwater level change.
According to the variation characteristics of stable hydrogen and oxygen isotopes in the thick clay soil profile and factors influencing soil moisture migration, the thick cohesive soil profile in the test site can be divided into three zones (Figure 7): the shallow mixing zone, the steady zone, and the deep mixing zone [44].In the first zone, named the shallow mixing zone, the liquid soil water migrated downward, and the gaseous soil water migrated upward in the form of water vapor.The vaporization front was located at a depth of 30 cm.The frequent alternation of rainfallevaporation caused the δD and δ 18 O values of the soil water below the vaporization front to vary greatly, showing the characteristics of vacillation.The depth of the bottom of the shallow mixed zone was affected by evaporation, with some differences in different periods, but the maximum depth was not deeper than the maximum evaporation depth of 200 cm [13].
The second zone, named the steady zone, ranged from the bottom of the shallow mixing zone to the top of the capillary saturation zone without evaporation.Soil water content was not affected by evaporation but by precipitation infiltration, which caused stable effective recharge to soil layers below by means of piston flow.Peaks and troughs of δD and δ 18 O values can indicate climate change characteristics during the historical infiltration period to a certain extent, but this indication is only applicable to soil layers with small changes in soil texture.
The third zone, named the deep mixing zone, ranged from the top of the capillary saturation zone to 1500 cm.Only piston flow existed in the whole zone.Groundwater and soil water above recharged together to this zone when the groundwater level increased, whereas the soil water of this zone recharged to groundwater.
The thick cohesive soil was recharged by rainfall infiltration, and the land use type was bare land with no crops at the field test site.Soil water was excreted into the atmosphere in the form of evaporation from within a depth of 200 cm [13].The change in cohesive soil moisture in this depth range was affected by both precipitation and evaporation.The change in cohesive soil moisture in the vadose zone at the depth range of 200-500 cm was only affected by infiltration from cohesive soil layers above 200 cm by piston flow.Due to the high clay content and low permeability of cohesive soil in the depth range of 440-570 cm, there was a certain retardation effect on the precipitation infiltration.Due to this retardation effect and the stable water flux and infiltration rate of the soil, the layers at the depth range of 500-1100 cm were named the unsaturated-saturated variable zone.The depth of the capillary saturation zone roof was affected by the microconfined groundwater level.The cohesive soil water content in this zone was affected by infiltration downward from the upper soil layer and infiltration upward from the lower soil layer by means of piston flow in the saturated zone caused by groundwater level elevation.There was a stable saturated zone of micro-confined groundwater in the range of 1100-1500 cm; the cohesive soil water content was affected by infiltration downward from the upper soil layer and infiltration upward from the lower soil layer by means of piston flow.

Mathematic Model
The land use type of the field test site was bare land, and the water absorption rate of the root system (S) was set to 0. One-dimensional vertical flow was mainly considered for vertical groundwater recharge.The Richards equation used to describe the vertical percolation in a vadose zone is as follows.
where h is pressure head (L); c(h) is the specific moisture capacity (L −1 ); K(h) is the unsaturated hydraulic conductivity (LT −1 ); z is the vertical coordinate (L) with zero referring to the ground and upwards being positive; t is time (T).The Van Genuchten-Mualem soil hydraulic functions were used to describe the unsaturated hydraulic properties in the models.

Initial Conditions and Boundary Conditions
The upper boundary was the atmospheric boundary condition with the surface layer defined by the evaporation and precipitation at the soil surface.The lower boundary condition was the variable pressure head defined by the micro-groundwater pressure head at the depth of 1500 cm.
Soil water potential data of the monitoring layers were used as the initial condition of the simulating profile.The upper boundary condition was controlled by data (precipitation, water surface evaporation, net solar radiation, humidity, air temperature, wind speed, and wind direction) collected by the meteorological monitoring system at the field test site.The lower boundary condition was controlled by the pressure head transferred from the auto-monitored groundwater level in the field test site.

Spatial and Temporal Discretization
The 0-1500 cm simulated soil profile was divided into nine layers and four mass balance zones (0-200 cm, 200-500 cm, 500-1100 cm, and 1100-1500 cm) according to the soil properties and the hydrogeological conceptual model of the field test site.The soil profile was divided into 300 units at equal intervals of 5 cm, and had 301 nodes and eight observation points (20 cm, 90 cm, 200 cm, 250 cm, 500 cm, 600 cm, 1100 cm, and 1500 cm).
The simulated period of 855 days ranged from 1 July 2018 to 1 November 2020.Time discretization was used in the simulation, and the interval of the time discretization was gradually adjusted according to the number of iterations of the convergence.The initial time interval was set to 1 × 10 −5 day; the maximum time step was set to 1 day.

Calibration and Validation
We divided the soil profiles into nine layers, and the basic parameters for soil water migration simulation include θs, θr, α, n, Ks, and l.Soil particle size distribution data and bulk density data (Table 1) of these nine layers were analyzed in the laboratory.We estimated θr, Ks, and the empirical shape parameters n and α using the H3 mode (Sand, Silt, Clay and Bulk Density; SSCBD for short) of Rosetta Lite (Salinity Laboratory of the USDA, Riverside, CA, USA) based on the neural network embedded in HYDRUS-1D from the data analyzed in the laboratory (Table 1).θr, Ks, n, and α were fitted parameters; the soil's mechanical composition (sand, silt, and clay), θs, and ρb were fixed parameters with specific values derived from measurements.
The initial values of the fitted parameters (Table 4) for each layer were estimated using the Neural Network Prediction embedded in the model, and then, the parameters were fitted by fitting the observations measured for this layer in the field test site.We compared observed field measurements with the results of the HYDRUS-1D simulations using the mean absolute errors (MAEs) and root mean square errors (RMSEs) shown in Table 5 [18,22].The average value of the simulation data is very close to the average value of observed data at each layer.The values of MAE and RMSE of these four layers are small and close to 0, which means the accuracies of model calibration and validation are good.The observed field measurements and simulation results of the four cohesive layers in the three mass balance zones were roughly the same (Figure 9).The difference between simulation data and observed data (|Oi − Pi|max) was larger at the depth of 0.2 m than at 0.9 m, followed by that at 6.0 m and 2.5 m.The numerical model can generally simulate the variation in soil water content in cohesive soil.The parameters of the Van Genuchten-Mualem equation were finally determined after multiple hydraulic parameter adjustments.The optimization results are presented in Table 6.The thick cohesive soil profile was divided into four mass balance zones according to the vertical hydrogeological conceptual model of the field test site (Figure 8).The depth range of the first mass balance zone was 0-200 cm, which belonged to the vadose zone.The change in cohesive soil moisture in this depth range was affected by precipitation and evaporation.Most of the water infiltrated into the vadose zone from precipitation returned back to the atmosphere by evaporation; only a small percentage of water was effective and recharged the cohesive soil moisture in the vadose zone [45].The depth range of the second mass balance zone was 200-500 cm, which belonged to the vadose zone.Water flux at the depth of 200 cm was the effective infiltration recharge from precipitation to soil water in the vadose zone.Water flux at the depth of 500 cm was the sum of the effective recharge from precipitation infiltration and the water release of cohesive soil in the vadose zone at the depth range of 200-500 cm.The depth range of the third mass balance zone was 500-1100 cm, also called the unsaturatedsaturated variable zone.It can be considered that the micro-confined groundwater had been recharged when it arrived at this zone, and the water flux at the bottom of this mass balance zone was the actual recharge to the micro-confined groundwater when the soil moisture migrated to the bottom of this zone.The depth range of the fourth mass balance zone was 1100-1500 cm, which was a stable saturated zone of micro-confined groundwater.Water flux at the bottom of this mass balance zone was the exchange capacity of soil water-groundwater.
Surface runoff could be formed when the rainfall was greater than 25 mm/d, which was mainly formed in the wet season from June to July.The cumulative rainfall was 1684 mm during two hydrological years, from 1 July 2018 to 30 June 2020.The cumulative surface runoff which discharged out to surface water was 272 mm (16.15%).The cumulative rainfall infiltration was 1382 mm (82.07%), but most of the rainfall infiltration returned back to the atmosphere by evaporation, having accumulated to 1190 mm, accounting for 70.67% (Figure 10).The cumulative water flux of the clay layer at the depth of 200 cm reflected the effective precipitation infiltration recharge that could be totally absorbed in the cohesive soil and finally recharged to the groundwater after removing the water amount of evaporation.Due to the low permeability coefficient of cohesive soil in the test site and the long recharge path, there was an obvious lag effect on precipitation infiltration recharge to the micro-confined groundwater aquifer.We can see from Figure 8 that the soil moisture in the vadose zone at the depth range of 200-500 cm was only affected by infiltration from cohesive soil layers above 200 cm by piston flow.So, the cumulative water flux of the cohesive soil layer at the depth of 200 cm can be regarded as an effective recharge of precipitation for the vadose zone.The cumulative water flux at the depth of 200 cm found through numerical simulation was 234 mm, so the effective recharge of precipitation to the vadose zone was 234 mm, accounting for 13.89% of the cumulative rainfall amount (1684 mm) during two hydrological years.The precipitation infiltration recharge coefficient we obtained is 0.1389.Compared with the experienced value (0.10-0.15) of the precipitation infiltration recharge coefficient (Rc) of cohesive soil in Anhui listed (Table 7) in the Handbook of Hydrogeology [46], we found that an approximate result was given by our model.This comparison confirmed the correctness of the simulations.Hengshui and Luancheng are located in the North China Plain.The recharge coefficient (Rc) of cohesive soil in the Jianghan Plain we obtained is lower than that in the North China Plain; this was mainly caused by its lower permeability.The constant soil water saturation zone comprised cohesive soil layers below the depth of 1100 cm during the simulation period.The cumulative water flux of the cohesive soil layer at the depth of 1100 cm coincided with the cumulative water flux at the lower boundary (1500 cm), which means that the cumulative water flux of the cohesive soil layer at the depth of 1100 cm could be considered as the actual recharge from the cohesive soil layer to micro-confined groundwater.During the simulation period of the two hydrological years, the cumulative water flux of the cohesive soil layer at the depth of 1100 cm was 349 mm, so the recharge of the thick cohesive soil moisture to the microconfined groundwater was 349 mm, which was greater than the cumulative water flux of the cohesive soil layer at the depth of 200 cm.The difference between the recharge values of the soil layers at these two depths indicated that thick cohesive soil released water to recharge groundwater in the dry year in the field test site.
From the cumulative water flux at different depths of the thick cohesive soil in the test site shown in Figure 10, it can be found that with an increase in depth, the cumulative water flux of the cohesive soil layer at different depths shows an increasing trend and a periodic approaching phenomenon in terms of years, indicating that the thick cohesive soil water in the test site migrates downward, generally, to recharge microconfined groundwater.
The effective precipitation infiltration recharge decreased with less precipitation in dry years (2018 and 2019), and the same result was found with respect to the pressure head of the micro-confined aquifer in the field test site.The capillary saturation zone moved down, and the thickness of the vadose zone increased.Soil water of cohesive soil layers below 200 cm was released to charge deeper soil layers and the groundwater aquifer, which caused the total water storage of the thick cohesive soil profile to decrease.The effective precipitation infiltration recharge increased with more precipitation in wet years (2020), and the same result was found with respect to the pressure head of the micro-confined aquifer in the field test site.The capillary saturation zone moved up and the thickness of the vadose zone decreased in wet years.With the effects of "soil water of upper layer moves down and groundwater level moves up", the soil water content of the cohesive soil layers and total water storage of the thick cohesive soil profile increased in general.
Taking the depth of 500 cm as the boundary, the cumulative water flux of the cohesive soil layer in the upper vadose zone was mainly affected by precipitation and evaporation, and the cumulative water flux of the deeper cohesive soil layer was mainly affected by micro-confined groundwater level changes (Figure 10).The annual exchange capacity of soil water and groundwater was controlled by the micro-confined groundwater pressure head to a great extent.The exchange capacity decreased as the micro-confined groundwater pressure head increased, leading to a decrease in the actual recharge rate of the groundwater.On the contrary, the actual recharge of the groundwater increased.The trend variation in the cumulative water flux of cohesive soil at different depths clearly revealed the water storage adjustment effect and water migration rules of the thick cohesive soil profile.

Conclusions
(1) The upward movement of the capillary zone roof caused by a rise in the groundwater level can increase the soil volume water content of layers below the depth of 5.0 m.The silty clay layer at the depths from 4.4 m to 5.7 m has a certain retardation effect on the cohesive soil water infiltration, which reduced the soil moisture content's change in the cohesive soil layer above the depth of 5.0 m caused by groundwater level changes.The soil profile of the field test site could be mainly divided into zones: the sensitive zone of rainfall infiltration (0-1.4 m), the buffer zone of rainfall infiltration (1.4-3.5 m), the migration zone of rainfall infiltration (3.5-5.0 m), and the rainfall infiltration and groundwater level co-influenced zone (5.0-15.0m).The ef-fects of micro-confined groundwater level dynamic changes were taken into consideration when dividing into zones.(2) We divided the dynamic change in soil volume water content into a "dehydration period" and a "water absorption period".Soil layers above the depth of 4.5 m continued infiltrating to lose water from the surface to the bottom in the unsaturated zone during the dehydration period without rainfall.The moisture decrease in the soil layers below 4.5 m started in advance in the 1.4-4.5 m range in the "dehydration period", affected by a groundwater pressure decrease.The soil moisture begins to infiltrate and migrate downward after reaching a positive water potential gradient, and then soil water continues moving downward through the whole profile gradually, layer by layer.The soil moisture response time of different layers shows a lag from the surface to the bottom in thick cohesive soil.Soil water in the unsaturated zone can infiltrate gradually into the saturated zone and finally recharge the groundwater under multiple rainfalls' accumulations.(3) The thick cohesive soil profile in the field test site can be divided into three zones.
The first zone, named the shallow mixed zone, was affected by precipitation and evaporation.The depth of its bottom was affected by evaporation, but it was never deeper than 200 cm in any period.The second zone, which had no evaporation, named the steady zone, ranged from the bottom of the shallow mixing zone to the top of the capillary saturation zone.Soil water content was only affected by precipitation infiltration, which caused stable effective recharge to the soil layers below by means of piston flow.The third zone, with only piston flow, named the deep mixing zone, ranged from the top of the capillary saturation zone to 1500 cm and was introduced to form the new conceptual model of the vertical movement of cohesive soil water when taking micro-confined groundwater level dynamic changes into consideration.Soil water above this zone and groundwater recharged together to the deep mixing zone when the groundwater level increased due to rainfall, and the soil water of this zone recharged to the groundwater.This research introduced a new mode as a reference for groundwater recharge studies.(4) The cumulative water flux at the depth of 200 cm was the effective infiltration recharge from precipitation to soil water in the vadose zone, which can be totally absorbed in cohesive soil and can finally recharge to groundwater.The precipitation infiltration recharge coefficient we obtained through simulation is 0.1389.The cumulative water flux of the cohesive soil layer at the depth of 1100 cm (349 mm) could be considered as the actual recharge from the cohesive soil layer to microconfined groundwater, which is greater than that of the depth of 200 cm (234 mm).The difference was not only because of a significant lag effect caused by a low saturated permeability coefficient and the long recharge path of the cohesive soil, but also by the exchange capacity, controlled by micro-confined groundwater pressure to a great extent.The exchange capacity decreased as the micro-confined groundwater pressure head increased, leading to a decrease in the rate of the actual recharge of groundwater.On the contrary, the actual recharge of groundwater increased.(5) The actual recharge to groundwater is the water flux exchange at the roof of the stable saturated zone.The change in the total soil moisture storage in the cohesive soil layers between the layer at the depth of 200 cm and the roof of the stable saturated zone is a non-negligible part and should be taken into consideration when evaluating the recharge of groundwater.Because of the buffering effect of the thick cohesive soil in the field test site, the water amount recharged to the groundwater in wet years was less than that in dry years, in general.On the other hand, although the permeability of thick cohesive soil is low, the recharge coefficient we obtained through the simulation is 0.1389, which means pollutants can infiltrate through soil layers from the surface to cause groundwater pollution.The cumulative effect of multiple rainfalls will accelerate the vertical transport process of water and pollu-tants in cohesive soil with high water content.In order to prevent groundwater pollution, pollution sources' control measures should be developed in areas covered with thick cohesive soil.

Figure 1 .
Figure 1.Monitoring and data collection systems of the field test site.The number "①-⑧" in the north-eastern part of this figure are the numbers of horizontal holes for instrument installation.The CS257 and CS650 probes were installed inside holes "⑤-⑥" at the east direction of WELL 3.

Figure 2 .
Figure 2. (a) The dynamic changes in rainfall, groundwater level, and soil volume water content of monitoring layers in profile 3E; (b) the dynamic changes in rainfall, groundwater level, and soil total water potential of monitoring layers in profile 3E.

Figure 3 .
Figure 3.The dynamic change diagram of rainfall, soil volume water content, and capillary zone roof elevation in profile 3E.Light blue dashed line and light orange dashed line represent the soil layer elevations at the depths of 6.0 m and 5.0 m.

Figure 4 .
Figure 4. Transition times of soil water at monitoring layers and groundwater response times to rain in the field test site.The vertical dashed line in color referring to the soil volume water content dynamic data of the cohesive soil layer stands for the water absorption period and dehydration period transition time of this layer, whereas the vertical solid line stands for the dehydration period and water absorption period transition time.GWL stands for groundwater level.

Figure 5 .
Figure 5. Soil volume water content and total water potential of thick cohesive layers are shown in the test site from January to July 2020.The left side is the soil water volume content's dynamic change in monitoring layers on the 15th of each month from January to July 2020, and the right side is the total water potential.

Figure 6 .
Figure 6.Vertical variations of stable isotope values (δD, δ 18 O) and soil volume water content in the test site on 5 July 2020.The red dashed lines are dividing lines of different zones.

Figure 7 .
Figure 7. Conceptual model of vertical movement of cohesive soil water in the field test site (modified from Yuan et al., 2012).The blue dashed line stands for the micro-confined groundwater level with certain fluctuations within 530 cm.The deeper the blue color is, the higher was the soil mois-

Figure 8 .
Figure 8. Vertical hydrogeological conceptual model of the field test site.

Figure 9 .
Figure 9.Comparison between observed field measurements and simulation result of 4 critical layers in 3 mass balance zones.

Figure 10 .
Figure 10.Water migration simulation results of critical depths of mass balance zones.

Table 1 .
Geological structure and soil texture information of test site.

Table 2 .
Soil, meteorology, and groundwater monitoring systems information.

Table 3 .
Information about transition time, response time, and lag time of soil water and groundwater.

Table 4 .
Initial hydraulic parameters of cohesive soil layers we obtained via the H3 mode (SSCBD mode) of Rosetta Lite for simulation.

Table 5 .
The accuracy assessment factors of model calibration and validation.O is the average value of observed data; P is the average value of simulation data; Oi is the observed data; Pi is the simulation data.

Table 6 .
Calibrated hydraulic parameters of cohesive soil layers used in soil water migration simulation.

Table 7 .
Comparison of recharge coefficients.