Next Article in Journal
Influence of Pre-Hydrolysis on Sewage Treatment in an Up-Flow Anaerobic Sludge BLANKET (UASB) Reactor: A Review
Next Article in Special Issue
Effects of Biochar Application and Irrigation Methods on Soil Temperature in Farmland
Previous Article in Journal
Numerical Simulation of Phosphorus Release with Sediment Suspension under Hydrodynamic Condition in Mochou Lake, China
Previous Article in Special Issue
Experimental and Simulation Investigation on the Kinetic Energy Dissipation Rate of a Fixed Spray-Plate Sprinkler
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Modeling Two-Dimensional Infiltration with Constant and Time-Variable Water Depth

1
Programa de Maestría y Doctorado en Ingeniería, Universidad Nacional Autónoma de Mexico, Ciudad de México 04510, México
2
Comisión Nacional de Agua. Avenida Insurgentes Sur 2416, Copilco El Bajo, Ciudad de México C.P. 04340, México
3
Instituto Mexicano de Tecnología del Agua. Coordinación de Riego y Drenaje. Paseo Cuauhnáhuac 62550 Col. Progreso Jiutepec, Morelos C.P.8532, México
*
Author to whom correspondence should be addressed.
Water 2019, 11(2), 371; https://doi.org/10.3390/w11020371
Submission received: 11 January 2019 / Revised: 10 February 2019 / Accepted: 19 February 2019 / Published: 21 February 2019
(This article belongs to the Special Issue Precision Agriculture and Irrigation)

Abstract

:
Water infiltration is simulated by obtaining the time infiltrated depth evolution and humidity profiles with the numerical solution of the two-dimensional Richards’ equation. The contact time hypothesis is accepted in this study and used to apply a unique form on time of the water depth evolution in the solution domain (furrow), as boundary condition. The specific form of such evolution in time was obtained from results reported in the literature based on the internal numerical full coupling of the Saint-Venant and Richards’ equations in border irrigation. Moreover, the equivalent hydraulic area between the border and the furrow was achieved by scaling the values of water depth. The analysis was made for three contrasting soil textures, and the comparison was done by computing the root mean square error (RMSE) indicator. The comparison was performed from the selection of five finite element meshes with different densities to discretize the solution domain of the two-dimensional Richards’ equation, combined with several time steps. Finally, a comparison was made between infiltrated depth evolution calculated with a constant water depth in the furrow to the one proposed in this work, finding important differences between both approaches. To expand the scope of this study and for a fuller exploration of the subject, the results were compared with results obtained by applying the HYDRUS-2D software. The results confirm that it is important to consider an internal full coupling of the Saint-Venant and Richards´ equations to improve furrow irrigation simulations.

Graphical Abstract

1. Introduction

Richards’ equation [1] is used in most hydrological models to describe groundwater flow in porous media. Richards’ equation is derived from the combination of a continuity equation and velocity ranges obtained with Darcy’s law; the three-dimensional form without water extraction by plants is:
C ( ψ ) ψ t = · [ K ( ψ ) ψ ] dK d ψ ψ z
where ψ is the water pressure potential into soil, expressed as the height of an equivalent water column (L) (positive in saturated zones and negative in unsaturated zones of soil); C ( ψ ) = d θ / d ψ is the humidity specific capacity of soil; θ = θ ( ψ ) is the water volume per unit volume of soil or water content; (L3·L−3) is a ψ function, known as the humidity characteristic curve or water retention curve; K = K ( ψ ) is hydraulic conductivity (L·T−1), in a saturated soil as a function of pressure potential; gravitational potential is assigned to spatial coordinate z positively oriented down (L); = ( / x , / y , / z ) is the gradient operator; x and y are spatial coordinates (L), and t is time (T).
In the literature, studies on infiltration modeling in surface irrigation can be found; for example, Liu et al. [2] proposed a coupled model in which surface water flow and solute transport are described using the zero-inertia equation and the average cross-sectional convection–dispersion equation, respectively, while the two-dimensional Richards’ equation and the convection–dispersion equation are used to simulate water flow and solute transport in soils, respectively. Dong et al. [3] developed a new numerical methodology based on the physical-based model of surface irrigation and the Monte Carlo simulation method to improve the modeling accuracy of surface irrigation performance at a field scale. Border irrigation was simulated by coupling the Saint-Venant equations and a one-dimensional Richards’ equation. Duncan et al. [4] presented a mathematical model that describes the movement of water and solutes in a ridge and furrow geometry. They focus on the effects of two physical processes: root water uptake and pond formation in the furrows. Special attention was paid to the physical description of ponding as an effect of transient rain events. Fan et al. [5] used vertical line source irrigation as a water-saving irrigation method in order to enhance water and nutrient delivery to the root zone. Therefore, reduction on soil evaporation and improvement on water and nutrient efficiency was achieved. To identify influencing factors, they performed computer simulations using the HYDRUS-2D software. Brunetti et al. [6] studied a hybrid Finite Volume–Finite Element (FV–FE) model that describes the coupled surface–subsurface flow processes occurring during furrow irrigation and fertigation. Furman [7] discussed some research gaps, led by the need to include vertical momentum transfer and to expand the use of fully coupled models toward surface irrigation applications. Banti et al. [8] developed a model for the simulation of furrow irrigation advance based on the Saint-Venant equations for one-dimensional surface flow and the two-dimensional Richards’ equation for porous media flow. Bautista et al. [9] proposed a methodology for estimating furrow infiltration under time-variable ponding depth. The methodology approximated the solution to the two-dimensional Richards’ equation, and was a modification of a procedure that was originally proposed for computing infiltration under constant ponding depth. Bautista et al. [10] expanded the analysis of a proposed furrow infiltration formulation based on an approximate solution to the two-dimensional Richards’ equation. The approach calculates two-dimensional infiltration flux as the sum of one-dimensional infiltration and a second term labeled the edge effect. Dong et al. [11] presented a hybrid coupling strategy by distinguishing the mass conservation and momentum conservation of surface flow and subsurface flow systems to reduce the iteration times of the momentum conservation equation and increase the simulation accuracy of coupled models. They proposed a coupling strategy, where the mass conservation component of surface flow model and subsurface flow model are iteratively coupled at the same time step to obtain the convergence value of surface flow depth, and then the momentum conservation component of surface flow model was externally coupled based on the convergence value of both the surface flow depth and infiltration rate to update the surface flow velocity. In the hybrid coupling procedure, infiltration was used as the common internal boundary condition. Soroush et al. [12] developed a furrow irrigation model based on the slow-change/slow-flow routing equation, which is an approximate reduced form of the Saint-Venant equations to a single equation with a single variable, the upstream volume of water. Downstream-propagating disturbances show that the rate of change of upstream inflow is small, with no limit on the Froude number, so it can be used for all slopes and with all common end conditions. To calculate resistance to flow, a composite model in terms of almost any boundary roughness was proposed. Infiltration is assumed to follow the Kostiakov formula. Wöhling et al. [13] presented a process based on a seasonal furrow irrigation model which describes the interacting one-dimensional surface–two-dimensional subsurface flow and crop growth during a whole growing period. The irrigation advance model was based on an analytical solution of the zero-inertia surface flow equations and was iteratively coupled with the two-dimensional subsurface flow model HYDRUS-2D. Wöhling et al. [14] showed a proposed furrow advance phase model (FAPS) which further develops the concepts of a previous study. An analytical zero-inertia surface flow model was iteratively coupled with the two-dimensional (2D) subsurface water transport model HYDRUS-2D. Wöhling et al. [15] proposed a model that overcomes the restrictions of traditional furrow irrigation models by rigorously describing the subsurface flow at computational nodes using the physically based two-dimensional (2D) model HYDRUS-2D. Surface flow is portrayed by an analytical zero-inertia model. Saucedo et al. [16] coupled one-dimensional Saint-Venant equations and two-dimensional Richards’ equation in border irrigation modeling.
Most of the above cited studies apply simplified models of either surface or subsurface water flow. The aim for this study is to note the disadvantage of using water depth average in two-dimensional infiltration simulations, and moreover to compare the infiltrated depth evolution in two-dimensional infiltration, in two cases: (a) considering a constant water depth in the furrow, and (b) considering a time–water depth evolution adapted as reported in literature in relation to a numerical full coupling of Saint-Venant equations and Richards’ equation in border irrigation. Otherwise, the importance of possible differences are analyzed comparing them to those calculated only by the spatial-temporal discretization of the two-dimensional Richards’ equation. This study focuses on phenomena taking place in the furrow cross-section, not along the furrow.

2. Materials and Methods

Irrigation is a three-dimensional phenomenon, therefore, water flux into the soil should be described as shown in Equation (1). This study accepts the contact time hypothesis, which implies that water has a predominantly vertical movement in the soil, i.e., the infiltration law is unique along the furrow. For the above it is possible to use the two-dimensional Richards’ equation form:
C ( ψ ) ψ t   =   x [ K ( ψ ) ψ x ]   +   z [ K ( ψ ) ( ψ z 1 ) ]
This equation can be solved on a solution domain as shown in Figure 1.
The function to define solution domain geometry on its upper border is [17]:
z   =   P s 2 [ 1 cos ( π 2 x L ) ]
Spatial pressure distribution as the initial condition necessary in obtaining two-dimensional Richards’ equation solution is:
ψ   =   ψ o ( x , z )
Boundary conditions are: EF ¯ Dirichlet frontier, prescribed potential according to the water depth in furrow center, AB ¯ , CD ¯ , DE ¯ and FA ¯ Neumann frontiers with zero flux, and BC ¯ frontier with unitary gradient:
ψ   =   h     ( z L / 2   z ) ,   x   EF ¯ ,   z AB ¯ ,   t > 0
K ( ψ ) ( ψ z ) z   =   0 ,   x   FA ¯ ,   z AB ¯ ,   t   >   0
K ( ψ ) ( ψ     z ) z   =   0 ,   x   DE ¯ ,   z CD ¯ ,   t   >   0
K ( ψ ) ( ψ ) x   =   0 ,   x = 0 ,   z AB ¯ , t > 0
K ( ψ ) ( ψ ) x = 0 ,   x = L ,   z CD ¯ ,   t > 0
( ψ     z ) z   =   1 ,   x   BC ¯ ,   z = P ,   t > 0
Richards’ equation solution is required to represent the hydrodynamic properties of soil, expressing the pressure potential ( ψ ) as a function of water content ( θ ) and hydraulic conductivity ( K ) as a function of θ .
Simulation of water flow requires a soil hydraulic characterization, as noted by Fuentes et al. [18], with a combination of characteristics from Fujita [19] and Parlange et al. [20] used in theoretical studies. Experimental studies may require the combination of the retention curve proposed by van Genuchten [21], with the restriction by Burdine [22], and the hydraulic conductivity curve proposed by Brooks and Corey [23] as they meet the integral properties of infiltration, and facilitate the identification of their parameters.
The retention curve is [21]:
θ ( ψ ) θ r θ s θ r =   [ 1 +   ( ψ ψ d ) n ] m
where ψ d is the water pressure characteristic value into the soil; m and n are empirical parameters related by Burdine’s restriction [22]: m = 1 – 2/n, with 0 < m < 1 and n > 2 ; θ s is the water content at effective saturation of soil, and θ r is the residual water content.
The hydraulic conductivity curve is [23]:
K ( θ ) = K s ( θ θ r θ s θ r ) η
where η is an empirically positive parameter. Cumulative infiltration per unit furrow length is:
A I ( x , t ) = P m q I ( x ,   y ,   z ,   t ) dP m
where q I ( x ,   y ,   z ,   t ) is the infiltration rate per unit at the furrow’s surface and P m is the wet perimeter of the furrow.

2.1. Numerical Solution with Finite Element Method

The numerical two-dimensional Richards’ equation solution used is made of finite elements in space, and implicit finite differences in time. The solution of the algebraic equations system that results from the application of the finite element method is effected using the preconditioned conjugate gradient method [24], which has been adapted to use a free from zeros computational storage matrix.
Multiplying Equation (2) by a weight function ( ν ) , the weak expression form of the two-dimensional Richards’ equation is obtained through integration using Green’s theorem on solution domain ( R ) limited by the border ( ) :
R C ( ψ ) ψ t ν dR + R K ( ψ ) [ ψ z ν z + ψ x ν x ] dR = R K ( ψ ) ν z dR + Γ q ν d Γ
where Γ is a fraction of under the Neumann’s condition with prescribed flux q . The Equation (2) solution is assumed as a linear combination of base functions φ i ( x , z ) defined in relation to the Kronecker’s delta function and applied to each particular node:
ψ ˆ n ( x , z , t ) = j = 1 n a j ( t ) φ j ( x , z )
where a j (t) are coefficients to be determined and n is the number of nodes where the finite element solution is obtained. It is replaced in the first weak form of Richards’ equation (Equation (14)), considering the following: (i) Weight functions are considered equal to the base functions ( φ ) corresponding to the interior nodes. (ii) A linear variation of the hydrodynamic characteristics is assumed on each element expressing it through the functions of form, i.e., C ˆ = φ g C g and K ˆ = φ g K g . (iii) A concentrated mass system is used in order to obtain a diagonal matrix and to improve the stability of the scheme [25,26]. Temporal derivative is approximated by an implicit scheme in finite differences and is obtained:
[ M k + 1 Δ t + K k + 1 ] a k + 1 = [ M k + 1 Δ t ] a k + B k + 1 + Q k + 1
where the matrices are calculated as indicated below when using linear base functions:
M kj = j = 1 n [ C g R φ ¯ g φ ¯ j φ ¯ k dR ] = e C j Δ 3
K kj = j = 1 n [ K g R φ g ( φ j x φ k x + φ j z φ k z dR ) ] = e K ¯ 4 Δ ( m j m k + p j p k )
B k = K g R φ g φ k z dR = e K ¯ 2 p k
Q k = Γ q φ ¯ k d Γ = e q L j 2
In previous equations, φ ¯ represents the so-called concentrated mass functions defined as unitary functions on a barycentric region corresponding to a specific node and zero in the rest of the domain. Δ is the element area, k ¯ is the element conductivity calculated as the arithmetic average of the conductivities obtained in each of its corners (as a consequence of the form adopted for the base functions), C j is the specific capacity estimated at the node j, L j is the border length corresponding to each node under Neumann’s condition, m and p are geometric factors defined according to the base functions: m i = z j z k and p i = x j x k where the subscripts i, j and k correspond to the corners of the triangular element and run on their three possible sequenced permutations.

2.2. Solution Domain Characteristics

The cross section of the furrow had the following dimensions: a width of 100 cm, and a depth of 150 cm; in the upper zone of the domain, the shape of the furrow obtained with Equation (3) was introduced.
To solve two-dimensional Richards’ equation numerically five finite element meshes with different densities were built (Figure 2). Characteristics of each mesh are shown in Table 1.
Time steps discretization proposed are Δ t = 1.0 , Δ t = 5.0 , Δ t = 10.0 , Δ t = 30.0 and Δ t = 60.0 s.
Selection of three different soil textures was made to cover a spectrum that involved different conditions. The contrasting chosen soil textures ensure that other soil textures have intermediate behaviors among the analyzed ones. Hydrodynamic characteristics of the three contrasting soil types are shown in Table 2. Residual saturation was assumed equal to zero in accordance with Fuentes et al. [18]. Total soil porosity was assumed as saturation water content ( ϕ ) , determined based on values provided by Rawls et al. [27]. For determined values of m and η parameters, a granulometric curve was built for each soil, based in sand, silt and clay percentages of triangle textures of USDA (U.S. Department of Agriculture), following Fuentes’ procedure [18]. Water content necessary to assign Richards’ equation initial condition was determined considering usable humidity of each soil type, assuming 50% humidity before irrigation. Usable humidity was determined by subtracting water content corresponding to field capacity and permanent wilting point.

2.3. Border Condition in the Furrow

To solve two-dimensional Richards’ equation numerically, it is necessary to know time–water depth evolution in furrow. Accepting the contact time hypothesis [28] allows proposing a unique form of time–water depth evolution in furrow throughout its length as a first approximation. The internal numerical coupling of the Saint-Venant and Richards’ equations for border irrigation was done by using a Eulerian–Lagrangian method based on [29], where the water depth over the soil was obtained by the numerical solution of the full Saint-Venant equations for mass and momentum conservation, while the infiltrated water depth was calculated by the finite element solution of the two-dimensional Richards’ equation. For a given border length, the Christiansen Uniformity Coefficient (CUC) was evaluated for a varied water intake flow. It was found that the CUC was very sensitive to the water intake flow; the water intake flow that corresponded to a maximum CUC was denominated optimal intake flow. The process was repeated for several lengths of a border and was found that a basically lineal relationship exist between the optimal intake flow and the border length. This last process was done for ten soil types of the texture triangle of USDA. Table 3 presents these values relative to three soil types, including optimal intake flow and the corresponding irrigation time. If use is made of a hydraulic resistance flow law of potential form [30], then the normal water depth for a border is obtained, which is transformed to an equivalent normal water depth for a furrow. Figure 3a–c exposes the evolution of equivalent furrow water depth at the beginning of the furrow, scaled from the evolution of border water depth at the beginning of the border obtained by the application of the described numerical model.
Through the internal numerical coupling of the Saint-Venant and Richards’ equations in border irrigation, Saucedo et al. [16] report values of optimal flow and corresponding irrigation time, for a slope of 0.002 m/m, indicated in Table 3. Considering normal water depth (perpendicular to the soil surface), a per unitary border width is obtained for a furrow corresponding to water depth.
Stability criteria was assumed per Banti et al. [8], and therefore a discretization was made in a space smaller than that used in the referred validation process, and different time steps mostly less than or equal to 30.0 s, so that the relative error maximum level was similar to that estimated in the referenced work.

3. Results and Discussion

The two-dimensional Richards’ equation solution was obtained applying finite element method combined with a temporal approximation based on implicit finite differences. The above was done by computational programming in C++ language. The model was applied in the solution domain to calculate the infiltrated depth, considering three contrasting soil textures. Results obtained considering different spatial-temporal discretizations were compared. Five meshes of finite element and five time steps were used. Results were evaluated using the root mean squared error (RMSE) [9]:
RMSE = 1 n i = 1 n ( M i m i ) 2
where n is the number of infiltrated depth data pairs, M i is the infiltrated depth obtained using the denser finite element mesh for a given time and m i is the infiltrated depth obtained with a different mesh from M i for the same time step.
Figure 4, Figure 5 and Figure 6 show the results of infiltrated depth for three different soil textures, with two denser finite element meshes and a time step of 1.0 s.
Figure 7, Figure 8 and Figure 9 present the volumetric soil water content (%). Typical irrigation times were used, and results presented consider a time step of 1.0 s, with spatial discretization according to a denser finite element mesh.
Figure 10a–c shows the values of the RMSE. These correspond to the results of infiltrated depth obtained with mesh 5 combined with meshes 1, 2, 3 and 4. With that combination the maximum error results; any other combination of meshes to determine the RMSE had a lower value than those presented below.
As presented above, results were compared considering five time steps as reference. Figure 11 shows the values of the RMSE indicator that compare the infiltrated depth evolution results, considering as standard the densest mesh combined with the shortest time step.
Instability of infiltrated depth evolution results for sandy loam soil were evident when time steps were greater than 10.0 s. This result is similar to studies reported in the literature [31,32].
Another aim was to compare results shown in Figure 3, Figure 4 and Figure 5 that were obtained when considering the infiltration model relative to the normal water depth indicated in Table 3 (constant water depth in time). The result is shown in Figure 12, Figure 13 and Figure 14. Mesh 5 and a time step of 1.0 s were used.
Relevant results are shown in Figure 12, Figure 13 and Figure 14. They make evident the convenience of considering time–water depth evolution in infiltration simulations. In fact, values of RMSE indicator shown in Figure 12, Figure 13 and Figure 14 are greater than values showed in Figure 10a–c.
To expand the scope of this study and for a fuller exploration of the subject, the results were compared with results obtained by applying the HYDRUS-2D program [33] (Version 2.0, PC-Progress, Prague, Czech Republic). The HYDRUS program numerically solves the two-dimensional Richards’ equation for saturated–unsaturated water flow. HYDRUS implements the soil-hydraulic functions of [21], using the statistical pore-size distribution model of [34] to obtain a predictive equation for the unsaturated hydraulic conductivity function in terms of soil water retention parameters. The expressions of [21] are given as follows:
θ ( ψ ) = { θ r + θ s θ r [ 1 + | α ψ | n ] m ψ < 0 θ s ψ 0
K ( ψ ) = K s S e l [ 1 ( 1 S e 1 / m ) m ] 2
where
m = 1 1 / n   , n > 1
The above equations contain six independent parameters: θ r , θ s , α , n , K s and l . The pore-connectivity parameter l in the hydraulic conductivity function was estimated [34] to be about 0.5 as an average for many soils. S e is the degree of saturation.
Results shown in Figure 12, Figure 13 and Figure 14 were compared with those obtained applying the normal water depth indicated in Table 3 (constant water depth in time) in the HYDRUS-2D program. Parameters used to implement the soil-hydraulic functions of van Genuchten–Mualem to obtain a predictive equation for unsaturated hydraulic conductivity function are shown in Table 4.
Results are shown in Figure 15, Figure 16 and Figure 17. A similar mesh to mesh 5 is used in the HYDRUS-2D program, as well as a time step of 1.0 s.
From the results shown in Figure 15, Figure 16 and Figure 17 it can be observed that the difference obtained with the HYDRUS-2D program and the model presented in this study were due to the different models used to obtain the humidity retention curve. It is possible to match the humidity retention curves of both models by adjusting the “m” parameter, but this is not enough to match the hydraulic conductivity functions. In addition, the difference between both models was due to spatial-temporal discretization, and the enforcement of a variable and constant water depth in each case.

4. Conclusions

By obtaining the time infiltrated depth evolution and humidity profiles, with a numerical solution of the two-dimensional Richards’ equation, the water infiltration during two-dimensional infiltration was simulated. The contact time hypothesis has been used to apply a unique form on time evolution of the water depth in the furrow throughout its length, as boundary condition. The specific form of such evolution in time was obtained from results reported in the literature based on the internal numerical coupling of the Saint-Venant and Richards’ equations in border irrigation. Moreover, the equivalent hydraulic area between the border and the furrow was achieved by scaling the values of water depth. The model was evaluated using three contrasting soil types and the comparison was done by computing the RMSE indicator (root mean square error). The comparison was performed from the selection of five finite element meshes with different densities, to discretize the solution domain of the two-dimensional Richards’ equation, combined with several time steps. Finally, a comparison was made between infiltrated depth evolution calculated with a constant water depth in the furrow to the one proposed in this work, finding important differences between both approaches. To expand the scope of this study and for a fuller exploration of the subject, the results are compared with results obtained by applying the HYDRUS-2D program. The results confirm that it is important to consider an internal full coupling of the Saint-Venant and Richards’ equations to improve furrow irrigation simulations. In other words, future efforts should be directed to the internal numerical full coupling of the Saint-Venant and two-dimensional Richards’ equations to adequately describe furrow irrigation.

Author Contributions

V.C. performed most of the analysis and numerical simulations showed in the study, he also contributed to the manuscript preparation; H.S. contributed to the design of the study and discussions of the results; C.F. revised the methodology, participated in the discussion of results and contributed to the final manuscript.

Funding

This research received no external funding.

Acknowledgments

The authors express their gratitude to the National Autonomous University of Mexico (UNAM), to the Irrigation and Drainage Coordination of the Mexican Institute of Water Technology (IMTA), and to the National Water Commission (CONAGUA) for their support in this investigation.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Richards, L.A. Capillary conduction of liquids through porous mediums. Physics 1931, 1, 318–333. [Google Scholar] [CrossRef]
  2. Liu, K.; Huang, G.; Xu, X.; Xiong, Y.; Huang, Q.; Šimůnek, J. A Coupled Model for Simulating Water Flow and Solute Transport in Furrow Irrigation. Agric. Water Manag. 2019, 213, 792–802. [Google Scholar] [CrossRef]
  3. Dong, Q.G.; Zhang, S.; Bai, M.; Xu, D.; Feng, H. Modeling the Effects of Spatial Variability of Irrigation Parameters on Border Irrigation Performance at a Field Scale. Water 2018, 10, 1770. [Google Scholar] [CrossRef]
  4. Duncan, S.J.; Daly, K.R.; Sweeney, P.; Roose, T. Mathematical modelling of water and solute movement in ridged versus flat planting systems. Eur. J. Soil Sci. 2018, 69, 967–979. [Google Scholar] [CrossRef]
  5. Fan, Y.; Huang, N.; Gong, J.; Shao, X.; Zhang, J.; Zhao, T. A Simplified Infiltration Model Predicting Cumulative Infiltration during Vertical Line Source Irrigation. Water 2018, 10, 89. [Google Scholar] [CrossRef]
  6. Brunetti, G.; Šimůnek, J.; Bautista, E. A Hybrid Finite Volume-Finite Element Model for the Numerical Analysis of Furrow Irrigation and Fertigation. Comput. Electron. Agric. 2018, 150, 312–327. [Google Scholar] [CrossRef]
  7. Furman, A. Modeling Coupled Surface–Subsurface Flow Processes: A Review. Vadose Zone J. 2008, 7, 741. [Google Scholar] [CrossRef]
  8. Banti, M.; Zissis, T.; Anastasiadou-Partheniou, E. Furrow Irrigation Advance Simulation using a Surface-Subsurface Interaction Model. J. Irrig. Drain. Eng. 2011, 137, 304–314. [Google Scholar] [CrossRef]
  9. Bautista, E.; Warrick, A.W.; Schlegel, J.L.; Thorp, K.R.; Hunsaker, D.J. Approximate furrow infiltration model for time-variable ponding depth. J. Irrig. Drain. Eng. 2016, 142, 04016045. [Google Scholar] [CrossRef]
  10. Bautista, E.; Warrick, A.W.; Strelkoff, T.S. New Results for an Approximate Method for Calculating Two-Dimensional Furrow Infiltration. J. Irrig. Drain. Eng. 2014, 140, 10. [Google Scholar] [CrossRef]
  11. Dong, Q.; Xu, D.; Zhang, S.; Meijian, B.; Yinong, L. A Hybrid Coupled Model of Surface and Subsurface Flow for Surface Irrigation. J. Hydrol. 2013, 500, 62–74. [Google Scholar] [CrossRef]
  12. Soroush, F.; Fenton, J.D.; Mostafazadeh-Fard, B.; Mousavi, S.F.; Abbasi, F. Simulation of furrow irrigation using the Slow-change/slow-flow equation. Agric. Water Manag. 2013, 116, 160–174. [Google Scholar] [CrossRef]
  13. Wöhling, T.; GSchmitz, H. Physically based coupled model simulating 1D surface-2D subsurface flow and plant water uptake in irrigation furrows. I. Model development. J. Irrig. Drain. Eng. 2007, 133, 6. [Google Scholar] [CrossRef]
  14. Wöhling, T.; Fröhner, A.; Schmitz, G.H.; Liedl, R. Efficient Solution of the Coupled One-Dimensional surface-Two Dimensional Subsurface flow during Furrow Irrigation Advance. J. Irrig. Drain. Eng. 2006, 132, 4. [Google Scholar] [CrossRef]
  15. Wöhling, T.; Singh, R.; Schmitz, G.H. Physically Based Modeling of Interacting Surface-Subsurface Flow during Furrow Irrigation Advance. J. Irrig. Drain. Eng. 2004, 130, 349–356. [Google Scholar]
  16. Saucedo, H.; Zavala, M.; Fuentes, C.; Castanedo, V. Optimal flow model for plot irrigation. Tecnología y Ciencias del Agua 2013, 4, 135–148. [Google Scholar]
  17. López Avendaño, J.E.; Palacios, O.L.; Fuentes, C.; Rendón, L.; García, N. Two-dimensional analysis infiltration in furrow irrigation. Agrociencias 1997, 31, 259–269. [Google Scholar]
  18. Fuentes, C.; Haverkamp, R.; Parlange, J.Y. Parameter constraints on closed–form soil-water relationships. J. Hydrol. 1992, 134, 117–142. [Google Scholar] [CrossRef]
  19. Fujita, H. The exact pattern of a concentration-dependent diffusion in a semi-infinite medium, part I. Text. Res. J. 1952, 22, 757–761. [Google Scholar] [CrossRef]
  20. Parlange, J.Y.; Lisle, I.; Braddock, R.D.; Smith, R.E. The three parameter infiltration equation. Soil Sci. 1982, 133, 337–341. [Google Scholar] [CrossRef]
  21. Van Genuchten, M.T. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci. Soc. Am. J. 1980, 44, 892–898. [Google Scholar] [CrossRef]
  22. Burdine, N.T. Relative permeability calculation from size distributions data. Trans. AIME 1953, 198, 171–199. [Google Scholar] [CrossRef]
  23. Brooks, R.H.; Corey, A.T. Hydraulic Properties of Porous Media; Colorado State University: Fort Collins, CO, USA, 1994. [Google Scholar]
  24. Noor, K.A.; Peters, J.M. Preconditioned Conjugate Gradient Technique for the analysis of symmetric structures. Int. J. Numer. Methods Eng. 1987, 24, 2057–2070. [Google Scholar] [CrossRef]
  25. Neumann, S.P. Saturated-unsaturated seepage by finite elements. J. Hydraul. Div. 1973, 12, 2233–2250. [Google Scholar]
  26. Mori, M. The Finite Element Method and Its Applications; Macmillan: New York, NY, USA, 1983. [Google Scholar]
  27. Rawls, W.J.; Brakensiek, D.L. Estimating Soil Water Retention from Soil Properties. Journal of Irrigation and Drainage. Div. Am. Soc. Civil Eng. 1982, 108, 166–171. [Google Scholar]
  28. Saucedo, H.; Fuentes, C.; Zavala, M. El sistema de ecuaciones de Saint-Venant y Richards del riego por gravedad: 3. Verificación numérica de la hipótesis del tiempo de contacto en el riego por melgas. Ingeniería Hidráulica en México 2006, 4, 135–143. [Google Scholar]
  29. Strelkoff, T.; Katopodes, N. Border irrigation hydraulics with zero inertia. J. Irrig. Drain. Div. 1977, 3, 325–342. [Google Scholar]
  30. Fuentes, C.; de León, B.; Saucedo, H.; Parlange, J.Y. Saint-Venant and Richards’ equations system of surface irrigation: The potential law of hydraulic resistance. Ingeniería Hidráulica en México 2004, 19, 65–77. [Google Scholar]
  31. Selker, J.S.; Steenhuis, T.S.; Parlange, J.Y. Wetting front instability in homogeneous Sandy soils under continuous infiltration. Soil Sci. Soc. Am. J. 1991, 56, 1346–1350. [Google Scholar] [CrossRef]
  32. Hendrickx, J.M.H.; Tzung-Mow, Y. Prediction of wetting front stability in dry field soils using soil and precipitation data. Geoderma 1996, 70, 265–280. [Google Scholar] [CrossRef]
  33. Šimůnek, J.; van Genuchten, M.T.; Šejna, M. HYDRUS Technical Manual 2.0. Available online: http://www.pc-progress.com/downloads/Pgm_Hydrus3D2/HYDRUS3D%20Technical%20Manual.pdf (accessed on 20 February 2019).
  34. Mualem, Y. A new model for predicting the hydraulic conductivity of unsaturated porous media. Water Resour. Res. 1976, 12, 513–522. [Google Scholar] [CrossRef]
Figure 1. Solution domain for two-dimensional Richards’ equation. Z is the vertical spatial coordinate, X is the horizontal spatial coordinate, P is the soil depth profile, PS is the furrow depth and h is the water depth in furrow center. A, B, C, D, E and F are frontier points of solution domain used to define border conditions.
Figure 1. Solution domain for two-dimensional Richards’ equation. Z is the vertical spatial coordinate, X is the horizontal spatial coordinate, P is the soil depth profile, PS is the furrow depth and h is the water depth in furrow center. A, B, C, D, E and F are frontier points of solution domain used to define border conditions.
Water 11 00371 g001
Figure 2. Finite element mesh for two-dimensional Richards’ equation numerical solution (mesh No 1).
Figure 2. Finite element mesh for two-dimensional Richards’ equation numerical solution (mesh No 1).
Water 11 00371 g002
Figure 3. (a) Approximate water depth evolution form for sandy loam soil. (b) Approximate water depth evolution form for silt loam soils. (c) Approximate water depth evolution form for clay loam soil.
Figure 3. (a) Approximate water depth evolution form for sandy loam soil. (b) Approximate water depth evolution form for silt loam soils. (c) Approximate water depth evolution form for clay loam soil.
Water 11 00371 g003aWater 11 00371 g003b
Figure 4. Infiltrated depth evolution for sandy loam soil. Typical irrigation time of 1.2 h. RMSE = 0.160 cm.
Figure 4. Infiltrated depth evolution for sandy loam soil. Typical irrigation time of 1.2 h. RMSE = 0.160 cm.
Water 11 00371 g004
Figure 5. Infiltrated depth evolution for silt loam soil. Typical irrigation time of 6.2 h. RMSE = 0.063 cm.
Figure 5. Infiltrated depth evolution for silt loam soil. Typical irrigation time of 6.2 h. RMSE = 0.063 cm.
Water 11 00371 g005
Figure 6. Infiltrated depth evolution for clay loam soil. Typical irrigation time of 31.4 h. RMSE = 0.019 cm.
Figure 6. Infiltrated depth evolution for clay loam soil. Typical irrigation time of 31.4 h. RMSE = 0.019 cm.
Water 11 00371 g006
Figure 7. Humidity distribution in sandy loam soil, using a typical irrigation time of 1.2 h. θ s = 0.450 cm3/cm3 (100% saturation).
Figure 7. Humidity distribution in sandy loam soil, using a typical irrigation time of 1.2 h. θ s = 0.450 cm3/cm3 (100% saturation).
Water 11 00371 g007
Figure 8. Humidity distribution in silt loam soil, using a typical irrigation time of 6.2 h. θ s = 0.525 cm3/cm3 (100% saturation).
Figure 8. Humidity distribution in silt loam soil, using a typical irrigation time of 6.2 h. θ s = 0.525 cm3/cm3 (100% saturation).
Water 11 00371 g008
Figure 9. Humidity distribution in clay loam soil, using a typical irrigation time of 31.4 h θ s = 0.475 cm3/cm3 (100% saturation).
Figure 9. Humidity distribution in clay loam soil, using a typical irrigation time of 31.4 h θ s = 0.475 cm3/cm3 (100% saturation).
Water 11 00371 g009
Figure 10. RMSE values comparing infiltrated depth for the three soil textures and five finite element meshes, (a) time steps of 1.0 and 5.0 s; (b) time steps of 10.0 and 30.0 s; (c) time steps of 60.0 s.
Figure 10. RMSE values comparing infiltrated depth for the three soil textures and five finite element meshes, (a) time steps of 1.0 and 5.0 s; (b) time steps of 10.0 and 30.0 s; (c) time steps of 60.0 s.
Water 11 00371 g010aWater 11 00371 g010b
Figure 11. Values of RMSE for different time steps.
Figure 11. Values of RMSE for different time steps.
Water 11 00371 g011
Figure 12. Infiltrated depth evolution for sandy loam soil. The dotted line corresponds to results considering a constant water depth in time (normal water depth). The continuous line corresponds to results considering a variable water depth (proposed in this study). Typical irrigation time of 1.2 h. RMSE = 0.460 cm.
Figure 12. Infiltrated depth evolution for sandy loam soil. The dotted line corresponds to results considering a constant water depth in time (normal water depth). The continuous line corresponds to results considering a variable water depth (proposed in this study). Typical irrigation time of 1.2 h. RMSE = 0.460 cm.
Water 11 00371 g012
Figure 13. Infiltrated depth evolution for silt loam soil. The dotted line corresponds to results considering a constant water depth in time (normal water depth). The continuous line corresponds to results considering a variable water depth (proposed in this study). Typical irrigation time of 6.2 h. RMSE = 0.292 cm.
Figure 13. Infiltrated depth evolution for silt loam soil. The dotted line corresponds to results considering a constant water depth in time (normal water depth). The continuous line corresponds to results considering a variable water depth (proposed in this study). Typical irrigation time of 6.2 h. RMSE = 0.292 cm.
Water 11 00371 g013
Figure 14. Infiltrated depth evolution for clay loam soil. The dotted line corresponds to results considering a constant water depth in time (normal water depth). The continuous line corresponds to results considering a variable water depth (proposed in this study). Typical irrigation time of 31.4 h. RMSE = 0.267 cm.
Figure 14. Infiltrated depth evolution for clay loam soil. The dotted line corresponds to results considering a constant water depth in time (normal water depth). The continuous line corresponds to results considering a variable water depth (proposed in this study). Typical irrigation time of 31.4 h. RMSE = 0.267 cm.
Water 11 00371 g014
Figure 15. Infiltrated depth evolution for sandy loam soil. Typical irrigation time of 1.2 h. RMSE = 0.721 cm was obtained during comparison of results obtained with a time–water depth evolution and results obtained with the HYDRUS-2D program.
Figure 15. Infiltrated depth evolution for sandy loam soil. Typical irrigation time of 1.2 h. RMSE = 0.721 cm was obtained during comparison of results obtained with a time–water depth evolution and results obtained with the HYDRUS-2D program.
Water 11 00371 g015
Figure 16. Infiltrated depth evolution for silt loam soil. Typical irrigation time of 6.2 h. RMSE = 0.516 cm was obtained during comparison of results obtained with a time–water depth evolution and results obtained with the HYDRUS-2D program.
Figure 16. Infiltrated depth evolution for silt loam soil. Typical irrigation time of 6.2 h. RMSE = 0.516 cm was obtained during comparison of results obtained with a time–water depth evolution and results obtained with the HYDRUS-2D program.
Water 11 00371 g016
Figure 17. Infiltrated depth evolution for clay loam soil. Typical irrigation time of 31.4 h. RMSE = 0.365 cm was obtained during comparison of results obtained with a time–water depth evolution and results obtained with the HYDRUS-2D program.
Figure 17. Infiltrated depth evolution for clay loam soil. Typical irrigation time of 31.4 h. RMSE = 0.365 cm was obtained during comparison of results obtained with a time–water depth evolution and results obtained with the HYDRUS-2D program.
Water 11 00371 g017
Table 1. Characteristics of finite element meshes.
Table 1. Characteristics of finite element meshes.
MeshElementsNodes Maximum   Increase   Δ x   ( cm ) Maximum   Increase   Δ z   ( cm )
157533410.010.0
214017645.05.0
3812241972.02.0
432,01416,2721.01.0
5127,73264,3870.50.5
Table 2. Hydrodynamic characteristics for three soil textures.
Table 2. Hydrodynamic characteristics for three soil textures.
Soil θ r   ( cm 3 / cm 3 ) θ s   ( cm 3 / cm 3 ) ψ d   ( cm ) K s   ( m / s ) m = 1 2 / n η
Sandy loam0.00.4509.521.390 × 10−40.100413.62
Silt loam0.00.52529.350.0167 × 10−40.116512.01
Clay loam0.00.47534.150.0042 × 10−40.071419.30
Table 3. Indicators to obtain time–water depth evolution form, reported by Saucedo et al. [16].
Table 3. Indicators to obtain time–water depth evolution form, reported by Saucedo et al. [16].
Soil
Sandy LoamSilt LoamClay Loam
Qopt (lps/m)0.024760.004460.00088
Irrigation time (h)1.26.231.4
Normal water depth (border irrigation) (cm)2.420.980.50
Normal water depth (furrow irrigation) (cm)8.924.993.24
Table 4. Hydrodynamic characteristics for three soil textures used in HYDRUS-2D.
Table 4. Hydrodynamic characteristics for three soil textures used in HYDRUS-2D.
Soil θ s   ( cm 3 / cm 3 ) K s   ( m / s ) m = 1 1 / n n (n > 1)l
Sandy loam0.4501.390 × 10−40.10041.11160.5
Silt loam0.5250.0167 × 10−40.11651.13180.5
Clay loam0.4750.0042 × 10−40.07141.07690.5

Share and Cite

MDPI and ACS Style

Castanedo, V.; Saucedo, H.; Fuentes, C. Modeling Two-Dimensional Infiltration with Constant and Time-Variable Water Depth. Water 2019, 11, 371. https://doi.org/10.3390/w11020371

AMA Style

Castanedo V, Saucedo H, Fuentes C. Modeling Two-Dimensional Infiltration with Constant and Time-Variable Water Depth. Water. 2019; 11(2):371. https://doi.org/10.3390/w11020371

Chicago/Turabian Style

Castanedo, Vladimir, Heber Saucedo, and Carlos Fuentes. 2019. "Modeling Two-Dimensional Infiltration with Constant and Time-Variable Water Depth" Water 11, no. 2: 371. https://doi.org/10.3390/w11020371

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop