Newton's Shear Flow Applied to Infiltration and Drainage in Permeable Media

The paper argues that universal approaches to infiltration and drainage in permeable media that pivot around capillarity and that led to dual porosity, non-equilibrium, or preferential flow need to be replaced by a dual process approach. One process has to account for relatively fast infiltration and drainage based on Newton's shear flow, while the other one is responsible for storage and relatively slow redistribution of soil water by focusing on capillarity. Already Schumacher (1864) postulated two separate processes. However, Buckingham's (1907) and Richards' (1931) apparent universal capillary-based approach to flow and storage of water in soils dominated. The paper introduces the basics of Newton's shear flow in permeable media. It presents experimental support for the four presumptions of (i) sharp wetting shock fronts; (ii) that move with constant velocities; (iii) atmospheric pressure prevails behind the wetting shock front; (iv) laminar flow. It further discusses the scale tolerance of the approach, its relationship to Darcy's (1856) law, and its extension to solute transport.


Introduction
Infiltration is the transgression of liquid water from above the surface of the permeable lithosphere to its interior, while drainage refers to liquid water leaving some of its bulk.
Infiltration and drainage still bear unsolved problems. For instance, Blöschl et al. (2019), in a most thorough, exhaustive, and detailed survey among hundreds of active hydrologists, compiled 23 Unsolved Problems in Hydrology (UPHs). The 7 th UPH asks "Why is most flow [in the unsaturated lithosphere P.G.] preferential across multiple scales and how does such behaviour co-evolve with the critical zone?" The critical zone in hydrology delineates the upper most layer of the lithosphere that is in direct contact with the atmosphere. It typically carries the terrestrial ecosystems and thus simultaneously provides water, air, nutrients, and mechanical support to roots of most terrestrial plant communities. In general, the critical zone is congruent with soil. This contribution presents a solution of the 7 th UPH.
The second section reviews the evolution of infiltration concepts in partially saturated soils since the second half of the 19 th century. The next one summarizes Newton's shear flow applied to flow in permeable media, while the other ones provide support to, and various applications of the approach.

Review of infiltration concepts
In the mid-19th century, there was an increasing interest in flows in saturated soils and similarly permeable media. Hagen, a German hydraulic engineer, and Poiseuille (1846), a French physiologist, independently analyzed laminar flow in thin capillary tubes. Darcy (1856), in the quest of designing a filtration system for the city of Dijon, empirically developed the concept of hydraulic conductivity as proportionality factor of flow's linear dependence on the pressure gradient. Dupuit (1863) expanded Darcy's law to two dimensions as perpendicular and radial flow between two parallel drainage ditches and towards a groundwater well, respectively. Schumacher (1864), a German agronomist, was probably the first who considered capillarity as the cause for simultaneous flows of water and gas in partially water-saturated soils.
He qualitatively compared the rise of wetting fronts in soil columns with the rise of water in capillary-sized glass tubes, and concluded that the wetting fronts rise higher but slower in finer textured soils compared with coarser materials. He also infiltrated water in columns of undisturbed soil and found that infiltration fronts progressed much faster than the rising wetting fronts. He suggested two separate processes for the two flow types: (i) slower capillary rise and (ii) faster infiltration, however, without further dwelling on infiltration. Lawes et al.(1882) concluded from the chemical composition of the drain from large lysimeters at the Rothamsted Research Station that "The drainage water of a soil may thus be of two kinds (1) of rainwater that has passed with but little change in composition down the open channels of the soil; or (2) of the water discharged from the pores of a saturated soil." Lawes et al. (1882) prioritized two separate flow paths to explain the observations. During the second half of the 19 th century irrigation agriculture spread in semi-arid areas and so increased the demand for better understanding of the soil-water regime. Buckingham (1907), working on an universal approach to the simultaneous storage and flow of water and air in soils, postulated the relationship between the capillary potential ψ (Pa) and the volumetric water content θ (m 3 m -3 ), also known as the water retention function, retention curve, or water release curve. The capillary potential follows from the Young (1805)-Laplace relationship, stating that the pressure difference between a liquid and the adjacent gas phase increases inversely proportional to the radius of the interface. In addition to the specific weight of the soil water, Buckingham (1907) introduced the spatial gradient of ψ as the other major driving force, thus allowing for the redistribution of soil water in all directions, evaporation across the soil surface, transpiration via roots, and capillary rise from perched water including groundwater tables. In analogy to Fourier's (1822) and Ohm's (1825) laws for heat flow and electrical current, and Darcy's (1856) law for water flow in saturated porous media, Buckingham (1907) also proposed the hydraulic conductivity for flow in unsaturated porous media as function of either K(θ) or K(ψ) (m s -1 ). According to , the British meteorologist Richardson (1922) was most likely the first who introduced a diffusion type of K-ψ-θ-relationship in the quest of quantifying water exchange between the atmosphere and the soil as lower boundary of the meteorological system. A second-order partial differential expression became necessary because ψ depends on θ, and both their temporal variations on flow, while flow is driven by the gradient of ψ. The race was on to the experimental determination of the K-ψ-θ-relationships.
For instance, Gardner et al. (1922) applied plates and blocks of fired clay with water-saturated pores fine enough to hydraulically connect the capillary bound water within soil samples with systems outside them. Richards(1931) applied the technique to the construction of tensiometers that directly measure ψ within an approximate range of 0 > ψ >≈ -80 kPa (ψ = 0 corresponds to the atmospheric pressure as reference). With the pressure plate apparatus he measured ψ-θrelationships and determined hydraulic conductivity K(ψ or θ). Similar to Richardson (1922), he presented a diffusion-type approach to the transient water flow in unsaturated soils.
Numerous analytical procedures evolved for solving the well-known Richards (1931) equation. Van Genuchten (1980, for instance, developed a closed form of K-ψ-θ-relationships that provide the base for the many hues of HYDRUS, a numerical simulation packages dealing with flow and storage of water and solutes in unsaturated soils (e.g. Simunek et al., 2008). Veihmeyer (1927) investigated water storage in soils in the quest of scheduling irrigation schemes. He proposed the water contents at the field capacity FC and at the permanent wilting point PWP as upper and lower thresholds of plant-available soil water, where FC gets established a couple of days after a soil was saturated under exclusion of evaporation (also referred to as 'drainable or gravitational soil water'). Various methods appeared on how to establish PWP that is accepted today at -15 bars. It became unavoidable that concepts based on Buckingham's (1907) fundamental and seminal work contradicted with practical and fieldoriented research. Veihmeyer (1954), for instance, stated "Since the distinction between capillary and other 'kinds' of water in soils cannot be made with exactness, obviously a term such as non-capillary porosity cannot be defined precisely since by definition it is determined by the amount of 'capillary' water in the soils".
Progress in field instrumentation and computing techniques allowed for producing and processing large data sets including the numerical solution of Richards' (1931) equation. In the late 1970s, the development increasingly unveiled substantial discrepancies between measurements and the numerous approaches to water movement in unsaturated soils based on Richards' (1931) capillarity-dominated theory. Particularly disturbing were observations on wetting fronts advancing much faster than expected from the Richards approach. Concepts like macropore flow (e.g., Beven and Germann, 1982) and flow at non-equilibrium with respect to the ψ-θ-relationship appeared. Jarvis et al. (2016) summarized as preferential all the flows in unsaturated porous media not obeying Richards' (1931) equation. See also Morbidelli et al. (2018) for a recent review on infiltration approaches. Beven (2018) argued that, for about a century, the hardly questioned preference given to capillarity denied recognition of concepts considering flow along macropores, pipes, and cracks. Indeed, there is an increasing number of contributions focusing on the dimensions and shapes of flow paths, their 3-d imaging, and trials to derive flows from them (e.g., Abu Najm et al., 2019). However, there is hardly an approach capable of applying the wealth of information about the paths to the quantification of flow.
Ignoring Veihmeyer's (1954) warning, the attraction of research on flow paths is so dominant that, for instance, Jarvis et al. (2016) flatly denied the applicability of Hagen-Poiseuille concepts to flow in soils. (See Germann, 2017, andJarvis et al. 2017). Moreover, advanced techniques of infiltration with non-Newtonian fluids led so far just to the description of path structures rather than more directly to the flow process (Atalah and Abou Najm, 2018). Wide-spread research in the types, dimensions, and shapes of 'macropores' and their apparent relationships to flow and transport mostly pivot around Richards (1931) equation that is numerically applied to either macropore-/ micropore-domains or by modelling flow and transport in the macropore domain with separate rules yet still maintaining a Richards-type approach in the micropores.
Both approaches allow for due exchange of flow and transport between the two domains.
Imaging procedures visualize flow in 2-d and 3-d in voids as narrow as some 10 µm, rising hope that the wealth of information gained at the hydro-dynamic scale will eventually lead to macroscopic models at the soil profile scale of meters (see, for instance, Jarvis et al., 2016).
Thus, Beven's (2018) denial of progress in infiltration research is here carried a step further.
The obsession with pores, channels, flow paths, and their connectivity, tortuosity, and necks actually retards research progress towards more general infiltration that should be based on hydro-mechanical principles as the 7 th Unresolved Problem in Hydrology demands.
A second thread, leading to the alternative infiltration approach presented here, is traced back to Schumacher's (1864) dual-processes. He suggested that infiltration follows rules, though unspecified at that time, that markedly differ from the capillary rise out of water tables.
Moreover, the alternative approach should be based on the same principles as Hagen-Poiseuille's (1846) and Darcy's (1856) laws, thus closing the gap of one to two orders of magnitude of hydraulic conductivity between saturated flow and flow close to saturation (Germann and Beven, 1981a).
In his quest of demonstrating the benefits of forests and reforestations on controlling floods and debris flows from steep catchments in the Swiss Alps and Pre-Alps, Burger (1922) measured in situ the time lapses Δt100 for the infiltration of 100 mm of water into soil columns of the same length. In the laboratory, he determined the air capacity AC (m 3 m -3 ) of undisturbed samples taken near the infiltration measurements, where AC is the difference of the specific water volume after standardized drainage on a gravel bed and complete saturation. Germann and Beven (1981b) found an encouraging coefficient of determination of r 2 = 0.77 when correlating via a Hagen-Poiseuille (1846) approach 76 pairs of Δt100and AC-values. Following Lawes et al. (1882), who distinguished between fast and slow drainage, Germann (1986) assessed the arrival times of precipitation fronts in the Coshocton lysimeters. Accordingly, rains of 10 (mm/d) sufficed to initiate or increase drainage flow within 24 hours at the 2.4-m depth if the volumetric water content in the upper 1.0 m of the soil was at or above 0.3 (m 3 m -3 ). The observations result in wetting front velocities greater than 3 x 10 -5 (m s -1 ). Beven and Germann (1981) modelled flow in tubes and planar cracks, and proposed kinematic wave theory according to Lighthill and Whitham (1955) as analytical approach to Newton's shear flow. Germann (1985) applied the theory successfully to data from an infiltration-drainage experiment carried out on a block of polyester consolidated coarse sand. The paper is considered a precursor of the following section that treats infiltration and drainage in permeable media as exclusively gravity driven and viscosity controlled, while capillarity may adsorb water from flow to the sessile parts of the system.

Theory a) Basic relationships
The approach is laid out at the hydro-mechanical scale of spatio-temporal process integration, allowing for its easy handling with analytical expressions, yet under strict observance of the balances of energy, momentum, and mass (i.e., the continuity requirements). The approach builds on four presumptions that are not necessarily common to soil hydrology: (i) infiltrating water forms a sharp wetting shock front; (ii) the wetting shock front moves with constant velocity; (iii) atmospheric pressure prevails in the mobile water between the wetting shock front and the surface; and (iv) flow is laminar (i.e.,

Reynolds numbers may not exceed values close to unity).
The interior of a permeable solid medium contains connected flow paths that are wide enough to let liquids pass across the volume considered. The definition purposefully avoids further specification of the flow paths' shapes and dimensions. Water supply to the surface is thought of a pulse P(qS, TB, TE), where qS (m s -1 ) is constant volume flux density from the pulse's beginning at TB to its ending at TE (both s). (The subscript S refers to the surface of the permeable medium). The pulse initiates a water content wave WCW of mobile water that is conceptualized as a film gliding down the paths of a permeable medium according to the rules of Newton's shear flow. The parameters film thickness F (m) and specific contact length L (m m -2 ) per unit cross-sectional area A (m 2 ) of the medium specify a WCW. Regardless of the thickness of F, atmospheric pressure prevails in the film. Figure 1 illustrates the concept. A WCW supposedly runs along the flow paths while forming a discontinuous and sharp wetting shock front at zW(t). The WCW partially fills the upper part of the medium within 0 ≤ z ≤ zW(t) with the mobile water content w(z,t) (m 3 m -3 ), where w < εθante with porosity ε and antecedent soil moisture θante, both (m 3 m -3 ). The lower part z > zW(t) remains at θante.
The coordinate z (m) originates at the surface and points positively down.
Newton (1729) where η (≈ 10 -6 m 2 s -1 ) is the temperature dependent kinematic viscosity of water, ρ (1000 kg m -3 ) is the water's density, v(f) (m s -1 ) is the velocity of the lamina at f in the verticaldown direction, and dv(f)/df is the velocity gradient in the horizontal direction. The expression (N). Integration of Eq.
[2] from the SWI, where v(0) = 0 (the non-slip condition), to f yields the parabolic velocity profile from the SWI to f as (3) The differential volume flux density at f is (m s -1 ). Its integration from the SWI at f = 0 to the air-water interface AWI at f = F produces the volume flux density of the film as (m 3 s -1 ), while the volume of mobile water per unit volume of the permeable medium from the surface to zW(t) amounts to (m 3 m -3 ). The constant velocity of the wetting shock front follows from the volume balance amounting to while the position of the wetting shock front as function of time becomes The terms relating to velocity depend exclusively on F 2 , Eqs. [3,7,8], while those relating to mobile water and its volume flux density also on L, Eqs. [5,6].
Under consideration of zW(t), L expresses the specific vertical contact area of the WCW per unit volume of the permeable medium as the locus where momentum, heat, capillary potential, water, solutes, and particles get exchanged between the WCW and the sessile parts of the medium.
Equations [3 to 8] hold during infiltration i.e., TB ≤ t ≤ TE. Input ends abruptly at TE and at z = 0 i.e., qS → 0, when and where the WCW collapses from f = F to f = 0. All the rear ends of the laminae are released at once at z = 0. They move downwards with v(f) according to Eq.
[3]. The outermost lamina moves the fastest with the celerity of the draining front as Celerity refers to the velocity of a flow property change. The slower moving wetting front intercepts the faster draining front at time TI (s), that follows from setting vW x (TI-TB) = cD x (TI-TE), as ( ) Thus, TI is an exclusive expression of the pulse duration. The wetting front intercepts the draining front at depth The rear ends of all the other laminae move with decreased celerities. According to Eq.
[9], the celerity cRE (f) of the rear end of a lamina between 0 < f < F is Rearranging the last two terms in Eq.
[12] and solving for f leads to the temporal position of the film thickness as ( ) Multiplication of FRE (z,t) with L leads to the spatio-temporal mobile water content of the WCW as trailing wave according to ( ) After TI and beyond ZI the draining front disappears and vW(z,t) decreases with time and depth. However, the shape of the WCW remains according to Eq. [14] over the depth range of 0 ≤ z ≤ zW(t). The volume balance of the WCW amounts to where VWCW (m) represents the water volume of the WCW that has infiltrated during TE-TB.
[15] for zW(t) leads to the temporal position of the wetting shock front as The first derivative of Eq.
[16] produces the wetting shock front velocity as Re ≤ 1 strictly defines laminar flow; however, depending on the application, Re > 1 might be tolerable, yet within an undisclosed range.
The following paragraph b) provides cuts of the WCW in the z-w(z,t)-plane of Fig. 2 while paragraph c) introduces cuts in the t-w(z,t)-plane.

b) Variation of a WCW with depth
The following inspects the spatial variation w(z,τ) of the WCW's mobile water content during the three intervals of (i) TB ≤ τ1 ≤ TE, (ii) TE ≤ τ2 ≤ TI , and (iii) TI ≤ τ3 < ∞ .
(i) TB ≤ τ1 ≤ TE: Position of the wetting shock front, mobile water content, and volume flux density are Thus, piston flow occurs during infiltration, TB ≤ τ1 ≤ TE.
(ii) TE ≤ τ2 ≤ TI: The position of the wetting shock front is the same as in Eq.
The mobile water content remains constant at w between zD(τ2)≤ z ≤ zW (τ2)

c) Variation of a WCW with time
The following inspects the time series of a WCW's mobile water content w(ζ,t) at the three depth ranges of (i) 0 ≤ ζ1 < ZI , (ii) ζ2 = ZI , and (iii) ζ3 ≥ ZI .
(i) 0 ≤ ζ1 < ZI: The arrival times of the wetting shock and draining fronts at ζ1 are while the mobile water content assumes the following values during the respective time intervals: (ii) ζ2 = ZI : At depth of front interception and after t ≥ TI the mobile water content becomes and the mobile water content as a function of time becomes

d) Routing pulse series
From mass balance requirement follows the celerity of an abrupt pulse increase from P1 to P2 with q2 > q1 and w2 > w1, as

Experimental support for the four presumptions
This section experimentally supports the four presumptions that provide the base for Eqs. [1 to 36]: (i) infiltrating water forms a sharp wetting shock front; (ii) the wetting shock front moves with constant v; (iii) atmospheric pressure prevails in the WCW; and (iv) flow is laminar (i.e., Tank (alHagrey et al., 1999) produced the data that support presumptions (i) to (iv).     infiltration (Flammer et al., 2001), Germann (2018a) concluded that atmospheric pressure prevails in the WCW between the wetting shock front and the soil surface at least during TB ≤ t ≤ TI.

e) Laminar flow, Reynolds number: presumption (iv)
From the application of Eq.

a) Coherence of the approach
The parameters F and L suffice to treat infiltration and drainage with Newton's shear flow approach, Eqs.
[1] to [20]. In principle, time series of either θ(Z,t) or q(Z,t) prmits calibration of the two parameters. Both procedures are introduced, using the data presented in Fig. 6. [20], and integrating the resulting expression from tW(Z) to t > tD(Z), yielding where Z refers to the depth of drainage flow at 2.0 (m). The specific contact area is the only factor left for matching Eq.
[37] to the data that resulted in L(q) = 3.3 x 10 4 (m 2 m -3 ), comfortably lying within the range of L(w). This demonstrates the coherence of Newton's shear flow approach to infiltration and drainage. Dubois (1991)

b) Scale tolerance
depths of 100, 150, and 200 (mm). Fluxes in each layer followed from Newton's shear flow approach. The flux differences from layer to layer deviated utmost by 19% from the corresponding water content changes in the volumes between the layers (Germann, 2014).
Dubois ' (1991) observations across 1800 (m) of crystalline rocks of the Mont Blanc massif and the water balance calculations of finger flow in the sand box of Hincapié and Germann (2010) at the scale of millimeters hint at the spatio-temporal tolerance of Newton's shear flow that may advance the approach to an attractive tool, for instance, for the study of infiltration into groundwater systems.

c) Preferential and retarded tracer breakthrough
Preferential flow in soil hydrology is frequently associated with enhanced and accelerated solute and pollutant breakthrough (e.g., Larsbo et al., 2014). However, Bogner and Germann (2019) reported considerable delays of tracer breakthrough compared with the first arrival of the wetting shock fronts at the bottoms of soil columns with heights of 0.4 (m). They referred to the phenomenon as 'pushing out old water' that is well known in catchment hydrology. They statistically explained 81% of the observed delay variations with combinations of L and F when applying Newton's shear flow to the data. Tracer exchange on large L from thin F of the WCW may be even faster than presumed 'preferential' tracer transport. Under consideration of the mechanistic parameters F and L, Newton's shear flow provides for a novel tool for the unambiguous investigation of tracer transport and exchange i.e., accelerated as well as decelerated breakthrough.

d) Gravity vs. capillarity
Schumacher (1864) suggested a two-process approach to water flow and storage in partially saturated permeable media. While he recognized capillarity as responsible for the water's rise, and probably also its contribution to water redistribution in soil columns, he left open the mechanism behind infiltration. This paper concentrates on infiltration that is completely gravity-driven and viscosity-controlled, yet allowing for water abstraction due to capillarity from the mobile to the immobile part of the permeable system. Concentrating on gravity and viscosity liberates infiltration and drainage from the omnipresence of capillarity in soil hydrology with the benefit of avoiding the difficult definitions of non-equilibrium flow and the separation of macropores from the remaining pores. With respect to capillarity, the relative contribution of gravity to flow varies according to cos(α), where α (°) is the angle of deviation from the vertical. Thus, at cos(0°) = 1, as in the cases presented above, gravity's contribution is at maximum; it reduces to cos(90°) = cos(270°) = 0, while it completely opposes capillarity at cos(180°) = -1. Darcy's (1856) law mutates to an extension of unsaturated vertical shear flow. From Eq.

e) Shear flow and Darcy's law
[5] follows: θante + w < ε and Δp/Δz = ρ g: where θante (m 3 m -3 ) is the antecedent volumetric water content, ε (m 3 m -3 ) is porosity, Δp/Δz (Pa m -1 ) is the pressure gradient, ρ (=1000 kg m -3 ) is the density of water, and µ = ρ η (Pa s) is dynamic viscosity. At saturation we get: θante + w = ε and Δp/Δz = ρ g: Darcy's law states that q  p/z i.e., volume flux density is a linear function of the flowdriving gradient with the proportionality factor Ksat. In view of the various dimensionalities of w  (L 1 ,F 1 ), v  ( L 0 , F 2 ), and q  (L 1 , F 3 ), linearity seems only possible if Fsat and Lsat remain constant and independent from p in the transition from gravity-driven to pressure-driven shear flow at saturation i.e., in the transition from Eq. [38] to Eq. [40]. This elaboration supports the linearity of Darcy's law, but it is not its independent proof. As a consequence, w = q/v also remains constant. Further, if θante + w = ε, dLsat/dp = 0 and dFsat/dp = 0 then follows the hypothesis that (Fsat x Lsat) represent (F x L)max leading to Ksat. However, other combinations of (F x L) in unsaturated media are feasible that may lead to q > qsat = Ksat. This unproven speculation opens an unexpected view on shear flow, that is in stark contrast to Richards (1931) capillary flow, where a priori Ksat > K(θ or ψ). See Germann and Karlen (2016) for further discussion.

f) Water abstraction from the WCW
Pressure in the WCW is atmospheric while ψ < 0 typically prevails ahead of it. Therefore, water is abstracted from the WCW onto L. Abstraction is usually completed during short periods as the θ(Z,t)-series in Fig. 6 demonstrate. The amount of abstraction shows in the difference between θend and θante.

Conclusions
Newton's shear flow provides for a cohesive approach to infiltration and drainage in permeable media, and no a priori decisions on pore properties are required. So far, the approach is in its descriptive mode, capable of quantifying infiltration and drainage with the two parameters film thickness F and specific contact area L. However, the analytical expressions facilitate the development of predictive model applications such as to groundwater recharge and to the transport of solutes and particles. Advances are expected from research, among other topics, in the relationships of F and L with antecedent soil moisture, intensity of infiltration, and hydraulic conductivity Ksat.
Finally, Newton's shear flow seems to have solved the 7 th Unsolved Problem in Hydrology (Blöschl et al., 2019) that asks "Why is most flow preferential across multiple scales and how does such behaviour co-evolve with the critical zone?". However, Newton's shear flow as the solution of the 7 th UPH did not evolve from the suggested dual-porosity perspective but from a hydro-mechanical point of view that does neither require preferential flow nor coevolution of flow paths.