A Percolation ‐ Based Approach to Scaling Infiltration and Evapotranspiration

Optimal flow paths obtained from percolation theory provide a powerful tool that can be used to characterize properties associated with flow such as soil hydraulic conductivity, as well as other properties influenced by flow connectivity and topology. A recently proposed scaling theory for vegetation growth appeals to the tortuosity of optimal paths from percolation theory to define the spatio-temporal scaling of the root radial extent (or, equivalently, plant height). Root radial extent measures the maximum horizontal distance between a plant shoot and the root tips. We apply here the same scaling relationship to unsteady (horizontal) flow associated with plant transpiration. The pore-scale travel time is generated from the maximum flow rate under saturated conditions and a typical pore size. At the field-scale, the characteristic time is interpreted as the growing season duration, and the characteristic length is derived from the measured evapotranspiration in that period. We show that the two scaling results are equivalent, and they are each in accord with observed vegetation growth limits, as well as with actual limiting transpiration values. While the conceptual approach addresses transpiration, most accessed data are for evapotranspiration. The equivalence of the two scaling approaches suggests that, if horizontal flow is the dominant pathway in plant transpiration, horizontal unsteady flow follows the same scaling relationship as root growth. Then, we propose a corresponding scaling relationship to vertical infiltration, a hypothesis which is amenable to testing using infiltration results of Sharma and co-authors. This alternate treatment of unsteady vertical flow may be an effective alternative to the commonly applied method based on the diffusion of water over a continuum as governed by Richards’ equation.


Introduction
How can one predict the complex time-dependence of the cumulative infiltration into a soil?Why do estimates of evapotranspiration (ET) at different spatial scales differ?What do infiltration and ET have in common and how do they differ?An initial attempt to answer these questions within the context of percolation theory is given here.
When rain falls upon an unsaturated porous medium such as dry soil, some runs off, some evaporates from the terrestrial surface, and some infiltrates into the medium.Of the latter, a fraction is taken up by the plants through the roots and evaporated from the stomata (a process termed transpiration), and the remainder continues downward to the water table (net infiltration or deep percolation).While the total return of water to the atmosphere includes evaporation and transpiration, known as evapotranspiration (quantities often measured together), our focus is on the plant pathway.Both transpiration and infiltration are instrumental in the drawdown of CO 2 from the atmosphere, the former through photosynthesis, the latter through chemical weathering.In each case, water is critical for the chemistry: in transpiration, for transporting the CO 2 from the atmosphere, plant roots and decaying organic material; in infiltration, as a reagent and as an agent to bring nutrients to plant roots, and thus facilitate nutrient uptake.
The fact that water carries reactive solutes or nutrients with it is relevant to the processes of weathering and vegetation growth, among others.The mathematical treatment of reactive solute transport is a possible source of techniques to address unsteady flow, since the majority of solute transport appears to be non-Gaussian1 [1].Such non-Gaussian transport is generally inconsistent with the Gaussian transport predicted by the Advection-Dispersion Equation (ADE), a formulation which is echoed in the treatment of unsteady flow using the Richards' equation [2].Thus, particularly if infiltration experiments do not prove consistent with the Richards' equation, there may be a motivation to address such discrepancies within the framework of non-Gaussian transport theories.While a wide range of mathematical techniques have been applied to generate non-Gaussian transport, including the Continuous Time Random Walk, or CTRW (see recent review in [3]), to our knowledge, the only technique that has also been adapted to the problem of variable infiltration is the fractional derivative, familiar from Fractional Advection Dispersion Equation (FADE) applications [4]; see [5] for application of the FADE to infiltration.
Scaling laws for solute transport based on percolation concepts, originally proposed by [6] and followed up by [7,8], have recently been applied to the growth and development of plants [9,10], as well as for chemical weathering (see e.g., [11,12]).These scaling laws treat the advective spreading of solutes into a medium that is originally fully filled with pure water, or is in the process of saturation.It is the effective "compressibility" of solutes (changing the volume of a plume as it spreads) that distinguishes solute transport from water flow under steady flow conditions, and leads to the significant non-Gaussian characteristics of the transport [11].In unsaturated, unsteady flow, water spreads (imbibes) into a medium that is originally filled with air, in contrast to steady flow, for which the medium water content does not change.In both cases, the mass advected into the medium is presumed to preferentially follow paths of optimal flow as described in percolation theory, making a predictive approach based on analogy attractive.Importantly, both unsteady flow and solute transport are treated traditionally (in continuum mechanics) as diffusive processes [13]; however, advection of solutes into such a medium is distinctively different from diffusion [12].The analogy between water infiltration and solute transport in terms of the concept of compressibility thus suggests that the same scaling framework that works for chemical weathering and plant growth can be applied towards unsteady water infiltration.
Therefore, the main objectives of this study are (1) to develop theoretical scaling laws for infiltration and transpiration processes in soils under unsteady flow conditions (e.g., the optimal transport distances and the most likely transport time, respectively); and (2) to compare the proposed models with experimental data.Below, we review the scaling theory of solute transport, followed by a brief discussion of scaling arguments linking vegetation growth [9] and soil formation [14].

Theory
In the following, we overview the classical theory of unsteady flow or solute transport only, as their results can be contrasted with the predictions of percolation theory.Since percolation theory (see e.g., [15,16]) is critical to the hypotheses tested in this manuscript, we turn to that subject first.Percolation theory comes in the three variants: site, bond and continuum (corresponding to volume fraction).Following [17,18], Hunt [19] has advocated using percolation formulas adapted to continuum percolation problems for predicting flow and transport phenomena of porous media.The reason is that, with the wide range of pore sizes and local coordination numbers, it is difficult to place such media on any kind of regular grid.Nevertheless, Sahimi [20] showed that wetting, in which water must occupy an entire pore, is a site percolation problem in principle, while drying, for which water must exit through a pore neck, should be conceptualized as a bond percolation problem.Either is more properly treated as an invasion percolation problem than as a random percolation problem while the fluid, either air or water, is incompressible at the pore scale.Such pore-scale incompressibility allows fluid to become trapped in pores rather than displaced, when a second fluid invades.Thus, wetting is invasion site percolation with trapping, drying is invasion bond percolation with trapping, and fully saturated conditions are termed random percolation.This classification scheme is important if one seeks relevant information on various percolation scaling predictions, as we do below.
One useful application from percolation theory is critical path analysis (CPA) developed to find dominant paths of charge or fluid flow through heterogeneous media [19,21,22].The concept is most easily understood in a network model constructed of pore bodies connected through pore throats (see e.g., [23]), for which the difficulty of flow is localized along the connections between pores and discretized accordingly.Within the percolation theory framework, sites and bonds are equivalent to bodies and throats in a network of pores, respectively.The strategy of CPA is to apply the percolation threshold-the minimum fraction required to form the sample-spanning cluster and let the system percolate-to the spread of local resistance values in order to find the most resistive element on the most highly conductive (or permeable) interconnected path through the system.Sahimi [16] observed that solute transport along such optimal paths as described in percolation theory should follow the same scaling laws as solute transport through actual percolation clusters.Here, we extend that conjecture to unsteady flow as well.

Percolation Scaling Theory for Solute Transport
Percolation theory provides a promising framework that has been widely applied to study flow and transport in interconnected networks and media.Percolation theory combines microscopic characteristics of a medium with its topological properties to determine macroscopic quantities, such as conductivity, permeability, and dispersion.A full theoretical treatment of the scaling of solute transport is beyond the scope of the present work, as we are interested in distances and times such as minimum and optimal transport distances and the most likely transport time as a function of system size, rather than arrival time distributions.The most likely transport time has been found to generate transport velocities and fluxes reasonably accurately for up to five orders of magnitude of time [12].For a comprehensive review of percolation theory applications on fluid flow and solute transport, see [16,24,25].
Percolation theory applies most simply to porous media for which flow or conduction along bonds (pore throats) between nodes (pore bodies) is either possible or impossible [26].Here, it delivers relatively simple scaling laws for, e.g., conduction or flow properties.Such binary media can be represented by networks with resistances of value 1 or infinity, but are not characteristic of typical natural porous media, since flow is strongly pore-size dependent and pore size typically spans several orders of magnitude.Instead of bond percolation, continuum percolation has been used [17][18][19] to address flow in more realistic pore networks, with irregular lattices made of pores of different sizes and shapes.In media with a wide range of pore sizes, critical path analysis [21,22] utilizes the percolation threshold value to generate a percolating sub-network composed of all resistances less than or equal to a "critical" value, a technique which can provide effective flow parameters such as conductivity or permeability [19].Thus, a sub-network composed of all the resistances less than or equal to a "critical" value just percolates.The importance in the present context is that Sahimi argued already more than 20 years ago [16] that solute transport in real porous media with wide ranges of pore sizes should also follow percolation scaling results because of the concentration of the flow along such critical paths.Because critical path analysis is applicable to virtually all natural soils [27], percolation theory is equally relevant to solute transport; thus, we begin by listing solute transport results from percolation theory.Lee et al. [6] addressed the spatiotemporal scaling of solutes transported across a two-dimensional system at the percolation threshold.The authors found that for a line source of solute in a system of length x, the minimum path length, L, for particles exiting on the opposite face scales with x as L ∼ x D min , with a most likely time of t ∼ x D b , where D min and D b are the chemical paths exponent and the fractal dimensionality of the percolation backbone, respectively.However, when the medium is highly heterogeneous, e.g., with a wide range of pore sizes, the length L of the optimal paths that provide the most rapid transport across the system is proportional to the system length to a different exponent, D opt , rather than to D min [28,29].D opt is the so-called optimal paths exponent, and is determined for the heterogeneous limit [30].The values of the exponents determined by [30] are mostly accurate to five significant figures, but the optimal paths exponent only to three.The most important distinctions in porous media that affect the exponents D b and D min relate to the connectivity of the paths and the specific percolation conditions: (i) imbibition (wetting) conditions correspond to site invasion percolation with trapping; (ii) drying conditions to invasion percolation with trapping; and (iii) saturated flow to random percolation [20].Following [6] the minimum path length L is (using D min from [30]): where x 0 is the bond length, and 2D and 3D refer to two and three dimensional systems.In the case of porous media, x 0 corresponds to the inter-pore distance or a typical grain size.Numerical prefactors of order unity are less well constrained, hence omitted from Equations ( 1) and ( 2).The value of the exponent D opt in 2D (for all conditions of saturation), 1.21, is similar to the value of D min for 2D imbibition or drainage, 1.21, and the optimal path lengths L scale the same as the shortest paths for imbibition or drainage, L = x 0 (x/x 0 ) 1.21 .In 3D D opt = 1.46.Similarly, the most likely travel time t is [30] We note that the exponent for 3D imbibition, 1.861, differs from the value for saturated conditions, 1.87 (Equation ( 3)).As we discuss in Section 6, such a small discrepancy actually leads to detectable differences in scaling of the cumulative infiltration.Although the conditions we address are (1) for vertical infiltration, wetting, 3D flow; and (2) for transpiration, desaturating 2D flow, we include the entire suite of values of these exponents in case conditions addressed in other investigations are different.In particular, swelling soils may develop mudcracks on drying, providing preferential flow topologies for imbibition that are constrained to two-dimensional surfaces.
In Equations ( 1)-( 4), t 0 is the time for fluid to traverse a distance x 0 , such that v 0 = x 0 /t 0 is the fluid flow rate at the pore scale.The contrast between the exponent 1.87 for the most likely time and 1.37 for the shortest distance (for 3D saturated conditions) helps explain why paths that do not appear to be extraordinarily tortuous (20 m in a 1 m column) can still retard the arrival of solutes by an enormous factor [11] (over 1000).Since we are interested in (1) vegetation growth and transpiration in 2D (optimal paths exponent); and (2) unsteady flow in 3D (for the wetting conditions in Philip's infiltration theory), we select here the following exponent values (from Equations ( 1)-( 4)): (i) D opt (2D) =1.21 for vegetation growth, and (ii) D b (3D imbibition) = 1.861 for vertical infiltration.These two exponents are essentially the same as the two values assumed relevant [9,12] for soil formation (1.87) and vegetation growth (1.21), although as a result, our analogy to soil formation and chemical weathering is not quite precise.
In [12], for example, soil was assumed to form under either saturated or saturating conditions, but only the exponent 1.87 was used.Given all the uncertainties in soil formation processes operating at time scales of up to tens of millions of years, ignoring the possibility that the appropriate exponent might, for some fraction of the time, be 1.861, seemed justified.However, in unsteady infiltration experiments, of much shorter duration and with much greater control, the exponent is constrained conceptually to be 1.861.

Experimental Data
Experimental data used for the analysis in this study are from three databases: (1) Falster et al. [31]; (2) Hunt [9] and (3) Sharma et al. [32].The Biometric and Allometric Database (BAAD), in [31] includes allometric data that consist of height and age for plantation species as well as from plants growing under natural conditions.The data points from plants grown in natural conditions, i.e., 5539 out of 6650 are used to compare with the percolation theory scaling relationship.The data given in [9], the Hunt dataset, were collected from well over 100 sources, but the natural plant data were selected to include only faster growing species, notably the three sequoia species, many eucalypts, eastern cottonwoods, douglas firs, tropical dipterocarps, etc. Crop data were scarce, as most studies address total biomass, while soil data corresponded mainly to depths to the Bw horizon.
The Sharma et al. database [32] consists of experimentally measured cumulative infiltration via double-ring infiltrometers at 26 sites on a watershed, called R-5, close to Chickasha, Oklahoma.R-5 is native grassland dominated by silt loam soil texture with a slope of near 3% and an area of 9.6 ha.

Relevance of Percolation Scaling to Vegetation Growth and Soil Formation
The primary purposes in this section are: (1) to show the basis of the pore-scale link between vegetation growth and soil formation; (2) provide a basis for upscaling the pore-scale water flow rate controlling vegetation growth to growing season transpiration rates, as addressed in the next section.This provides both a context and a method for the following topics.
Hunt [9,14] showed that the relationship with D b = 1.87 (for 3D flow connectivity) and x 0 a typical grain diameter predicts soil depths over an extensive range of time scales (a decade to 100 Ma), where v 0 = x 0 /t 0 was taken to be the pore-scale infiltration rate.Equation ( 5) provides a soil formation rate of R s = dx/dt = 0.53v 0 (x/x 0 ) −0.87 or 0.53v 0 (t/t 0 ) −0.46 .Furthermore, ref. [12] showed that the related scaling of solute transport rates (from the velocity of the centroid of the solute spatial distribution) could predict soil formation rates at time scales of up to 60 Ma, and that that the predicted proportionality of R s to v 0 is in accordance with the observed four orders of magnitude variability in R s when precipitation also varied over four orders of magnitude (from the Atacama desert to the New Zealand Alps).Thus, rate of flow of water through the surface medium is responsible for much of the variability in soil formation rates.
Considering vegetation linear extent (root radial extent, RRE, or plant height), Hunt and Manzoni [9,10] showed that with D opt = 1.21 (2D) and x 0 taken as a typical pore diameter (~10, the geometric mean plant xylem diameter, Watt [33] hold for time scales ranging from minutes to 100 kyr.In this case, x is either the plant height, or the root radial extent, which have been shown to be closely equivalent over time scales from a growing season to about 40 years [10,[34][35][36][37]. Once again, v 0 = x 0 /t 0 is the pore-scale flow rate, ranging here from about 24 cm/year to 20 m/year.This is a similar, though narrower range of flow rates than those relevant for soils under typical conditions i.e., about 2.0 cm/year to 20 m/year.The soil depths and the linear extent of natural vegetation computed from Equations ( 5) and ( 6) are given in Figure 1.The vegetation data incorporate a wide range of data sets: at time scales from minutes to days, these are root tip extension rates, for weeks to hundreds of years, they are plant heights (which are a proxy for root radial extent over most of this range), and for thousands of years up, the measurements are subsurface extents of clonal organisms, such as Posidonia oceanica and various tree species, or fungi.We note that the leaf area index could also be used to characterize plant growth [38], instead of plant height.At the longest time scales, the soil data are mostly paleosols and laterites.At the shortest time scales, the soil data include primarily soil formation after disturbances, such as landslides, glacial retreat, or tree-throw.For typical pore and particle size of 10 µm and 30 µm, respectively, the presented pore scale velocities bound the data.We point out that since larger pores have a smaller retention capacity as they are easily drained, their contribution to plant growth is smaller [33].We will apply essentially the same Equations ( 5) and ( 6) to vertical infiltration and horizontal flow associated with transpiration, respectively.

Potential Relationship with Unsteady Flow
At spatial and temporal scales larger than the (pore) scales considered above, we hypothesize that the ratio of a transpiration depth to a transpiration time scale can take the place of the pore-scale flow velocity v 0 = x 0 /t 0 .Thus, x 0 /t 0 → x g /t g , where x g is the transpiration during a growing season (measured typically as a depth) and t g is the length of the growing season.Using this hypothesis, we rewrite the scaling for vegetation growth in Equation ( 6) as Furthermore, using total ET depths x g (instead of transpiration data, which is much less abundant) that range from 0.02 m (the Namibian desert at the dry end of the scale [39]) to 1.65 m (tropical rainforests and savannahs [40]) for a growing season of t g = 0.5 year, provides x = 0.02 m (t/180 days) 0.83 and x = 1.65 m (t/180 days) 0.83 , for a minimum and a maximum transpiration, or, equivalently, by Equation (7), for a minimum and maximum growing season plant height, respectively.As shown in Figure 2, these values of x (t g ) bound the world's plant heights [31] at a time corresponding to the length of the growing season, meaning that Equation (7) generates essentially identical scaling predictions of plant growth rates as depicted in Figure 1.
The data of [31] Falster et al. ( 2015) is a compilation performed by many authors over many years of dozens of studies worldwide.Note that regnans minimum and maximum predictions are generated from the range in adult Eucalyptus regnans heights (a factor 22) across a precipitation gradient [41] together with the demonstration [10], that all Eucalyptus growth rates in the study of [42] followed the scaling of Equation ( 6), but with variations in the t 0 values consistent with a gradation of flow rates across their study sites.Such a close correspondence between predicted upper and lower vegetation growth limits and the observed range of actual vegetation growth rates in Figure 2 (the independent meta-data set of [31]) implies the possibility of applying an analogous upscaling to the horizontal soil moisture transport associated with the plant transpiration.Further, it also suggests a close relationship between transpiration "depths" and plant heights, which is also indicated by examining record breaking crop heights [9].For amaranth, sunflowers, hemp, and corn, these record heights all are approximately 10 m in a growing season, and for tomatoes, 20 m in a year.Such rates correspond fairly closely to some classic literature source values for the saturated hydraulic conductivity, K S .In particular, a geometric mean K S value for soils is roughly 1 µm/s [43,44] and typical saturated subsurface flow rates are of the same order of magnitude [45].An amount of 1 µm/s corresponds to a flow distance of about 32 m in 1 year.Thus, such high transpiration rates required for growth rates on the order of 20 m/year would require the existence of very nearly saturated conditions much of the time in a typical soil.This approximate correspondence suggests that the soil saturated hydraulic conductivity may generate a rough maximum for plant growth rates, if the water (and nutrient) supply is adequate.We note that in the above we consider saturation as optimal conditions in terms of conducting water and nutrients; unsaturated conditions are optimal in terms of aerobic root metabolism.

Figure 1.
Vegetation growth and soil formation as a function of time.Scaling predictions for vegetation growth, both natural (green) and intensively managed (blue), using maximum pore scale flow rates of 20 m/year, and minimum flow rates of 0.24 m/year.We assume pore and particle sizes of 10 μm and 30 μm, respectively, in accord with typical silt particle sizes and the usual assumption that pore diameters are about 30% that of particles [46].The vegetation data are divided into crops [9], for which nitrogen and other nutrients were supplied; "fast plants" [9,10] and BAAD [31].The sources for the soil database are referenced in [9] and [10].The same maximum flow rate prediction (marked "Flow"), but using vertical solute transport scaling ("Physical Transport"), generates the upper bound on soil production as on plant growth ("Biological Transport") using the horizontal optimal paths scaling."Phys Min" and "Bio Min" are the corresponding predictions at the given lower flow rates.The lower bound for soil formation is consistent with a smaller flow velocity than for vegetation growth, in accord with the fact that soil formation, albeit at a very slow rate, can occur in regions where precipitation is insufficient to support plant growth (e.g., the Atacama desert). -

Log [x (m)]
Figure 1.Vegetation and soil formation as a function of time.Scaling predictions for vegetation growth, both natural (green) and intensively managed (blue), using maximum pore scale flow rates of 20 m/year, and minimum flow rates of 0.24 m/year.We assume pore and particle sizes of 10 µm and 30 µm, respectively, in accord with typical silt particle sizes and the usual assumption that pore diameters are about 30% that of particles [46].The vegetation data are divided into crops [9], for which nitrogen and other nutrients were supplied; "fast plants" [9,10] and BAAD [31].The sources for the soil database are referenced in [9] and [10].The same maximum flow rate prediction (marked "Flow"), but using vertical solute transport scaling ("Physical Transport"), generates the upper bound on soil production as on plant growth ("Biological Transport") using the horizontal optimal paths scaling."Phys Min" and "Bio Min" are the corresponding predictions at the given lower flow rates.The lower bound for soil formation is consistent with a smaller flow velocity than for vegetation growth, in accord with the fact that soil formation, albeit at a very slow rate, can occur in regions where precipitation is insufficient to support plant growth (e.g., the Atacama desert).[31] is compared with our predicted upper and lower limits on growth rates using the pore-scale flow rates discussed in the text."Regnans maximum" is the maximum plant height of E. regnans computed from Equation ( 7) with 1 µm/s flow rate [9], a value which corresponds to the prediction from the analogous equation at larger scales using a corresponding ET of ~1650 mm in a growing season."Regnans minimum" is generated from "regnans maximum" by multiplying by a factor 1/20, in general accord with the variability in precipitation to pan evaporation ratio determined in [41]."Overall minimum" is obtained by multiplying "regnans maximum" by the ratio of smallest to largest ET values in the BAAD, 1/82.5.The hydraulic limit is postulated as a cavitation limit [47].Note that the minimum and maximum scaling relationships could be equally expressed in terms of the limiting growing season transpiration values, 20 mm at the dry end of the spectrum [39] and 1650 mm at the wet end [40], which are plotted on the graph at the time of six months, a typical growing season.The correspondence between plant height and growing season transpiration suggests the viability of growth models in terms of growing season transpiration values.

Relevance of Percolation Scaling to Infiltration
Philip infiltration theory [48,49] represents cumulative infiltration I as having two terms, one a steady-state component, proportional to the time, the second a transient component proportional to the square root of time.Existence of higher order terms is also acknowledged, although these terms are mostly neglected.The two-parameter Philip's equation is Here, S is sorptivity, and A is a constant coefficient that depends on soil properties and initial water content.In the soil physics literature, in Equation ( 8) the first term represents the effect of gravity forces and the second term is due to sorptive forces.In Philip's theory [50], A is proportional to the saturated hydraulic conductivity, K S , (1/3 ≤ A/K S ≤ 2/3), while S is proportional to the square root of K S , suggesting a dependence of S ~A0.5 , as noted (and also investigated) by [32].However, as we show below, our treatment yields a different proportionality, S ~A0.77 .
We propose that the transient term is not governed by diffusion, but by the same scaling relationship that describes solute transport under saturated (or saturating) conditions.This transient term is suggested to correspond conceptually to the solute that is carried along with the water, with the steady-state infiltration term thus representing the advecting fluid.The cumulative infiltration is, therefore, given by (9) in which the first term (x 0 /t 0 )t represents the steady state infiltration, and the second term (x 0 /t 1/D b 0 )t 1/D b the scaling power law from percolation theory.Note that in Equation ( 9) we have explicitly used the value of D b = 1.861 for conditions of imbibition, rather than 1.87, appropriate for saturated conditions.Although the two resulting exponents (1/D b ) of 0.5348 (saturated) and 0.537 (imbibition) are quite similar, since the medium is unsaturated by definition, one should use the latter, as in Equation ( 9).
Consider our Equation ( 9) where t 0 = x 0 /v 0 , and v 0 is proportional to the saturated hydraulic conductivity, K S .Most formulations (see e.g., [51]) of K S lead to a proportionality to x 0 2 (as does the critical path analysis of [52]), such that x 0 /v 0 is proportional to x 0 −1 .This makes A proportional to x 0 2 (or K S ) but S proportional to x 0 1.53 , or K S 0.77 , suggesting a different proportionality S ≈ A 0.77 , than that in Philip's theory, S ≈ A 0.5 .
Is there an analysis that can confirm the relevance of the exponents, 1 and 0.5?Could it distinguish between 0.5 and 0.54 (0.537 rounded up) in the transient term?The analysis proposed by [32] can reasonably confirm the approximate form of the equation, but, although the analysis does suggest that the transient-term exponent is not precisely 0.5, it cannot establish that it is actually 0.54.However, Figure 1b in [53] indicates that a much better collapse to a single curve results when infiltration data are scaled according to xt −0.55 than for xt −0.5 , meaning that an exponent value of 0.537 is a better approximation than 0.5.Although the values of the exponents in [30] are exceptionally well constrained, it is known that this percolation scaling is precise only in the limit of large system sizes [15].Systems may be considered large at upwards of 100 bonds on a side, which for typical soil particle sizes in the tens of microns means length scales of millimeters and up [25].
In the first step of their analysis (which does not address the specific form of A and S), Sharma et al. [32] consider the functions and These functions transform Equation ( 8) to β = τ + τ 0.5 , independent of the specific values of A and S, as stated by [32], but generate a slightly different result if I is given by Equation (9).In particular, we have Perfect scaling, as predicted by Philip infiltration, is no longer quite developed, and a weak source of scatter is introduced (the factor t 0 0.04 ), generally in conformance with the results.However, in order to determine with certainty whether the residual scatter can be attributed to the additional time factor produced by our particular theoretical description, it would be necessary to know the individual values of t 0 .Unfortunately, these values are not obtainable from the published data.It is clear, however, that the deviation must be due to a discrepancy in the power, because the formulation of the scaling functions [32] eliminates any scatter due to variability in the coefficients.For a first estimate of how much uncertainty is introduced by the factor t 0 0.04 , note that t 0 is proportional to τ −1.06 .From this, we estimate the variability in t 0 using the largest and smallest τ values in the data set, 0.271 and 5.8 × 10 −4 , respectively, providing minimum and maximum multiplicative factors of 1.05 and 1.32.Multiplication of the second term by these values gives the lower and upper bounds, respectively, on the theoretical prediction in Figure 3.It should be noted finally that additional terms in the infiltration treatment due to Philip would also generate deviations in scaling.Interestingly, however, our lowest order theoretical treatment appears to generate just the right magnitude of scaling deviations to be compatible with experiment.An important point relates to the limits of applicability of percolation techniques.For very short time scales, or equivalently, very small systems, fluctuations associated with the statistical variability introduced by small system sizes will produce a wide variability in results, and, typically, overestimate the infiltration (Figure 3), while underestimating the mean [19].In accord with these observations, we reproduce in Figure 4 a bi-logarithmic plot of the data and predictions given in Figure 3.The results are in accord with the statistical interpretation, and with previous numerical results [19] for the finite-size scaling of the hydraulic conductivity, K, which similarly had a few values of K that are much higher than predicted, and a larger number that are much lower (Figure 4).This general conclusion is also in accord with the treatment of solute transport [54], for which the simple scaling result underestimated the variability in solute dispersion at shorter distance scales, requiring calculation of the entire solute arrival time distribution for quantitative agreement with the experiment.If infiltration is analogous, the corresponding treatment of unsteady flow could introduce pore-size distribution effects into the short time (or distance) predictions.
Finally, to confirm the accuracy of our above postulate, we test in Figure 5 the prediction in [32] regarding the relationship between A and S. Fitting of the data into a power-law form yields the relationship S ~A0.727 , which is in reasonably good agreement with our prediction in the discussion following Equation (8) (0.77, 6% difference), but deviates significantly (by 45%) from the value predicted by the classical Philip's infiltration theory, 0.5.
Water 2017, 9, 104 10 of 14 finally that additional terms in the infiltration treatment due to Philip would also generate deviations in scaling.Interestingly, however, our lowest order theoretical treatment appears to generate just the right magnitude of scaling deviations to be compatible with experiment.An important point relates to the limits of applicability of percolation techniques.For very short time scales, or equivalently, very small systems, fluctuations associated with the statistical variability introduced by small system sizes will produce a wide variability in results, and, typically, overestimate the infiltration (Figure 3), while underestimating the mean [19].In accord with these observations, we reproduce in Figure 4 a bi-logarithmic plot of the data and predictions given in Figure 3.The results are in accord with the statistical interpretation, and with previous numerical results [19] for the finite-size scaling of the hydraulic conductivity, K, which similarly had a few values of K that are much higher than predicted, and a larger number that are much lower (Figure 4).This general conclusion is also in accord with the treatment of solute transport [54], for which the simple scaling result underestimated the variability in solute dispersion at shorter distance scales, requiring calculation of the entire solute arrival time distribution for quantitative agreement with the experiment.If infiltration is analogous, the corresponding treatment of unsteady flow could introduce pore-size distribution effects into the short time (or distance) predictions.
Finally, to confirm the accuracy of our above postulate, we test in Figure 5 the prediction in [32] regarding the relationship between A and S. Fitting of the data into a power-law form yields the relationship S ~ A 0.727 , which is in reasonably good agreement with our prediction in the discussion following Equation (8) (0.77, 6% difference), but deviates significantly (by 45%) from the value predicted by the classical Philip's infiltration theory, 0.5.10) vs. scaled time , calculated from Equation (11), for 26 experiments from [32] (digitized from their Figure 2b).The infiltration data represent field measurements using double-ring infiltrometers and 26 infiltration tests carried out on a 9.6-ha grassland watershed (called R-5) near Chickasha, Oklahoma.Given that t0 is proportional to τ −1.06 , we estimated the variability in t0 in Equation ( 12) using the largest (0.271) and smallest 5.8 × 10 -4 values of τ in the dataset.This provides minimum and maximum multiplicative factors of 1.05 and 1.32.Multiplication of the second term in Equation ( 12) by these values gives the lower and upper bounds, respectively for the theoretical predictions.The blue and red curves are the minimum and maximum expected from the theory in this study.Philip's infiltration theory allows no scatter, and would, to within the resolution of the graph, coincide with the blue curve.11), for 26 experiments from [32] (digitized from their Figure 2b).The infiltration data represent field measurements using double-ring infiltrometers and 26 infiltration tests carried out on a 9.6-ha grassland watershed (called R-5) near Chickasha, Oklahoma.Given that t 0 is proportional to τ −1.06 , we estimated the variability in t 0 in Equation ( 12) using the largest (0.271) and smallest 5.8 × 10 −4 values of τ in the dataset.This provides minimum and maximum multiplicative factors of 1.05 and 1.32.Multiplication of the second term in Equation ( 12) by these values gives the lower and upper bounds, respectively for the theoretical predictions.The blue and red curves are the minimum and maximum expected from the theory in this study.Philip's infiltration theory allows no scatter, and would, to within the resolution of the graph, coincide with the blue curve.Replotting the data and predictions given in Figure 3 bi-logarithmically (same color coding) demonstrates that the data agree with the statistical nature of percolation.In particular, as times and corresponding length scales increase relative to pore-scale times, the statistical constraints from percolation theory become more restrictive, and large fluctuations related to small sample size become less important.[32], digitized from their Figure 5, for the relationship between S and A, the coefficients of the transient and steady-state terms of the infiltration equation, respectively.Note that the exponent observed from the data, 0.727 (black line), is in much closer agreement with our prediction using percolation theory, 0.77 (in red), than with the value in Philip's infiltration theory, 0.5 (in blue).

Conclusions
Existing scaling predictions from percolation theory of vertical solute transport, chemical weathering and soil formation [9,11,12,14]  Replotting the data and predictions given in Figure 3 bi-logarithmically (same color coding) demonstrates that the data agree with the statistical nature of percolation.In particular, as times and corresponding length scales increase relative to pore-scale times, the statistical constraints from percolation theory become more restrictive, and large fluctuations related to small sample size become less important.
Water 2017, 9, 104 11 of 14 Replotting the data and predictions given in Figure 3 bi-logarithmically (same color coding) demonstrates that the data agree with the statistical nature of percolation.In particular, as times and corresponding length scales increase relative to pore-scale times, the statistical constraints from percolation theory become more restrictive, and large fluctuations related to small sample size become less important.[32], digitized from their Figure 5, for the relationship between S and A, the coefficients of the transient and steady-state terms of the infiltration equation, respectively.Note that the exponent observed from the data, 0.727 (black line), is in much closer agreement with our prediction using percolation theory, 0.77 (in red), than with the value in Philip's infiltration theory, 0.5 (in blue).

Conclusions
Existing scaling predictions from percolation theory of vertical solute transport, chemical weathering and soil formation [9,11,12,14] [32], digitized from their Figure 5, for the relationship between S and A, the coefficients of the transient and steady-state terms of the infiltration equation, respectively.Note that the exponent observed from the data, 0.727 (black line), is in much closer agreement with our prediction using percolation theory, 0.77 (in red), than with the value in Philip's infiltration theory, 0.5 (in blue).

Conclusions
Existing scaling predictions from percolation theory of vertical solute transport, chemical weathering and soil formation [9,11,12,14] invoke 3D imbibition (or saturated) fractal dimensionality of the backbone, relating time and distance of solute transport.Here, by comparison with an analysis of a series of experiments by Sharma et al. [32], it was shown that essentially the same scaling describes vertical infiltration.For horizontal infiltration, we adopt scaling arguments for the growth of vegetation based on the ability to find essential nutrients to the scaling of transpiration.As a consequence, we are able to express plant growth extent directly as the transpiration depth.
Overall, the above arguments support the use of percolation scaling formulations for unsteady flow in heterogeneous media such as soils, both horizontally (along layers), and vertically (through them).The implication of our approach, where percolation scaling is based on advection rather than diffusion, is that the Richards' equation (which takes the form of a diffusive process) is not an adequate description for scaling of unsteady flow.These two results for flow bear important implications for the formation of soil and the growth of the vegetation within the soil.

Figure 3 .
Figure 3. Scaled parameter , determined from Equation(10) vs. scaled time , calculated from Equation(11), for 26 experiments from[32] (digitized from their Figure2b).The infiltration data represent field measurements using double-ring infiltrometers and 26 infiltration tests carried out on a 9.6-ha grassland watershed (called R-5) near Chickasha, Oklahoma.Given that t0 is proportional to τ −1.06 , we estimated the variability in t0 in Equation (12) using the largest (0.271) and smallest 5.8 × 10 -4 values of τ in the dataset.This provides minimum and maximum multiplicative factors of 1.05 and 1.32.Multiplication of the second term in Equation (12) by these values gives the lower and upper bounds, respectively for the theoretical predictions.The blue and red curves are the minimum and maximum expected from the theory in this study.Philip's infiltration theory allows no scatter, and would, to within the resolution of the graph, coincide with the blue curve.

Figure 3 .
Figure 3. Scaled parameter β, determined from Equation(10) vs. scaled time τ, calculated from Equation(11), for 26 experiments from[32] (digitized from their Figure2b).The infiltration data represent field measurements using double-ring infiltrometers and 26 infiltration tests carried out on a 9.6-ha grassland watershed (called R-5) near Chickasha, Oklahoma.Given that t 0 is proportional to τ −1.06 , we estimated the variability in t 0 in Equation (12) using the largest (0.271) and smallest 5.8 × 10 −4 values of τ in the dataset.This provides minimum and maximum multiplicative factors of 1.05 and 1.32.Multiplication of the second term in Equation (12) by these values gives the lower and upper bounds, respectively for the theoretical predictions.The blue and red curves are the minimum and maximum expected from the theory in this study.Philip's infiltration theory allows no scatter, and would, to within the resolution of the graph, coincide with the blue curve.

Figure 4 .
Figure 4.Replotting the data and predictions given in Figure3bi-logarithmically (same color coding) demonstrates that the data agree with the statistical nature of percolation.In particular, as times and corresponding length scales increase relative to pore-scale times, the statistical constraints from percolation theory become more restrictive, and large fluctuations related to small sample size become less important.

Figure 5 .
Figure 5. Sorptivity S as a function of Philip model parameter A. Comparison of the prediction of percolation scaling (red) with the 26 experiments collected by[32], digitized from their Figure5, for the relationship between S and A, the coefficients of the transient and steady-state terms of the infiltration equation, respectively.Note that the exponent observed from the data, 0.727 (black line), is in much closer agreement with our prediction using percolation theory, 0.77 (in red), than with the value in Philip's infiltration theory, 0.5 (in blue).

Figure 4 .
Figure 4.Replotting the data and predictions given in Figure3bi-logarithmically (same color coding) demonstrates that the data agree with the statistical nature of percolation.In particular, as times and corresponding length scales increase relative to pore-scale times, the statistical constraints from percolation theory become more restrictive, and large fluctuations related to small sample size become less important.

Figure 4 .
Figure 4.Replotting the data and predictions given in Figure3bi-logarithmically (same color coding) demonstrates that the data agree with the statistical nature of percolation.In particular, as times and corresponding length scales increase relative to pore-scale times, the statistical constraints from percolation theory become more restrictive, and large fluctuations related to small sample size become less important.

Figure 5 .
Figure 5. Sorptivity S as a function of Philip model parameter A. Comparison of the prediction of percolation scaling (red) with the 26 experiments collected by[32], digitized from their Figure5, for the relationship between S and A, the coefficients of the transient and steady-state terms of the infiltration equation, respectively.Note that the exponent observed from the data, 0.727 (black line), is in much closer agreement with our prediction using percolation theory, 0.77 (in red), than with the value in Philip's infiltration theory, 0.5 (in blue).

Figure 5 .
Figure 5. Sorptivity S as a function of Philip model parameter Comparison of the prediction of percolation scaling (red) with the 26 experiments collected by[32], digitized from their Figure5, for the relationship between S and A, the coefficients of the transient and steady-state terms of the infiltration equation, respectively.Note that the exponent observed from the data, 0.727 (black line), is in much closer agreement with our prediction using percolation theory, 0.77 (in red), than with the value in Philip's infiltration theory, 0.5 (in blue).