Modeling of the Hydrological Processes in Caatinga and Pasture Areas in the Brazilian Semi-Arid

: The semi-arid regions of northeastern Brazil have historically suffered from water shortage. In this context, monitoring and modeling the soil moisture’s dynamics with hydrological models in natural (Caatinga) and degraded (Pasture) regions is of fundamental importance to understand the dynamics of hydrological processes. Therefore, this work aims to evaluate the hydraulic parameters in Caatinga and Pasture areas using the Hydrus-1D inverse method. Thus, ﬁve soil hydraulic models present in Hydrus-1D were used, allowing the comparison of the single-porosity model with more complex models, which consider the dual porosity and the hysteresis of the porous medium. The hydraulic models showed better adjustments in the Caatinga area (RMSE = 0.01–0.02, R 2 = 0.61–0.97) than in the Pasture area (RMSE = 0.01–0.03, R 2 = 0.61–0.90). Regarding the hydraulic parameters, for all models, the Pasture showed smaller saturated hydraulic conductivity and water content values of the mobile region than the Caatinga. This fact demonstrates the negative impact of compaction and change in natural vegetation in the Brazilian semi-arid. The dual-porosity model presented the best ﬁt to the data measured in the Pasture area. However, a single-porosity model could be considered representative of the Caatinga area. The results showed that Caatinga areas contribute to maintaining soil moisture and increasing the water storage in semi-arid regions. methodology, writing—original writing—review A.E.C.d.G.d.C.R.: methodology, investigation, writing—original draft preparation, writing—review and editing; A.C.D.A.: resources, supervision, project administration, funding. Highlights: Modeling soil water dynamics in Caatinga and Pasture areas; application of ﬁve different soil hydraulic models; Hydrus has a better ability to simulate water dynamics in the Caatinga area; the Ks and θ s parameters were lower in the Pasture area; the results demonstrated the need to test more than one hydraulic model for the study regions.


Introduction
The semi-arid region of northeastern Brazil has historically suffered from water shortage [1]. This area is characterized by the Caatinga biome, which has been suffering from anthropization processes. Particularly, in the Brazilian semi-arid region, the increase in the intensity of land use and the reduction of native vegetation cover causes the degradation of natural resources [2]. The removal of the caatinga, due to deforestation, for the implantation of pastures and agricultural crops can cause soil degradation. This fact contributes to changes in the dynamics of water transfer processes in the soil, plant and atmosphere system.
Pasture areas are generally exposed to intense animal trampling, causing soil compaction [3]. This contributes to a severe reduction in macroporosity, an increase in soil density, and reduction of water infiltration in the most superficial layers, negatively affecting the water dynamics in the hydrological cycle [4]. Therefore, one of the current challenges is the generation of data regarding the effect of land-use change on soil water dynamics. This information is scarce and necessary for conservationist management models of water resources.
In this work, a strategy for a better understanding of the impact resulting from landuse changes on hydrological processes was carried out through the continuous monitoring of soil moisture and other climatic variables in an area with natural vegetation of the and modeling the soil moisture dynamics with hydrological models in the Caatinga and Pasture areas is of fundamental importance to understand the dynamics of the hydrological processes in these environments. For this purpose, the National Observatory of Water and Carbon Dynamics in the Caatinga Biome (NOWCDCB) is a pioneering initiative that aims to gather a climate and eco-hydrological database in the semi-arid region of northeastern Brazil. Within this initiative, flux towers were installed to perform climatic and hydrodynamic monitoring in Caatinga and Pasture areas in the municipality of Serra Talhada (Pernambuco, Brazil). The AmeriFlux platform integrates all data monitored by these towers, which are identified as BR-CST (Caatinga) and BR-GST (Pasture) [31].
Thus, this work aims to analyze the effect of land use on the unsaturated soil hydraulic properties and hydrological process in the Caatinga (BR-CST) and Pasture (BR-GST) areas using the inverse method of Hydrus-1D. For this purpose, the following models were analyzed: single porosity (SP), dual porosity of Durner (DPD), dual porosity with mass transfer (DPTM), dual porosity with transfer of pressure (DPTP), and hysteresis in retention curve and hydraulic conductivity curve (HRCC). Figure 1 shows the towers of the National Observatory of Water and Carbon Dynamics in the Caatinga Biome (NOWCDCB). This project is located in Pasture (BR-GST) and Caatinga (BR-CST) areas, which are approximately 2.5 km apart. These areas are in the hydrographic basin of the Pajeú River in Serra Talhada, municipality of Pernambuco state, northeastern Brazil. The climate is classified according to Köppen as BSwh' (semi-arid climate), being characterized as warm and semi-arid, with summer rains concentrated between December and May (85%). The average annual precipitation is approximately 640 mm, and the average monthly air temperature ranges between 23.1 and 26.7 • C, with an annual average of 25.2 • C [32]. Chromic Luvisol soils, with depths varying between 0.40 and 0.90 m [33]. This Caatinga area is composed of shrubs (Mimosa hostilis, Mimosa verrucosa, and Croton sonderianus) and tree species (Anadenanthera macrocarpa, Spondias tuberosa, Caesalpinia pyramidalis, and Ziziphus joazeiro), and most of them have a height of approximately 8.0 m.

Study Area
The Pasture area (BR-GST), with 18 ha, is located in Lagoinha farm (7 • 56 50.53 S and 38 • 23 29.11 O, 450 m), also in the municipality of Serra Talhada. The area was planted with Urochloa mosambicensis, which is a C4 non-native grass, originating from the dry tropics of Africa [34]. The predominant soil is Chromic Luvisol [35], whose depths vary from 0.40 to 0.50 m [33]. Just as the Caatinga, during the rainy season, sheep and cattle graze this area. Then, during the dry season, it becomes highly exposed because of the lack of vegetation cover.

Monitoring of the Hydrological Data
The towers of the NOWCDCB project installed in the Pasture (BR-GST) and Caatinga (BR-CST) areas monitored soil moisture, precipitation, temperature, radiation, and wind velocity at a 2 m height. The data set was obtained in 2014, with data collection in the Caatinga starting in February.

Meteorological Data Monitoring
Historical series monitored by NOWCBC towers were used to measure atmospheric flows, water, energy, and CO 2 . The towers were equipped with sensors to measure the wind direction and speed (model 034, Campbell Scientific Inc., Logan, UT, USA), relative humidity, and air temperature (model HMP45C, Vaisala, Campbell Scientific Inc., Logan, UT, USA). Solar global radiation (Rg) was measured with a pyranometer (model LI-200X, LICOR Inc., Lincoln, NE, USA), and net radiation (Rn) was measured with a net radiometer (NRLITE, Kipp & Zonen, Delft, The Netherlands). Soil heat flux (G) was measured using two soil heat flux plates (model HFT3, REBS, Seattle, WA, USA) inserted 0.05 m below the soil surface. Two temperature sensors (model 108 L, Campbell Scientific Inc., Logan, UT, USA) were also located at 0.02 and 0.08 m below the soil surface to calculate the surface ground heat flux [36]. Total rainfall was measured with a tipping bucket rain gauge (model TE 525 WS-L, Texas Electronics, Dallas, TX, USA). Measurements from all sensors were recorded by a data logger (model CR10X, Campbell Scientific Inc., Logan, UT, USA) every 60 s. Mean/sum data were logged every 1800 s.
The daily reference evapotranspiration was calculated using the data acquired from the tower and applying the Penman-Monteith equation [37], which is expressed in Equation (1).
where ET o is the reference evapotranspiration (mm·day −1 ); ∆ is the slope of the vapor pressure curve at saturation versus air temperature (kPa • C −1 ); γ is the psychometric constant, 0.067 kPa • C −1 ; U 2 is the wind speed at 2 m height (m·s −1 ); e s is the vapor pressure at saturation (kPa); e a is the current vapor pressure (kPa); and T med is the average air temperature taken at 2 m ( • C).

Monitoring of Soil Water Content
The TDR CS616 probes (Campbell Scientific, Logan, Utah, USA), installed at 10, 20, 30, and 40 cm depths (Figure 2), were used to monitor soil moisture. The TDR sensors were calibrated based on the recommendations of Campbell's manual [8,38].  In this work, the cubic regression was used to adjust the soil moisture to pulse the sensors' time. Figure 3 shows the adjustment with the cubic model for both sites.

Soil Sample Properties
The soil physical properties were determined for both study areas. The disturbed soil samples were collected less than 5 m away from the moisture sensors, at layers of 5-15, 15-25, 25-35, and 35-45 cm. The particle size characterization of these samples was performed according to NBR 7181 [39]. This analysis was carried out in the Laboratory of Soil Physics at the Nuclear Energy Department (Departamento de Energia Nuclear-DEN) of the Federal University of Pernambuco (UFPE). Table 1 shows the soil particle size characterization of the soil and its textural classification. The bulk density (Ds) values were estimated according to [40] and [33] for the Caatinga and Pasture areas, respectively. All soil layers below the towers were homogeneous and classified as sandy loam. The determined particle density (ρ p ) was 2.42 and 2.50 g·cm −3 for the soil in the Caatinga and Pasture areas, respectively.  Table 2 presents, for each area, the average hydraulic parameters of the van Genuchten-Mualem [17] (VGM) equation obtained by pedotransfer functions of the Rosetta Lite program. These parameters were used as an initial estimation for the Hydrus-1D inverse method.

Modeling the Unsaturated Water Flow
The simulation of soil water dynamics was carried out with the Hydrus-1D. This program is capable of simulating water flow and solute transport in porous media with different saturation levels. The soil hydraulic models used were single porosity (SP), dual porosity of Durner (DPD), dual porosity with mass transfer (DPTM), dual porosity with transfer of pressure head (DPTP), and hysteresis in retention curve and hydraulic conductivity curve (HRCC).

Modeling with Single-Porosity (SP) Systems
The modeling with the SP model solves the Richards' equation (Equation (2)) using finite elements in the spatial discretization and finite differences in the temporal discretization [15].
∂θ ∂t where θ is the volumetric soil moisture (L 3 ·L −3 ), t is the time (T), z is the vertical coordinate (L), K(h) is the hydraulic conductivity, h is the matrix potential, and S(h) is the root water extraction term ((L 3 ·L −3 ) T −1 ). For the SP model, the van Genuchten-Mualem equation [17] was used to determine the relationships of θ(h) and K(h).
where θ r and θ s are residual and saturated water contents (L 3 ·L −3 ), respectively; K s is the saturated hydraulic conductivity (L·T −1 ); Se is the effective saturation; α (L −1 ), l, and n represent empirical shape parameters, with m = 1 − 1/n.

Single-Porosity (SP) Model with Hysteresis
The hysteresis phenomenon occurs when pressure head and soil water content do not have a single relationship, but different relationships in sorption and desorption. Sorption is the gradual increase in soil moisture by decreasing soil pressure. Desorption corresponds to the drying of the saturated soil, gradually increasing its suction. The Hydrus-1D code uses the empirical model introduced by [41] for hysteresis approaches. The HRCC model in the retention curve requires the main drying and main wetting curves. These two curves are described with Equation (6) using the parameter vectors θr d , θs d , θm d , α d , and n d and θr w , θs w , θm w , α w , and n w , respectively, where the subscripts d and w indicate wetting and drying, respectively (Equation (6)).
As specified in Equation (6), θ s and α are the only independent parameters describing hysteresis in the retention function.
According to the hysteresis model, the drying sweep curves are scaled from the main drying curve. In addition, the wetting sweep curves are determined from the main wetting curve. The scaling factors for the drying sweep curves can be obtained by adopting the hypothesis that the main drying curve is a reference curve in the scaling equation (Equation (7)).
making each scanning curve, θ(h), pass the point (θ ∆ , h ∆ ) characterizing the latest reversal from wetting to drying. Substituting this reversal point into the equation above, and assuming that θ r = θ r d , leads to Equation (8).
The scaling procedure results in a fictitious value of the parameter θs' for the drying scanning curve (Equation (9)).
in which the fictitious parameter θ r is now used. The scaling factor α θ for a given sweep curve can be obtained by replacing the inversion point (θ ∆ , h ∆ ) and the full saturation point (θ s , 0) into the equation above, and subtracting the two resulting equations to eliminate θ r to give Equation (10).

Dual-Porosity Model with Systems at Equilibrium
The modeling with dual-porosity systems at equilibrium solves Richards' equation (Equation (2)) with Durner's dual-porosity (DPD) model [20]. The DPD model divides the porous medium into two regions governed by the van Genuchten-Mualem relationship, which is a weighted average of the functions that describe the pore system (Equation (11)).
where w i is the percentage of each region, α i , n i , and m i represent empirical parameters, and Se i is the effective saturation in each region (i = 1 and 2).

Dual-Porosity Model with Mass Exchange
The dual-porosity mass transfer and pressure head transfer (DPTM and DPTP, respectively) models assume that water flow is restricted to interconnected pores (macropores) and that there is no flow in unconnected pores [42]. This conceptualization leads to dualporosity flow and transport models with two parts [43], dividing the liquid phase into mobile (θ m ) and immobile (θ im ) regions (Equation (13)).
The dual-porosity expression is based on a mixed formulation of Richards' equation [42], where the mobile and immobile water flow are described by a mass balance, as shown in Equations (14) and (15), respectively.
where S m and S im represent the drainage terms for mobile and immobile regions, respectively, and Γ w is the water transfer rate. The water between the fracture and matrix regions, Γ w , is assumed to be proportional to the difference in effective water contents of the two regions using the first-order rate equation (Equation (16)).
The rate of exchange of water between the fracture and matrix regions, Γ w , can also be assumed to be proportional to the difference in pressure heads between the two pore regions (Equation (17)).
where (T −1 ) and α ω (L −1 T −1 ) are first-order rate coefficients. Table 3 presents a summary of the necessary parameters to simulate the water flow in the soil's unsaturated zone for each model. Table 3. Soil hydraulic parameters used in each model SP, DPD, DPTM, DPTP, and HRCC.

. Atmospheric Boundary Conditions
Beer's law (Equations (18) and (19)) was used, separating the potential evaporation and transpiration fluxes, which partitions the solar radiation component of the energy budget via interception by the canopy [44].
where ET o , T p , and E p are potential evapotranspiration, transpiration, and evaporation fluxes (LT −1 ), respectively; LAI is the leaf area index; SCF is the soil cover fraction; and k is a constant governing the radiation extinction by the canopy as a function of sun angle, the distribution of plants, and the arrangement of leaves (between 0.5-0.75), which is 0.6 for the semi-arid Pernambuco [5]. Interception (I) was considered when the leaf area index (LAI) was entered, being defined by Equation (20) [45][46][47].
The LAI was calculated from the rectified equation generated by [48] in the municipalities of Serra Talhada and Petrolina (Pernambuco) (Equation (21)).
where NDVI is the normalized difference vegetation index, which, for this study, was determined by [49].

Root Water Uptake
The sink term, S, is defined as the volume of water extracted by the plant from a unit of representative volume of soil per unit of time due to plant water uptake, which is described by Equation (22) [50].
where the function α(h) is a prescribed dimensionless function of the soil water pressure head (0 ≤ α ≤ 1), and Sp the potential water uptake rate (T −1 ). For root water extraction, Hydrus-1D uses five pressure heads proposed by [50]. For both areas in this study, we used the pressure head established by [51] for the sandy loam soil, with woody vegetation for Caatinga and grassy for Pasture (Table 4).

Model Calibration and Validation
The calibration of the hydraulic parameters was performed by inverse modeling with the Hydrus-1D. The time frame used for calibration and validation was from the beginning of 2014 until 30 September and 01 October to 31 December, respectively. In the temporal discretization, the time unit was given in days, with the initial time considered as 0 (the day before the simulation) and the final time referring to the calibration and validation days mentioned above. Regarding the convergence criteria, Hydrus-1D solved Richards' nonlinear equation with a maximum of 10 iterations, allowing a tolerance of 0.05 cm and 0.001 for the potential and water content, respectively. Besides that, the soil profile was discretized using a mesh with 101 nodes.

Statistical Parameters Used to Evaluate Soil Hydraulic Models
Some statistical criteria helped the evaluation of the model's performance. For such analyses, we considered the root of the quadratic error of the mean (RMSE), which has an optimal value close to 0; Nash and Sutcliffe limited modeling efficiency (NSE), which limits the Nash coefficient between −1 (worst model performance) and 1 (best model performance); Kling-Gupta Efficiency coefficient (KGE), which has the same meaning as the NSE coefficient; determination coefficient (R 2 ), which indicates the degree of correlation between the independent variables and the dependent variable, ranging from 0 to 1, with an optimal value close to one; Pearson's correlation coefficient (r), which indicates the intensity of the linear association between the variables; the Willmott agreement index (d), which indicates the degree of precision of the equation, ranging from 0 to 1, with 1 being a perfect agreement; and the model's performance index (c), which allows analyzing the precision and accuracy of the results obtained through the product of the two coefficients, r and d, as proposed by [52]. fall. With a similar scenario, the BR-GST tower presented a wet season between February and May, representing 60% of the total annual rainfall. In addition, both towers also detected rainfall in November, representing approximately 20% of the annual rainfall. During rainy seasons, the Brazilian northeast is influenced by several meteorological systems, where the intertropical convergence zone (ZCIT) is the main one [54]. Meanwhile, the El Niño and Atlantic Dipole govern this region's dry season [55].   Based on the APAC data, it is observed that the wet season in the region begins in January and ends in April (65% of total annual rainfall). In 2014, the BR-CST tower identified a wet season between January and May, representing 68% of the total annual rainfall. With a similar scenario, the BR-GST tower presented a wet season between February and May, representing 60% of the total annual rainfall. In addition, both towers also detected rainfall in November, representing approximately 20% of the annual rainfall. During rainy seasons, the Brazilian northeast is influenced by several meteorological systems, where the intertropical convergence zone (ZCIT) is the main one [54]. Meanwhile, the El Niño and Atlantic Dipole govern this region's dry season [55]. Figure 5 presents the data series of volumetric moisture and daily precipitation observed in the towers. In the BR-CST tower ( Figure 5A), analyzing the soil moisture behavior as a function of rainfall, one verified an increase in soil moisture over all depths for rainfall greater than 20 mm. Meanwhile, lower rainfall only modified soil moisture in the first 10 cm of the profile.

Rainfall and Soil Moisture Dynamics
Moreover, one observed two large peaks of soil moisture, corresponding to rainfalls of 49 and 74 mm on 03/07 and 11/17. However, the water content did not reach the saturated water content. Due to the absence of rain, a slight moisture decay at a 10 cm depth was noticed. Using Hydrus-1D in the Caatinga area [12], also observed that the first 15 centimeters of soil retains most of the rainfall. He further noted an abrupt increase in soil surface moisture after a rainfall event. Figure 5 also indicates a similar relationship between precipitation and soil water content for both towers. Soil water content also changed after rainfall above 20 mm but with smaller peaks. In the dry period, soil water content remained constant at the highest depths. This effect occurs due to the low soil depth that limits deep percolation and favors runoff [7].  Figure 6 shows the inverse method adjustments for the soil hydraulic models SP, DPD, DPTM, DPTP, and HRCC at the BR-CST tower. In general, the models presented similar behavior. In the calibration step, the fits were better for the lower depths. At 30 and 40 cm depths, the model showed the same sensitivity to rain events as to 10 and 20 cm depths. However, the monitored soil moisture presented a higher variation at the surface.
The inverse method of one-dimensional water flow models finds difficulty in demonstrating the field spatial reality from data measured with TDR and tensiometers [56]. This type of failure occurs because the Hydrus-1D model is unable to account for the 3D horizontal water flow and lateral drainage that occur in the field.
In the validation step, the models reacted to the soil moisture peak of November and October, but they overestimated these peaks. In general, among the models, the DPTP simulated the lowest soil moisture values in rain events. Figure 7 shows the model adjustments for the Pasture area. They demonstrated sensitivity to rain events, but, overall, the DPTM model had the highest simulated soil moisture values. In the calibration, the models overestimated most rain events. Nevertheless, the values of the SP model were the closest to the measured data. In the validation step, the models also overestimated the soil moisture peaks. In the first peak, at a depth of 10 cm, only the HRCC model was able to simulate the moisture recession. This may prove that hysteresis must be accounted for. However, for the second and third peak, the DPTM model exhibited a better representation of moisture dynamics. At a depth of 20 cm, the SP model was the one that best represented soil moisture in the first rain event, while the HRCC model better represented the second and third rain events. Moreover, the DPD, DPTM, and DPTP models underestimated the dry season data at depths of 20 and 30 cm. However, at a depth of 40 cm, the DPTM model better represented the second rain event.  Table 5 shows the statistics of the data monitored by the BR-CST tower after being adjusted by the models. In the calibration step, the NSE coefficient indicated that only the DPTP model at a 30 cm depth was worse than using the average value of the data (NSE < 0). However, in the validation, at a depth of 20 cm, the NSE coefficient indicated that the SP and DPD models did not simulate the monitored data properly, and, at a depth of 40 cm, only the DPTM model achieved a satisfactory performance. The RMSE indicated errors ranging between 0.01 and 0.02 for the calibration and validation. The determination coefficient ranged from 0.61 to 0.88 in the calibration and from 0.86 to 0.97 in the validation. The r coefficient indicated very high (0.7-0.9) and almost perfect (>0.9) correlations for the models. Lastly, the performance index (c) indicated very good (0.76-0.85) rates for calibration and optimal (>0.85) values for validation. On the one hand, at a depth of 10 cm, the DPTM model presented the best fit for both calibration and validation. The only exception was the KGE in the validation, for which the DPTP model performance was better. On the other hand, at a depth of 20 and 30 cm, the SP and the DPD model showed better calibration performance. The NSE coefficient indicated the DPTP as the most representative model of the system in the validation step. However, at 40 cm depth, calibration suggests that the best model is the HRCC. However, in the validation, the best model was the DPTM.
Based on the statistics average associated with each depth, the best model, according to the NSE and KGE coefficients, was the DPTM. However, based on R 2 , r, and c, the best-adjusted model was the SP, indicating that this model can be used in Caatinga areas. Even though it is simpler than the other models, it fits as well as them or even better. Table 6 shows the statistics of the data adjusted to the models analyzed in the BR-GST tower. In the calibration, at a depth of 10 cm, the NSE coefficient indicated that only the SP and DPTP models represent the measured data properly. Conversely, R 2 indicated that these models had the worst performance. In the validation, only DPTP (at 20 cm depth) and DPTM (at 40 cm depth) models were not satisfactory according to the NSE coefficients. The R 2 for these models ranged from 0.61 to 0.90. Besides that, RMSE ranged from 0.01 to 0.03 in the calibration and from 0.01 to 0.02 in the validation. Overall, the correlations (r) were moderate (0.4-0.69) and strong (0.7-0.89) in both the calibration and the validation. Moreover, the performance index (c) was estimated for the models as sound (0.66-0.75) and ideal (0.76-0.85). An exception was the SP model, which presented poorly (0.51-0.60) and sound (0.66-0.75) adjustments [52].
The statistical indices indicate that, at a 10 cm depth, the DPD and HRCC models had the best performance for calibration and validation, respectively, while at depths of 20 and 40 cm, the ideal model was the DPTP. For the 30 cm depth, calibration and validation were better analyzed by the DPTP and DPTM models, respectively.
Based on the average of the KGE, R 2 , r, and c coefficients for each depth, the best adjustments correspond to the DPTP model. However, considering the NSE and RMSE statistics, the best fit was the SP model. This fact shows that, contrary to the Caatinga, the use of the SP model in Pasture areas is not indicated because it could lead to inconsistent estimations of soil moisture. This occurs due to a better structured and less dense soil in the Caatinga area, favoring the water redistribution process in the soil matrix. In contrast, in the Pasture area, the soil is poorly structured and undergoes greater compaction and more intense weathering processes. This leads to the formation of cracks that create a preferential path for water in the soil. The consequence is a change in the paradigm of soil moisture dynamics, making the SP model inappropriate and promoting the use of dual porosity approaches. Table 7 shows the Hydraulic parameters of the models adjusted by the inverse method in Hydrus-1D. The saturated hydraulic conductivity (Ks) ranged from 107.1 to 445.9 cm·d −1 for the Caatinga and 23 to 206 cm·d −1 for the Pasture, depending on the chosen model. Previous study in the Brazilian semi-arid has implemented Hydrus-1D inverse method for two experimental areas [8]. The authors obtained Ks values of 234.2 and 243.7 cm·d −1 , which are within the range observed in the field for the BR-CST tower. A research in the semiarid region of Tunisia, also using the inverse method, obtained Ks values between 28 and 200 cm·d −1 for different soil depths. The authors also observed an inverse correlation between soil compactions and Ks [57].
In the BR-CST tower, the HRCC model had the highest and lowest Ks values in its drying and wetting curves, respectively, where Ks d . Ks w−1 was 4.16. In the BR-GST tower, the HRCC d model had the highest Ks value, while the DPTM model had the lowest Ks value. The ratio between Ks d and Ks w was lower in the Pasture (about 1.36) than in the Caatinga (4.16), indicating higher hysteresis in the Caatinga hydraulic conductivity curve. In general, the BR-CST tower had the highest values of Ks, except for the HRCC model (in the wetting curve). This fact occurs due to the compaction and the land-use change suffered by the soil of the Pasture area. Furthermore, [5] also detected a lower Ks in the non-vegetated Caatinga area (180 cm·d −1 ) than in the vegetated Caatinga area (200 cm·d −1 ).
In an analysis of the hydrodynamic properties of Caatinga and Pasture soil through the Beerkan methodology, the Caatinga Ks was almost twice that of the Pasture [4]. The authors claimed that intense animal trampling caused compaction and reduced soil macroporosity, which increased soil density and reduced water infiltration. A comparative analysis of areas of preserved Caatinga, Caatinga under regeneration, cultivated land, and degraded Pasture was carried out in the municipality of Serra Talhada, state of Pernambuco [58]. They observed that preserved Caatinga had the highest water infiltration rate in the soil and the lowest water and soil loss. Meanwhile, [59] compared the change in physical water indicators over an area of 40 × 40 m after 21 days of intensive grazing. Through this study, they confirmed that soil compaction due to intense cattle trampling causes a significant reduction of Ks by almost 50%.
The θs m was higher in the Caatinga area. This happens due to the extensive use of Pasture soil and the cumulative effect of cattle trampling. These factors influence the increase in density and resistance to penetration, leading to a decrease in soil porosity in the first centimeters of the soil [60]. In contrast, θs im (DPTM and DPTP models) was higher in the area of the BR-GST tower. The θ sim corresponded to 31% (DPTM) and 44% (DPTP) of total saturated humidity in BR-CST. In BR-GST, θs im was equal to 57% (DPTM) and 68% (DPTP) of θs. This fact is associated with a higher fines percentage present in the BR-GST, which is about 1.5 times higher than in the BR-CST. The higher percentage of θs im than θs m is associated with a higher water retention capacity of the porous soil matrix. Furthermore, the introduction of the dual-porosity model improves the accuracy in describing the hydraulic behavior of the porous medium [29].
This performance is also evident in the w 2 parameter (DPD model), which is six times higher in BR-GST than in BR-CST, therefore indicating that the Pasture area holds two porous regions of more significant evidence. The θs (θs im + θs m for dual-porosity models) ranged from 0.3 (HRCC w ) to 0.458 cm 3 ·cm −3 (DPTM) in the BR-CST and between 0.2 (SP, DPD, and HRCC w ) and 0.35 cm 3 ·cm −3 (DPTM) in the BR-GST. A study in the semi-arid region of China used the DPTP model at depths between 0 and 50 cm for sandy loam soil, determining θs values within the range observed in the Caatinga (θ s = 0.438 cm 3 ·cm −3 ), with θs m = 0.43 cm 3 ·cm −3 and θ sim = 0.008 cm 3 ·cm −3 [61]. Using the SP model for sandy loam soils, another study obtained θ s values between 0.342 and 0.369 cm 3 ·cm −3 [8].
The n parameter ranged from 1.374 (HRCC d ) to 1.715 (DPD 2 ) in the BR-CST and from 1.321 (DPD 2 ) to 2.35 (HRCC d ) in the BR-GST. This parameter is associated with soil texture, being smaller for soils with finer textures. However, for the most part, the lowest values for this parameter were obtained in the BR-CST tower (except for the DPD model n 2 ), which has the lowest percentage of fines. Using the HYPROP evaporation method in conjunction with Hydrus-1D to determine hydraulic parameters, [62] applied the SP and DPD models in soils with different texture and pore sizes. The authors obtained values for the n parameter of 1.32, 1.83, 2.59, and 10 for soils with fine, medium, coarse, and very coarse textures, respectively, when using the SP model. For the DPTP model, the n im was 1.44 and 1.67 for BR-CST and BR-GST, respectively [62]. In another study n im was 1.402 (0-50 cm) and 1.278 (50-100 cm) [61]. These values differ from those obtained by [28], who determined n im equal to 6, using the DPTP model with the Hydrus-2D inverse method.
The α parameter ranged from 0.004 (DPTMim) to 0.066 cm −1 (DPD2) in the BR-CST and from 0.001 (SP) to 0.08 cm −1 (HRCC W ) in the BR-GST. Within this range are the values obtained by [5], who determined α equal to 0.028 cm −1 . Using the DPTM model [63] defined α equal to 0.023 cm −1 . However, using the DPD model, [63] measured α values equal to 0.2 and 0.002 cm −1 , while [60] found α values of 0.013 cm −1 and α im of 1.007 cm −1 using the DPTP model. In the HRCC model, the α w /α d ratio was 2.14 and 1.33 for the BR-CST and BR-GST towers, respectively. The obtained values are within the range determined by [64], who verified an average value for α w /α d of 2.24 ± 1.25. The authors further stated that the higher the α w /α d value, the higher the level of soil cohesion, indicating greater particle cohesion in the soil of the BR-CST tower than in the BR-GST, even with a lower percentage of fines. This fact may occur due to the higher level of organic matter in Caatinga vegetated soil. A study evaluated the impact of different sizes of soil aggregates on the hysteresis of their water retention characteristics [65]. According to the author, the size of the aggregates has a significant impact on the hysteresis effect of water retention in soils, and the negligence of the hysteresis effect may cause inaccuracy in soil water content estimation [64].
The BR-CST tower also exhibited lower values of the pore connectivity parameter (l). Besides that, for this area of Caatinga, the l parameter was negative for the SP, DPD, and HRCC models. Negative l parameters may occur because these models do not consider dual porosity with pressure or moisture exchange between the mobile and immobile regions of the soil. Moreover, according to [65,66], negative values can be understood as correction factors that cause a more gradual drop in the unsaturated hydraulic conductivity. Materials with fine pore distribution have commonly negative values of the l parameter, or l can even be interpreted as an empirical parameter that serves to restrict relative hydraulic conductivity [67]. Figure 8 shows accumulated infiltration (sum (Infil)), accumulated evaporation (sum (Evap)), accumulated root water uptake (sum (vRoot)), accumulated drainage flow (sum (vBot)), and variation of the stored volume (Volume) for soil hydraulic models (SP, DPD, DPTM, DPTP, and HRCC) in the Caatinga area. According to these results, the drainage flow is much smaller than the amount of total infiltrated water, favoring the evapotranspiration process. When comparing the models, evaporation ( Figure 8C) and root water uptake ( Figure 8D) showed different behaviors. DPTM had the highest evaporation and the lowest root water uptake, while the behavior of the DPTP model was precisely the opposite. The results of the other models for these variables were slightly different. Regarding the accumulated infiltration ( Figure 8A), the models had similar behavior. By fitting the model to water retention and water pressure head parameters at different depths, we obtained similar infiltration rates at the surface. Figure 9 shows the accumulated water flows simulated by soil hydraulic models (SP, DPD, DPTM, DPTP, and HRCC) in the Pasture area. Unlike the Caatinga, the lower drainage flow presented high magnitudes in the HRCC model ( Figure 9A), medium magnitudes in the SP and DPTP models, and small magnitudes in the DPTM and DPD models. The HRCC model was the one that estimated the lowest accumulated root water uptake, whereas the DPTM model obtained the highest estimation of this variable. Regarding the stored water height ( Figure 9E), the models presented similar sensitivity to rain events, with changes in the magnitudes of the stored water volume peaks. The DPTM model was the one with the highest volume peaks, followed by DPD, DPTP, HRCC, and SP, respectively. In the recession phase, the HRCC model had the most significant stored volume decay, while the SP model kept the volume approximately constant. According to [68], the process of soil water drainage occurs faster when hysteresis is considered.

Actual Evapotranspiration for the Caatinga and Pasture Areas
In order to ascertain the models' ability to estimate the actual evapotranspiration (Eta), the estimate of evaporation and transpiration provided by the models was compared with the measurement by the eddy covariance (EC) method. Figure 10 shows the daily ETa simulated by the models and measured by the EC in the Caatinga (BR-CST tower). The estimated models did not present significant divergences among them. It is possible to observe that the models present a good fit from month 08, responding to the Eta peaks and remaining constant in the low ETa period. However, in the period prior to month 08, the models presented oscillations where the ETa measured by the EC remained constant ( Figure 10). Table 8 shows the statistics related to the models' ability to simulate ETa in the Caatinga area. The best statistics were from the SP model. The NSE and KGE indices ranged from 0.04 to 0.09 and from 0.35 to 0.37, respectively, as they are greater than zero, indicating that the models were adequate. The RMSE ranged from 0.11 to 0.12 cm·d −1 and the R 2 ranged from 0.51 to 0.56. The correlation coefficients (r) were considered as strong (0.7-0.89) and with good performance indices (0.66-0.75).  Figure 11 shows the daily ETa simulated by the models and measured by the EC in the Pasture area. The models accurately simulated the dry period; however, they overestimated the ETa peaks. The ETa peaks estimated by the models were remarkably similar, except for the HRCC model. HRCC underestimated ETa in most parts of the simulated period. Figure 11. Comparison between the ETa estimated by the models (A) and the measurement by eddy covariance (B) in the Pasture area. Table 9 shows the statistics related to the models' ability to simulate ETa in the Pasture area. All models presented negative NSE and KGE coefficients, with the exception of the HRCC model, which obtained a KGE value equal to 0.26. The RMSE varied between 0.09 (HRCC) and 0.16 (DPD) cm·d −1 and the R 2 between 0.28 (HRCC) and 0.5 (DPTM). The correlation coefficients (r) were considered moderate (0.4-0.69) for the SP, DPD, DPTP, and HRCC models and with a strong correlation only for the HRCC model. The c index indicated that the HRCC, SP, DPD, and DPTP models had poor performance and the DPTM model had moderate performance. Thus, it is noted that there is a great impact on the physical paradigm of the soil due to the anthropogenic processes presented in the Pasture, leading to the need to separate the vadose zone into a mobile and immobile region. In addition, in tropical regions, there is a more accelerated pedogenetic or soil formation process, particularly in semi-arid regions. Higher temperatures and more intense water action in the short rainy periods cause phenomena such as retraction and expandability that influence the processes of water storage. This can be confirmed by the measured and modeled Eta. The ETa in the Caatinga was higher than in the Pasture area. These values may be associated with the depth of the roots and the interception of rainfall by vegetation.
The models in both areas responded to ETa peaks and remained constant during the low ETa period. The Hydrus-1D model showed a better simulation of ETa for the Caatinga area than the Pasture area.
As seen in the results of Tables 8 and 9, the SP model can be used to characterize the Caatinga area. However, the results indicate that an approach using dual-porosity models (DPTM and DPTP) is more appropriate in anthropized areas. However, the improvements achieved via the DPTM and DPTP methods shown in Table 9 did not demonstrate significance; therefore, the hypothesis of the single-porosity model could not be discarded.

Comparison of the Extreme Rainfall Event That Occurred in November by the SP Model in the Caatinga and the DPTP Model in the Pasture
SP and DPTP models provided the best adjustments for the Caatinga and Pasture areas, respectively. Therefore, the extreme precipitation event that took place between 16 and 18 November 2014 was simulated for both models on an hourly scale to compare the main flows in these areas. Thus, Figure 12 shows the fluxes of runoff, accumulated bottom discharge (sum (vBot)), accumulated infiltration (sum (Infil)), and soil water storage (Volume) for both areas. One noted a greater runoff ( Figure 12A) in the Pasture than in the Caatinga. In contrast, accumulated infiltration was greater in the Caatinga ( Figure 12B) than in the Pasture. A possible reason for this would be the well-structured soil in Caatinga, leading to a higher total porosity and water infiltration rate than a poorly aggregated and compacted Pasture soil [69]. The accumulated bottom discharge ( Figure 12C) behaves similarly in the Caatinga and Pasture areas. However, at 10 h, there is a stronger increase in the accumulated bottom discharge in the Caatinga area. Soil water storage ( Figure 12D) is sounder in the Caatinga area than in the Pasture area.
In general, this work presented a detailed study of the processes of water transfer in the soil-plant-atmosphere system in the Brazilian semi-arid region. The physical models of the soil took into account the dual porosity of the porous medium and the hysteresis in soil water properties in Caatinga (with shrub cover) and Pasture (with grass cover) areas.
It was possible to observe that in the Caatinga area, due to the soil being better structured, there is a redistribution of water in a smoother and more uniform way across the soil. Capillary flow is also noted and less influenced by the preferential flow. This led to a little divergence between the physical models analyzed; therefore, it is possible to use the single-porosity model to represent the soil in the Caatinga area. On the other hand, in the Pasture, there is a higher preferential runoff effect on cracks in the soil and more impact of water fractionation into mobile and immobile water contents. This fact makes it essential to analyze this porous medium in terms of dual porosity. This occurs due to the strong weathering process that the soil surface undergoes due to the absence of vegetation, higher insolation rates, and also animal trampling, which led to higher runoff rates and lower infiltration rates.
One can observe that the processes of redistribution and storage of water from the deepest layers occur faster and in greater volume in the Caatinga than in the Pasture ( Figure 12B). Furthermore, one can conclude that the preservation of Caatinga areas contributes to the maintenance of soil moisture ( Figure 12D) and, therefore, an increase in water storage in semi-arid regions.
In addition, future research can observe the effect of the estimation of evapotranspiration methods on the quality of the adjustments. In addition, the effect of heterogeneity using information from several layers can also be considered. Efforts to improve knowledge about soil structure in the Caatinga and Pasture areas have been carried out by the National Observatory of Water and Carbon Dynamics in the Caatinga Biome (NOWCDCB) with the application of ground penetrating radar and x-ray microtomography techniques. These evaluations will allow the improvement of simulation processes with more accurate information.

Conclusions
The present work evaluated the capacity of SP, DPD, DPTM, DPTP, and HRCC models to simulate the water dynamics in the Caatinga and Pasture areas. Therefore, the Hydrus-1D program was used to simulate water flows in the soil's vadose zone. Besides that, one also estimated the soil hydraulic properties in the Caatinga and Pasture areas using the inverse method of Hydrus-1D.
In general, the hydraulic models showed better adjustments in the Caatinga than in the Pasture. Thus, all the soil hydraulic models achieved a better adjustment for soils with greater moisture variability. Although all the adjusted models were able to respond to the observed moisture peaks, most of them overestimated these peaks. The models also simulated water content peaks in the dry season, which were not observed in the field. In general, Caatinga areas contributed to the maintenance of soil moisture and increase in water storage in semi-arid regions. The results also indicated a lower efficiency of the inverse method in periods with moderate rainfall.
The SP model presented the best adjustments for the Caatinga area, and the DPTP model presented the best adjustment for the Pasture area. However, the SP and DPTP models showed no significant difference, and the SP model can be used for both areas. That is, the implementation of the preferential flow hypothesis did not significantly improve the simulation.
The Pasture area presented smaller Ks and θs m than the Caatinga area in all the models, demonstrating the negative impact of compaction due to land-use change (conversion of Caatinga to Pasture areas) in the Brazilian semi-arid. On the other hand, θs im was higher in the Pasture area than in the Caatinga. This result shows that for soils with a sound fine percentage and/or greater compaction, it is essential to consider dual-porosity models.
Models in both areas respond to ETa peaks and remain constant during the low ETa period. The ETa was better simulated with the Hydrus-1D model for the Caatinga area than for the Pasture area. In addition, it was demonstrated that the HRCC model has great relevance in the estimation of the ETa.