# Effects of the Grain Size on Dynamic Capillary Pressure and the Modified Green-Ampt Model for Infiltration.

1. Introduction

Infiltration involves gas and liquid flows in porous media and occurs during precipitation or when liquid contaminants leak underground or onto the soil surface. This contributes to runoff generation, crop irrigation, and transport of nutrients and contaminants. Infiltration is a main process in subsurface hydrology and plays a crucial role in geohazards, such as landslides, flooding, and groundwater contamination. Richards' equation is the most general model to address such flows with spatially and temporally variable saturation [1, 2]. Nevertheless, computational time is an issue of numerically solving Richards' equation in a large-scale simulation. Therefore, the simple one-dimensional Green-Ampt model (GAM)  is an attractive alternative [4, 5]. Actually, the GAM has been widely incorporated in large-scale hydrological processes and erosion models  such as HEC-HMS , WEPP , SWAT , and ANSWERS  models. However, describing the transient behavior during early infiltration based on the GAM remains a challenging task.

The classical GAM for downward infiltration can be expressed as an ordinary differential equation (ODE) [3, 11]:

[[theta].sub.s] - [[theta].sub.i]/[K.sub.s] [dl.sub.f]/dt [l.sub.f] = [h.sub.w] + [s.sub.f] + [l.sub.f], (1)

where [[theta].sub.s] is the saturated water content, [[theta].sub.i] is the initial water content, [K.sub.s] is the saturated hydraulic conductivity, [l.sub.f] is the length of the porous medium with infiltrated liquid (the wetting front depth), [v.sub.f] = [dl.sub.f]/dt is the wetting front velocity, [s.sub.f] = [P.sub.c]/[rho]g is the equilibrium suction head, [P.sub.c] is capillary pressure, [rho] is the water density, g is the gravitational acceleration, and [h.sub.w] is the ponding height. In the model, a wetting front moves downward as a sharp moving boundary during downward infiltration. The porous media above and below the wetting front are assumed to have water saturated with [[theta].sub.s] and unsaturated with [[theta].sub.i], respectively. In addition, the air pressure in the porous medium is assumed to be constant and, therefore, the viscous pressure drop due to the air movement is neglected. The capillary pressure [P.sub.c] across the moving wetting front is commonly assumed to be constant and determined from the soil water retention curve [12, 13].

Nevertheless, during the past few decades, studies conducted under nonequilibrium conditions have indicated that the soil water retention curve may depend on the dynamics of water flow [14, 15]. Water content measured under transient flow conditions has been shown to be significantly different from that measured under static and steady-state conditions [16-21]. Considering the soil water retention curve based on thermodynamics, Hassanizadeh and Gray , postulated the existence of a dynamic component in the unsaturated water flow. This dynamic component depends on both flow dynamics and the process of either drainage or imbibition. In addition, studies conducting column experiments have shown that during downward infiltration, the water pressure head at the gas-water or liquid-liquid interfaces changes as the flow velocity changes [11, 23-28]. The aforementioned studies suggest that the capillary pressure under dynamic conditions can be less than that under static conditions during infiltration. Therefore, infiltration can be better described by a GAM with a nonequilibrium suction head [s.sub.f].

The modified Green-Ampt model (MGAM) was first developed in Hsu and Hilpert . They verified the model by the experimental data from upward infiltration , and downward infiltration . Hsu and Hilpert  showed that the MGAM is better than the classical GAM at describing both transient capillary rise and downward infiltration in dry sands. For downward infiltration, the MGAM can also be expressed as an ordinary differential equation :

[[theta].sub.s] - [[theta].sub.i]/K [dl.sub.f]/dt [l.sub.f] = [h.sub.w] + [s.sub.f] - [gamma]/d[rho]g [??] [([eta]/[gamma] [dl.sub.f]/dt).sup.[beta]] + [l.sub.f], (2)

where [gamma] is the interfacial tension, d is the grain size, [??] = [alpha][epsilon] ([epsilon], [[theta].sub.i]) is the lumped parameter related to a nondimensional function ([epsilon]) of porosity [epsilon] and [[theta].sub.i], [alpha] and [beta] are model parameters, and [eta] is the dynamic viscosity of water. The additional term ([gamma]/d[rho]g)[??][[([eta]/[gamma])([dl.sub.f]/dt)].sup.[beta]] is added into GAM to consider the dynamic effect of capillary pressure.

The physical concept of the MGAM and the additional parameters [??] and [beta] are based on the theory of the pore-scale dynamic contact angle. The upscaling of dynamic contact angle to the dynamic effect of capillary pressure has been derived and discussed in Hsu and Hilpert . Based on the studies of dynamic contact angle, [beta] is in the range of 1/2 to 1/3, but the value of a highly relies on the solid surface structure [31-34].

Hsu et al. and Pellichero et al. [11, 26] showed that the MGAM approach can avoid the initial unphysical infinite rate of infiltration and better describe the experimental results than can the classical GAM for downward infiltration at an early stage in both dry and prewetted sand columns with varied initial water contents. Indications are that the grain size and pore size distributions affect the dynamic effect on the soil water retention curve [20, 35-37]. DiCarlo  showed that when a saturation overshoot occurs during the downward infiltration, the length of the gravitational fingering tips varies with the grain size and distribution. Hilpert  proposed a theory based on velocity-dependent capillary pressure that correctly predicts the degree to which the saturation overshoot depends on both the liquid contents and grain size. Hsu and Hilpert  pointed out that the model parameters [??] and [beta] in MGAM might be related to the porosity and grain size. Based on (2), the dynamic effect on the suction head should be inversely proportional to the grain size. However, the effect of the grain size and pore size distribution of the sand column on the nonequilibrium suction head as well as the model parameters [??] and [beta] in (2) of the MGAM have not yet been systematically investigated. Therefore, we performed a series of infiltration experiment to systematically investigate the effect of the grain size of the sand column on the nonequilibrium suction head as well as the model parameters [??] and [beta].

In this study, the cumulative infiltration depths were measured during a series of downward infiltration experiments in dry sand columns with different grain sizes, subject to different constant ponding heights. The remainder of this paper is organized as follows. In Section 2, we describe our experimental method and preparation of the porous materials. In Section 3, we describe the results of the experiments, and in Section 4, we compare the predictions from classical GAM and MGAM as well as the values of fitting parameters. We also discuss the effects of the grain size on the dynamic capillary pressure and the fitting of the equilibrium suction head.

2. Experimental Setup

2.1. Experimental Material Properties and Grain Size Measurement. In this study, sands (glass beads) with four grain sizes (labeled B35, B60, B150, and B320) were used for our infiltration experiments. The sands had the same shape, with a particle density of 2.5 g/[cm.sup.3]. Because glass beads after conditioning hardly aggregate, we could exclude the influence of the aggregates on the observed dynamic capillary pressure. Moreover, the grain size distribution (GSD) of the sands was determined through sieve analysis (ASTM C136) of eight different mesh screens.

2.2. Downward Infiltration Experiment. The infiltration experiments were conducted using glass filter columns (inner diameter = 2.6 cm, depth = 60 and 30 cm, cross-section area = 5.3 [cm.sup.2]). Based on Wang et al. , the sizes of the fingering flow for downward infiltration in sands are mostly ranged from 3 to 15 cm. The fingering flow with a width of less than 3 cm was only observed by Glass et al. . In this study, our sand column with the inner diameter of 2.6 cm is small enough to avoid 2D fingering flows. All the sands were conditioned by being rinsed with deionized water to remove impurities and dried overnight in an oven at 100[degrees]C. We packed sands B35, B60, and B150 into the 60 cm column and sand B320 into the 30 cm column. To achieve uniform packing, the sands poured into the column were maintained at a constant 3 cm distance between the supply funnel and the top of the sand packing. In addition, a rubber hammer was used to tap the top of the sand in the column to obtain an even more homogeneous packing.

Figure 1 is a schematic representation of the experiment setup. A Mariotte's bottle was connected through Tygon tubing and a valve to the sand column to maintain the hydraulic head. The bottle was placed on the top of an analytical balance (Sartorius TE1502S) to record the weight at 0.2 s intervals and automatically send the data to a computer. When the valve was opened, the water in the Mariotte's bottle flowed into the column packed with dry sand, infiltrated the sand and reached the bottom of the column. The weight change measured by the analytical balance was used to calculate the cumulative infiltration and infiltration rate. All infiltration experiments in this study were conducted in triplicate.

2.3. Saturated Hydraulic Conductivity ([K.sub.s]) Measurement. Saturated hydraulic conductivity ([K.sub.s]) is an important parameter of the GAM, as it determines the flow rate of water in the saturated soil. [K.sub.s] in this study was determined by the constant-head method. The measurement procedure was the same as that for the infiltration experiment previously described. Following the infiltration experiments, the sand in the column was at saturation, and the flow rate gradually reached equilibrium. We measured the equilibrium infiltration rate, and [K.sub.s] was estimated based on Darcy's law:

[K.sub.s] = QL/Ah, (3)

where Q is the flow rate, L is the length of the sand column, A is the total cross-section area, and h is the constant hydraulic head. The flow rate value was obtained by calculating the weight change of water per unit time.

3. Results

3.1. Grain Size Distribution of Experimental Materials. Figure 2 shows the grain size distributions of the sands. Sand B35 was the coarsest in grain size, ranging from 2 to 7 x [10.sup.-2] cm and was followed in order by B60 (distributed between 1 and 6 x [10.sup.-2] cm), B150 (ranging from 5 x [10.sup.-3] to 2.3 x [10.sup.-2] cm), and B320 (which had the finest grain size and was distributed between 1 x 10-4 and 2 x [10.sup.-2] cm).

The curves in Figure 2 are plotted based on a grain size distribution equation provided by Fredlund et al. . We used the same equation to calculate the geometric mean grain size ([d.sub.50]). Table 1 shows the geometric mean grain sizes of sands B35 (4.51 x [10.sup.-2] cm), B60 (3.48 x [10.sup.-2] cm), B150 (1.38 x [10.sup.-2] cm), and B320 (1.3 x [10.sup.-3] cm).

The sands were uniform. Specifically, the uniformity coefficients of sands B35, B60, B150, and B320 were 1.33, 1.46, 1.55, and 3.85, respectively. In the discussion that follows, we used [d.sub.50] to quantify the effect of grain size on the infiltration model and dynamic capillary pressure.

3.2. Results of Infiltration Experiments. Table 1 shows the measured porosity, bulk density, and saturated hydraulic conductivity. In each of the 36 infiltration replicate experiments, we controlled the bulk density and porosity of the soil column. The average bulk density of sand B35 column was 1.36 g/[cm.sup.3] and the average porosity was 0.46[cm.sup.3]/[cm.sup.3]. The bulk density of sand B60 was 1.42-1.43 g/[cm.sup.3], B150 was 1.38-1.39 g/[cm.sup.3], and B320 was 1.33-1.36 g/[cm.sup.3]. The porosity of sands B60, B150, and B320 was 0.43, 0.45,

and 0.47 [cm.sup.3]/[cm.sup.3], respectively.

Under the assumption of GAM, the sand has a water content [theta] = [[theta].sub.s] behind the wetting front and an initial water content [[theta].sub.i] ahead of the front. Mobile water content ([DELTA][theta]) indicates the moisture change ([[theta].sub.s] - [[theta].sub.i]) at the wetting front and can be estimated from the total cumulative infiltration and total volume of the sand column, that is, [DELTA][theta] = [F.sub.total]/[V.sub.total] . In our study, the [DELTA][theta] of sand B35 column was obtained from 0.42 to 0.44 [cm.sup.3]/[cm.sup.3], B60 from 0.39 to 0.41 [cm.sup.3]/[cm.sup.3], B150 from 0.40 to 0.42[cm.sup.3]/[cm.sup.3], and B320 from 0.39 to 0.44 [cm.sup.3]/[cm.sup.3].

The results of the saturated hydraulic conductivity from additional experiments show that B35 had the largest value, followed by B60, B150, and B320, at 1.84-1.97 x [10.sup.-1], 6.42-6.53 x [10.sup.-2], 1.16-1.42 x [10.sup.-2], and 0.81.1 x [10.sup.-3] cm/s, respectively.

We used the analytical balance to measure the mass ([m.sub.w]) of water infiltrating the sand column. In addition, the cumulative infiltration (F) was calculated from the mass ([m.sub.w]), water density ([rho]), and soil column cross-sectional area (A), such that F = [m.sub.w]/[rho]A. We plotted the relationship between F and time (t) according to different water ponding heights and soil sizes. Figure 3 shows F versus t in our experiments. Each plotted point represents the average result of three replicate experiments. However, one of the infiltration experiments of soil B320 column at a water ponding height of 20 cm failed because of water leakage, a fact that was not mentioned in our subsequent discussion. Figures 3(a)-3(d) clearly show that F increased as the water ponding height increased.

3.3. Wetting Front Depth and Velocity versus Water Ponding Height and Grain Size. Wetting front depth ([l.sub.f]) can be calculated from F and [DELTA][theta], such that [l.sub.f] = F/[DELTA][theta]. Figure 4 shows [l.sub.f] versus t in our experiments. In general, the [l.sub.f] monotonously increased with t. In addition, water ponding height [h.sub.w] = 40 cm had the largest [l.sub.f] value at the same time, followed by [h.sub.w] = 20 cm and [h.sub.w] = 10 cm. For different grain sizes, the wetting front of B35 with the largest grain size and fastest infiltration rate spent 55-80.8 s to arrive [l.sub.f] = 60 cm, B60 spent 125-217.6 s, and B150 spent 597.2-980 s. Finally, the wetting front of B320 with the smallest grain size and slowest infiltration rate spent 2040-2430 s to arrive [l.sub.f] = 30 cm.

We divided [dl.sub.f] by dt to obtain the wetting front velocity ([v.sub.f]) and plotted the relationship of [v.sub.f] and t in Figure 5. To compare the wetting front velocities in different experiments with t, we normalized the X axis through t divided by the infiltration end time ([t.sub.ie]). Figure 5 shows [v.sub.f] versus standardized time (t/[t.sub.ie]). We can see that v f was the fastest at the beginning of infiltration, gradually decreased with t, and finally approached an equilibrium value. Regardless of the water ponding height, [v.sub.f] increased as the grain size increased.

The main purpose of this study is to examine the relationship between flow velocity and capillary pressure during transient infiltration. The combination of ponding depths and grain sizes can lead to a range of flow rates and front velocities for examining our hypothesis. In addition, our experimental results shown in Figures 4 and 5 confirm that our experimental design in accord with the basic physical assumptions of infiltration theory and the GAM. It implies that the fingering flow did not occur during our experiments.

3.4. Infiltration Simulation of GAM and MGAM. In this study, experimental data were modeled by GAM and MgAm. We used a MATLAB implementation of ODE solvers to provide a numerical solution of GAM (1) and MGAM (2) in downward infiltration conditions. GAM required only a parameter [s.sub.f] for simulation, whereas other parameters [K.sub.s] and [DELTA][theta] were inputted as fixed values (listed in Table 1). Parameters used by MGAM included [s.sub.f], [??], [beta], d, [K.sub.s], and [DELTA][theta], where [K.sub.s] and [DELTA][theta] were inputted as fixed values (listed in Table 1). Here, grain size (d) was represented by the geometric mean grain size ([d.sub.50]).

The fitting parameters were inversely determined from the measured data using a MATLAB nonlinear fitting function. The nonlinear fitting function uses a Levenberg-Marquardt nonlinear least squares algorithm to solve bound-constrained optimization problems . MATLAB nonlinear fitting function must specify the initial value of [s.sub.f], [??], [beta] and d. Based on the studies of Hsu et al.  and Pellichero et al. , we set [??] = 100 and [beta] = 0.3 as the initial values for model fitting. The initial value of d of sands B35, B60, B150, and B320 was [d.sub.50] =4.51 x [10.sup.-2], 3.48 x [10.sup.-2], 1.38 x [10.sup.-2], and 1.3 x [10.sup.-3] cm, respectively. The initial value of [s.sub.f] was obtained by

[s.sub.f] = 2[gamma] cos [phi]/r[rho]g, (4)

where [gamma] is interfacial tension (in which the surface tension of water and air is [gamma] = 7.2 x [10.sup.-2] N/m), g is gravity acceleration (g = 9.81 m/[s.sup.2]), p is water density ([rho] = 1000 [m.sup.3]), [phi] is contact angle, and r is effective pore radius. In this study, the contact angle was [phi] = 30[degrees]; this was the static contact angle of quartz and water suggested by Friedman . We used the empirical formula of Glover and Walker  to convert grain size (d) into effective pore size (r):

r = d/2 [square root of (8/3 x [1.5.sup.2]/8[[epsilon].sup.2x1.5])], (5)

where [epsilon] is porosity. The aforementioned initial values were substituted in (1) and (2) to simulate GAM and MGAM for downward infiltration, respectively.

The results of infiltration simulation of GAM and MGAM are shown in Figure 6, and the model fitting parameters are listed in Table 2. In Figure 6, the black circles are the wetting front depths, which were calculated from the measured values of [m.sub.w]. The remaining four curves show the simulation results of GAM and MGAM. The orange-dotted curve (GAM) was plotted by GAM with the calculated [s.sub.f] based on (4). The yellow-dashed curve (MGAM) was simulated by MGAM with [??] = 100, [beta] = 0.3, [s.sub.f] (which was calculated by (4)), and d = [d.sub.50]. The purple-dashed curve (optimized GAM) used the nonlinear fitting method to simulate downward infiltration by GAM with [s.sub.f] as a fitting parameter. Finally, the green-solid curve (optimized MGAM) was represented by MGAM with [??], [beta], [s.sub.f], and d as fitting parameters. Both the purple-dashed and green-solid curves with fitting parameters showed better simulation results; these curves could approach the black circles of the wetting front depth. The orange-dotted curve with a fixed [s.sub.f] was farthest from the black circles, followed by the yellow-dashed curve with fixed [??], [beta], sf, and d in the middle. This indicates that the performance of MGAM with multiparameters for simulating downward infiltration was superior to that of the classical GAM.

Regardless of the grain size of the sand, the optimized GAM and optimized MGAM with the fitting parameters had good simulation results. For GAM and MGAM with no fitting parameters, the smaller the grain size, the more the simulation results deviated from the black circles, particularly the sand B320.

Table 2 shows that two of the fitting suction heads [s.sub.f] of sand B35 based on GAM were negative. Because the grain size of the B35 was larger than that of other sands, the equilibrium [s.sub.f] was relatively small or even close to zero (Figure 7 and Table 2). Nevertheless, to have a negative suction head for a sand column packed with glass beads is still abnormal. It implies that the [s.sub.f] was underestimated by using classical GAM during a downward infiltration with high wetting front velocity. According to Chow et al. , the [s.sub.f] is distributed from 0.97 cm (0.0097 m) to 156.5 cm (1.565 m) from sand to clay. From the grain size distribution in Figure 2, the texture class of B320 can be classified as silt loam. Its [s.sub.f] value ranges from 2.92 cm (0.0292 m) to 95.39 cm (0.9539 m). However, theoretical values of sf based on (4) are over 2 m. The [s.sub.f] value for sand B320 is larger than that from Chowet al. , but it should still be in the range of silt loam or soils slightly smaller than silt loam for practical problems.

The performances of GAM and MGAM were evaluated using root mean square error (RMSE). For each simulation, the difference between the measured wetting front depth ([l.sub.f(m)]) and the calculated wetting front depth ([l.sub.f(p)]) was expressed in RMSE as:

RMSE = [square root of [([[summation].sup.N.sub.i=1] ([l.sub.f(m)] - [l.sub.f(p)]).sup.2]/N, (6)

where N is the number of observed samples. By statistically verifying RMSE, we could evaluate the reliability of the model simulation; when the RMSE was small, the model was more reliable. The results of RMSE for infiltration simulation of GAM and MGAM are listed in Table 2. In general, the RMSEs of most MGAMs (0.001-0.008) in this study were less than those of GAMs (0.002-0.011), indicating that the simulation results of MGAM were superior to those of GAM. Except for two sets of simulation results for B60 and B320, the RMSE of GAM was smaller than that of MGAM. Fortunately, the RMSE of these two sets of results were also very small, which meant that the performance of MGAM was also good.

4. Discussion

4.1. Relationship between Wetting Front Suction Head Sf) and Grain Size. The relationship between the suction head and grain size can be linked by using (4) and (5). Figure 7 shows the relationship between the grain sizes of the sand columns and the equilibrium suction heads estimated by GAM and MGAM. The solid line in Figure 7 shows the theoretical relationship based on (4) and (5) between the grain sizes of the sand columns and the equilibrium suction heads with the assumed effective contact angle between air, water, and the glass bead of 30[degrees]. Using the same equations, we determined the effective pore radius and contact angles by fitting the effective contact angle to the grain size and the determined equilibrium suction heads. The dashed lines in Figure 7 show the fitted relationship between grain sizes and the equilibrium suction heads estimated by GAM and MGAM with the fitted effective contact angle of 84[degrees] and 63[degrees], respectively.

The equilibrium suction heads estimated by the GAM were systematically lower than those estimated by the MGAM. Consequently, the effective contact angle estimated by the GAM was also systematically greater than that by the MGAM. The overestimated effective contact angle from GAM was because we neglected the dynamic effect when estimating the suction head. The effective contact angle estimated by GAM was the average dynamic contact angle during the infiltration process. Nevertheless, the contact angle (63[degrees]) estimated by MGAM was still greater than that (30[degrees]) measured in previous studies, even when the dynamic effect was considered. One explanation for this is related to the contact angle hysteresis: the contact angle no longer corresponds to the static equilibrium angle but is larger when the liquid is advancing and smaller when the liquid is receding . In general, during the infiltration, the water advances the air in the sand. Therefore, the contact angle during the infiltration should be in the form of an advancing contact angle, which is commonly greater than the equilibrium contact angle.

4.2. MGAM Fitting Parameters for Dynamic Capillary Pressure. Hsu and Hilpert  pointed out that the parameters [??] and [beta] in the MGAM are related to the dynamic effect on the capillary pressure, and their values might not be the same for different initial water contents, porosities, and grain sizes. Figure 8 shows two trends in which the fitted [??] is proportional to the porosity and inversely proportional to the grain size. From the empirical relationship proposed by Stauffer  relating the dynamic coefficient to measurable porous media and fluid properties, the equation suggests that the dynamic effects would be greater for fine-grained soil with high air entry pressure and high porosity. Camps-Roach eat al.  observed that the dynamic coefficient statistically increases as the grain size decreases. Therefore, our fitted values of the parameters were consistent with those of the aforementioned studies.

Limited numbers of Darcy-scale studies focus on determining the [beta] values in sand. Weitz et al.  yielded [beta] = 0.5 [+ or -] 0.1 for water displacing decane in a glass bead column. Hsu and Hilpert  tested a range of [beta] values from 1/5 to 1 based on the studies of pore-scale dynamic contact angles on flat surfaces and in circular capillary tubes, which can be described by cos [[phi].sub.eq] - cos [[phi].sub.dyn] = [??][Ca.sup.[beta]] [33, 34, 46]. The best [beta] value was approximately 0.31 for the experiments performed by Geiger and Durnford . Pellichero et al.  showed that [beta] = 0.3 resulted in the best fit in downward infiltration experiments with the same sand but varying ponding heights. The same [beta] value also worked well in similar experiments but with various initial water contents . In this study, [beta] was set as one of the fitting parameters. Figure 9 shows that the fitted [beta] was approximately 0.3, and the error bars overlapped for most grain sizes except for those of B320, which had the highest [beta] value (up to 0.34) with the smallest grain size. Nevertheless, the value of the fitted [beta] in this study agreed with the values derived or determined from the aforementioned Darcy and pore-scale studies. We confirmed again that 0.3 is a universal value for [beta] for Darcy's scale sand experiments.

5. Conclusions

We performed a series of downward infiltration experiments and compared the predictions derived from classical GAM and MGAM as well as the values of fitting parameters. In this study, cumulative infiltration monotonously increased with time, and a greater water ponding height showed more cumulative infiltration. This kind of result was also reflected in the wetting front depth. Furthermore, the sand column with smaller grain size had a lower infiltration rate and wetting front velocity, and a greater amount of time was spent on infiltration. Regarding the model simulation, the RMSEs of MGAM were less than those of GAM, indicating that the performance of MGAM with multiparameters for simulating downward infiltration was superior to that of the classical GAM.

Because the classical GAM ignores the dynamic effect on estimating the suction head, the effective contact angle estimated by GAM was systematically greater than that by MGAM, meaning that the equilibrium suction from GAM was systematically less than that from MGAM. Moreover, the contact angle is typically in the form of an advancing one during infiltration. Therefore, the equilibrium contact angle estimated by MGAM was greater than the intrinsic contact angle. In addition, the fitted [??] was proportional to the porosity and inversely proportional to the grain size in our study, indicating that the dynamic effect of capillary pressure increased with decreasing grain size and increasing porosity. Furthermore, the fitted [beta] was approximately 0.3, which is consistent with the values obtained from most of Darcy's and pore-scale studies. This study confirmed that 0.3 can be the universal value of [beta] for Darcy-scale sand experiments.

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

https://doi.org/10.1155/2018/8946948

Conflicts of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper.

Acknowledgments

The authors thank the editors and anonymous referees for their thoughtful comments and suggestions. The authors also want to thank the Ministry of Science and Technology (MOST), Taiwan, for financially supporting this research under Grant MOST 106-2628-M-002-009-MY3 and MOST 106-2923-E-009-001-MY3.

References

 D. A. Barry, J.-Y. Parlange, G. C. Sander, and M. Sivaplan, "A class of exact solutions for Richards' equation," Journal of Hydrology, vol. 142, no. 1-4, pp. 29-46, 1993.

 G. Sposito, "Recent advances associated with soil water in the unsaturated zone," Reviews of Geophysics, vol. 33, no. S2, pp. 1059-1065, 1995.

 W. H. Green and G. A. Ampt, "Studies on soil physics part I--the flow of air and water through soils," Journal of Agricultural Science, vol. 4, pp. 1-24, 1911.

 S. M. Hsu, C.-F. Ni, and P.-F. Hung, "Assessment of three infiltration formulas based on model fitting on Richards equation," Journal of Hydrologic Engineering, vol. 7, no. 5, pp. 373-379, 2002.

 G. Liu, J. R. Craig, and E. D. Soulis, "Applicability of the Green-Ampt infiltration model with shallow boundary conditions," Journal of Hydrologic Engineering, vol. 16, no. 3, pp. 266-273, 2011.

 A. Van den Putte, G. Govers, A. Leys, C. Langhans, W. Clymans, and J. Diels, "Estimating the parameters of the Green-Ampt infiltration equation from rainfall simulation data: why simpler is better," Journal of Hydrology, vol. 476, pp. 332-344, 2013.

 A. D. Feldman, "Hydrologic modeling system HEC-HMS: technical reference manual," in US Army Corps of Engineers, Hydrologic Engineering Center, 2000.

 J. M. Laflen, W. Elliot, D. Flanagan, C. Meyer, and M. Nearing, "WEPP-predicting water erosion using a process-based model," Journal of Soil and Water Conservation, vol. 52, no. 2, pp. 96-102, 1997.

 P. Tuppad, K. R. Douglas-Mankin, T. Lee, R. Srinivasan, and J. G. Arnold, "Soil and water assessment tool (SWAT) hydrologic/water quality model: extended capability and wider adoption," Transactions oftheASABE, vol. 54, no. 5, pp. 1677-1684, 2011.

 F. Bouraoui and T. A. Dillaha, "ANSWERS-2000: runoff and sediment transport model," Journal of Environmental Engineering, vol. 122, no. 6, pp. 493-502, 1996.

 S.-Y. Hsu, V. Huang, S. Woo Park, and M. Hilpert, "Water infiltration into prewetted porous media: dynamic capillary pressure and Green-Ampt modeling," Advances in Water Resources, vol. 106, pp. 60-67, 2017.

 V. T. Chow, D. R. Maidment, and L. W. Mays, Applied Hydrology, McGraw-Hill, 1988.

 S. P. Neuman, "Wetting front pressure head in the infiltration model of Green and Ampt," Water Resources Research, vol. 12, no. 3, pp. 564-566, 1976.

 R. S. Mokady and P. F. Low, "The tension-moisture content relationship under static and dynamic conditions1," Soil Science Society of America Journal, vol. 28, no. 4, 1964.

 G. C. Topp, A. Klute, and D. B. Peters, "Comparison of water content-pressure head data obtained by equilibrium, steady-state, and unsteady-state methods 1," Soil Science Society of America Journal, vol. 31, no. 3, p. 312, 1967.

 E. Diamantopoulos and W. Durner, "Dynamic nonequilibrium of water flow in porous media: a review," Vadose Zone Journal, vol. 11, no. 3, 2012.

 W. C. Lo, C. C. Yang, S. Y. Hsu, C. H. Chen, C. L. Yeh, and M. Hilpert, "The dynamic response of the water retention curve in unsaturated soils during drainage to acoustic excitations," Water Resources Research, vol. 53, no. 1, pp. 712-725, 2017.

 D. M. O'Carroll, T. J. Phelan, and L. M. Abriola, "Exploring dynamic effects in capillary pressure in multistep outflow experiments," Water Resources Research, vol. 41, no. 11, 2005.

 T. Sakaki, D. M. O'Carroll, and T. H. Illangasekare, "Direct quantification of dynamic effects in capillary pressure for drainage-wetting cycles," Vadose Zone Journal, vol. 9, no. 2, 2010.

 F. Stauffer, "Time dependence of the relations between capillary pressure, water content and conductivity during drainage of porous media," in IAHR Symp. on Scale Effects in Porous Media. Thessaloniki Greece 29 Aug.-1. IAHR, Madrid, Spain, 1978.

 L. Zhuang, S. M. Hassanizadeh, C. Z. Qin, and A. de Waal, "Experimental investigation of hysteretic dynamic capillarity effect in unsaturated flow," Water Resources Research, vol. 53, no. 11, pp. 9078-9088, 2017.

 S. M. Hassanizadeh and W. G. Gray, "Mechanics and thermodynamics of multiphase flow in porous media including interphase boundaries," Advances in Water Resources, vol. 13, no. 4, pp. 169-186, 1990.

 T. Annaka and S. Hanayama, "Dynamic water-entry pressure for initially dry glass beads and sea sand," Vadose Zone Journal, vol. 4, no. 1, pp. 127-133, 2005.

 D. A. DiCarlo, "Stability of gravity-driven multiphase flow in porous media: 40 years of advancements," Water Resources Research, vol. 49, no. 8, pp. 4531-4544, 2013.

 S. L. Geiger and D. S. Durnford, "Infiltration in homogeneous sands and a mechanistic model of unstable flow," Soil Science Society of America Journal, vol. 64, no. 2, 2000.

 E. Pellichero, R. Glantz, M. Burns, D. Mallick, S. Y. Hsu, and M. Hilpert, "Dynamic capillary pressure during water infiltration: experiments and Green-Ampt modeling," Water Resources Research, vol. 48, no. 6, 2012.

 N. Weisbrod, T. McGinnis, M. L. Rockhold, M. R. Niemet, and J. S. Selker, "Effective Darcy-scale contact angles in porous media imbibing solutions of various surface tensions," Water Resources Research, vol. 45, no. 4, 2009.

 D. A. Weitz, J. P. Stokes, R. C. Ball, and A. P. Kushnick, "Dynamic capillary pressure in porous media: origin of the viscous-fingering length scale," Physical Review Letters, vol. 59, no. 26, pp. 2967-2970, 1987.

 S. Y. Hsu and M. Hilpert, "Incorporation of dynamic capillary pressure into the Green-Ampt model for infiltration," Vadose Zone Journal, vol. 10, no. 2, 2011.

 T. Tabuchi, "Infiltration and capillarity in the particle packing," Records Land Reclam Res., vol. 19, pp. 1-121, 1971.

 P. G. D. Gennes, X. Hua, and P. Levinson, "Dynamics of wetting: local contact angles," Journal of Fluid Mechanics, vol. 212, no. 1, p. 55, 1990.

 P. Joos, P. van Remoortere, and M. Bracke, "The kinetics of wetting in a capillary," Journal of Colloid and Interface Science, vol. 136, no. 1, pp. 189-197, 1990.

 T. E. Mumley, C. J. Radke, and M. C. Williams, "Kinetics of liquid/liquid capillary rise," Journal of Colloid and Interface Science, vol. 109, no. 2, pp. 398-412, 1986.

 E. Schaffer and P.-z. Wong, "Contact line dynamics near the pinning threshold: a capillary rise and fall experiment," Physical Review E, vol. 61, no. 5, pp. 5257-5277, 2000.

 G. Camps-Roach, D. M. O'Carroll, T. A. Newson, T. Sakaki, and T. H. Illangasekare, "Experimental investigation of dynamic effects in capillary pressure: grain size dependency and upscaling," Water Resources Research, vol. 46, no. 8, 2010.

 D. A. DiCarlo, "Experimental measurements of saturation overshoot on infiltration," Water Resources Research, vol. 40, no. 4, 2004.

 D. Wildenschild, J. W. Hopmans, and J. Simunek, "Flow rate dependence of soil hydraulic characteristics," Soil Science Society of America Journal, vol. 65, no. 1, 2001.

 M. Hilpert, "Velocity-dependent capillary pressure in theory for variably-saturated liquid infiltration into porous media," Geophysical Research Letters, vol. 39, no. 6, 2012.

 Z. Wang, J. Feyen, and D. E. Elrick, "Prediction of fingering in porous media," Water Resources Research, vol. 34, no. 9, pp. 2183-2190, 1998.

 R. J. Glass, J. Y. Parlange, and T. S. Steenhuis, "Immiscible displacement in porous media: stability analysis of three-dimensional, axisymmetric disturbances with application to gravity-driven wetting front instability," Water Resources Research, vol. 27, no. 8, pp. 1947-1956, 1991.

 M. D. Fredlund, D. G. Fredlund, and G. W. Wilson, "An equation to represent grain-size distribution," Canadian Geotechnical Journal, vol. 37, no. 4, pp. 817-827, 2000.

 G. A. F. Seber and C. J. Wild, Nonlinear Regression, Wiley-Interscience, Hoboken, N.J, 2003.

 S. P. Friedman, "Dynamic contact angle explanation of flow rate-dependent saturation-pressure relationships during transient liquid flow in unsaturated porous media," Journal of Adhesion Science and Technology, vol. 13, no. 12, pp. 1495-1518, 1999.

 P. W. Glover and E. Walker, "Grain-size to effective pore-size transformation derived from electrokinetic theory," Geophysics, vol. 74, no. 1, pp. E17-E29, 2009.

 L. Makkonen, "A thermodynamic model of contact angle hysteresis," The Journal of Chemical Physics, vol. 147, no. 6, article 64703, 2017.

 L. H. Tanner, "The spreading of silicone oil drops on horizontal surfaces," Journal of Physics D, vol. 12, no. 9, pp. 1473-1484, 1979.

Yi-Zhih Tsai (iD), (1) Yu-Tung Liu, (2) Yung-Li Wang, (1) Liang-Cheng Chang, (3) and Shao-Yiu Hsu (1)

(1) Department of Bioenvironmental Systems Engineering, National Taiwan University, Taipei 10617, Taiwan

(2) Department of Water Resources Engineering and Conservation, Feng Chia University, Taichung 40724, Taiwan

(3) Department of Civil Engineering, National Chiao Tung University, Hsinchu 30050, Taiwan

Correspondence should be addressed to Shao-Yiu Hsu; syhsu@ntu.edu.tw

Received 20 April 2018; Accepted 12 July 2018; Published 26 August 2018

Caption: Figure 1: Schematic of the experimental setup for the downward infiltration experiments.

Caption: Figure 2: Grain size distribution curves of sands B35, B60, B150, and B320.

Caption: Figure 3: Cumulative infiltration versus time of the sands: (a) B35, (b) B60, (c) B150, and (d) B320.

Caption: Figure 4: Wetting front depth of sand versus time: (a) B35, (b) B60, (c) B150, and (d) B320. The column length of sands B35, B60, and B150 was 60 cm and of sand B320 was 30 cm.

Caption: Figure 5: Wetting front velocity of sands B35, B60, B150, and B320 under ponding height of: (a) [h.sub.w] = 10 cm, (b) [h.sub.w] = 20 cm, and (c) [h.sub.w] = 40 cm.

Caption: Figure 6: Infiltration simulation of GAM and MGAM of sands B35, B60, B150, and B320.

Caption: Figure 7: Wetting front suction head [s.sub.f] versus grain size. In all sands, [s.sub.f] from MGAM was larger than [s.sub.f] from GAM. The solid curve estimates [s.sub.f] from a static contact angle of 30[degrees], the dotted line represents the best contact angle (63[degrees]) with [s.sub.f] from MGAM, and the dashed-dotted line represents the contact angle = 84[degrees] with [s.sub.f] from GAM.

Caption: Figure 8: MGAM parameter [??] versus (a) grain size and (b) porosity.

Caption: Figure 9: MGAM parameter [beta] versus grain size. The fitted [beta] was approximately 0.3.
```Table 1: Geometric mean grain size, porosity, bulk density,
and saturated hydraulic conductivity for column experiments
with different ponding heights.

Geometric
mean grain        Water
size,           ponding             Bulk
[d.sub.50]       height,            density,
x [10.sup.-4]     [h.sub.w]        [[rho].sub.b]
Materials          (cm)            (cm)         (g x [cm.sup.-3])

10         1.36 [+ or -] 0.006
B35                 451             20         1.36 [+ or -] 0.002
40         1.36 [+ or -] 0.017

10         1.43 [+ or -] 0.003
B60                 348             20         1.42 [+ or -] 0.005
40         1.42 [+ or -] 0.008

10         1.38 [+ or -] 0.004
B150                138             20         1.38 [+ or -] 0.007
40         1.39 [+ or -] 0.000

10         1.36 [+ or -] 0.004
B320                13              20         1.34 [+ or -] 0.013
40         1.33 [+ or -] 0.006

Porosity,          Mobile water content,
[epsilon]              [DELTA][theta]
([cm.sup.3]              ([cm.sup.3]
Materials        x [cm.sup.-3])           x [cm.sup.-3])

0.46 [+ or -] 0.003      0.44 [+ or -] 0.011
B35            0.46 [+ or -] 0.001      0.42 [+ or -] 0.006
0.46 [+ or -] 0.007      0.43 [+ or -] 0.006

0.43 [+ or -] 0.001      0.41 [+ or -] 0.014
B60            0.43 [+ or -] 0.002      0.39 [+ or -] 0.003
0.43 [+ or -] 0.003      0.39 [+ or -] 0.003

0.45 [+ or -] 0.002      0.40 [+ or -] 0.001
B150           0.45 [+ or -] 0.003      0.40 [+ or -] 0.005
0.45 [+ or -] 0.000      0.42 [+ or -] 0.021

0.46 [+ or -] 0.002      0.39 [+ or -] 0.007
B320           0.47 [+ or -] 0.005      0.44 [+ or -] 0.040
0.47 [+ or -] 0.002      0.42 [+ or -] 0.036

Saturated hydraulic
conductivity,
[K.sub.s] x
[10.sup.-2]
Materials       (cm x [s.sup.-1])

19.68 [+ or -] 0.683
B35           18.35 [+ or -] 0.132
18.37 [+ or -] 1.591

6.53 [+ or -] 0.202
B60            6.42 [+ or -] 0.079
6.50 [+ or -] 0.393

1.21 [+ or -] 0.024
B150           1.16 [+ or -] 0.055
1.42 [+ or -] 0.391

0.08 [+ or -] 0.003
B320           0.11 [+ or -] 0.026
0.10 [+ or -] 0.028

Table 2: Model parameters and RMSE of GAM and MGAM for
downward infiltration simulation.

GAM
[h.sub.w]
Materials       (cm)        [S.sub.f] (m)     RMSE     [S.sub.f] (m)

10             0.011        0.007         0.105
B35              20            -0.009        0.011         0.094
40            -0.052        0.009         0.09

10             0.053        0.002         0.112
B60              20             0.025        0.005         0.133
40             0.015        0.007         0.135

10             0.134        0.008         0.326
B150             20             0.147        0.004         0.323
40             0.099        0.006         0.31

10             0.528        0.003         1.063
B320             20             0.467        0.002         1.007
40             0.279        0.009         1.496

MGAM

Materials      [??]     [beta]     d (x [10.sup.-4] cm)     RMSE

86.138     0.305             425             0.004
B35           93.316     0.297             456             0.008
100.741     0.283             442             0.006

85.391     0.339             361             0.003
B60           98.828      0.3              351             0.002
93.037     0.29              361             0.004

97.797     0.288             138             0.003
B150          98.817      0.3              138             0.001
98.585     0.291             137             0.002

95.888     0.335              16             0.003
B320          97.124     0.335              17             0.004
109.98     0.305              13             0.002
```
Title Annotation: Printer friendly Cite/link Email Feedback Research Article Tsai, Yi-Zhih; Liu, Yu-Tung; Wang, Yung-Li; Chang, Liang-Cheng; Hsu, Shao-Yiu Geofluids 9TAIW Jan 1, 2018 7096 Study on the Stability of the Coal Seam Floor above a Confined Aquifer Using the Structural System Reliability Method. Mechanism of Permeability Evolution for Reservoir Sandstone with Different Physical Properties. Porosity