2 Formulation of the Model

[page 509, §1]   
[509.1.1] Consider a onedimensional, homogeneous, isotropic and incompressible porous column filled with two immiscible Newtonian fluids. [509.1.2] In a one dimensional model transversal variations and column thickness are neglected. [509.1.3] Let SW=SW⁢x,t denote the saturation of wetting fluid (called water), and SO=SO⁢x,t the saturation of nonwetting fluid (called oil). [509.1.4] Here time is t≥0, and x∈0,L is the position in a column of length L. [509.1.5] Each of the two fluid phases is considered to consist of a continuous, mobile, percolating subphase, and of a discontinuous, isolated, trapped or nonpercolating subphase as discussed in [23, 28, 25, 24, 26]. [509.1.6] Following the notation of [26] the percolating phase of water is [page 510, §0]    indexed by i=1 and its nonpercolating phase is indexed by i=2. [510.0.1] The water saturation is then SW=S1+S2. [510.0.2] The percolating oil phase is indexed as i=3 and its nonpercolating phase by i=4. [510.0.3] The oil saturation is SO=S3+S4. [510.0.4] The volume fraction ϕi of phase i is defined as ϕi=ϕ⁢Si where the volume fraction of the pore space, also called porosity, is ϕ. [510.0.5] The volume fraction of the solid matrix is denoted as ϕ5=1-ϕ. [510.0.6] Volume conservation requires

to hold.

[510.1.1] The mass balance of fluid phase i can be expressed in differential form as

where ϱi⁢x,t,ϕi⁢x,t,vi⁢x,t are mass density, volume fraction and velocity of phase i as functions of position x∈S=0,L⊂R and time t∈R+. [510.1.2] Mi is the mass transfer rate from all other phases into phase i.

[510.2.1] The momentum balance is written as (i=1,2,3,4)

where Σi is the stress tensor in the ith phase, Fi is the body force per unit volume acting on the ith phase, mi is the momentum transfer into phase i from all the other phases, and Di/D⁢t=∂/∂⁡t+vi⁢∂/∂⁡x denotes the material derivative for phase i. [510.2.2] The stress tensor for the nonpercolating phases is defined as the momentum flux across surfaces in the threedimensional continuum (see [29] for more discussion).

[510.3.1] For a macroscopically homogeneous porous medium

is assumed. [510.3.2] Incompressible fluids are assumed so that their densities

are independent of x,t.

[510.4.1] The percolating and the nonpercolating phases are able to exchange mass through breakup and coalescence of droplets, ganglia and clusters. [510.4.2] The mass transfer rates must depend on rates of saturation change. [510.4.3] They are here assumed to be

[page 511, §0]    where η2,η4 are constants. [511.0.1] The parameters SW*,SO*, S2*,S4* are defined by

where SW⁢i,SO⁢r are limiting saturations for S2,S4. [511.0.2] In eq. (7) the shorthand ∂t=∂/∂⁡t is used and

denotes the Heaviside unit step function. [511.0.3] Equation (7) follows from the form used in [26] for small rates of saturation change. [511.0.4] The mass exchange depends on the sign of ∂t⁡SW. [511.0.5] The sign determines the type of process. [511.0.6] It can switch locally between drainage and imbibition. [511.0.7] This results in hysteresis. [511.0.8] The structure of the mass exchange term was chosen with hindsight such that theoretical results obtained in the residual decoupling approximation agree with experimental measurements of capillary pressure. [511.0.9] The mass exchange terms have recently been further generalized to reproduce not only capillary pressure, but also experimental capillary desaturation curves [27].

[511.1.1] Turning to the momentum balance, note first that the inertial term will not be neglected in this paper. [511.1.2] The stress tensor for the four phases are specified as

where P1 and P3 are the fluid pressures in the percolating phases. [511.1.3] The constants P2*,P4* and exponents γ,δ are constitutive parameters. [511.1.4] The assumptions for the nonpercolating phases reflect their modified pressure. [511.1.5] This phenomenon seems to have been observed in experiment [2, Fig. 2, p. 233]. ^a (This is a footnote:) ^aFig. 2 in [2] shows that the pressure measured by the pore pressure transducers PPT3 and PPT4 rebounds after the end of the infiltration, i.e. when the DNAPL has passed and the transducers measure the pressure of water and residual PCE [511.1.6] In applications, the parameters P2*,P4* and γ,δ are determined by measuring capillary pressure curves (see below). [511.1.7] The body forces are assumed to be given by gravity and capillarity. [511.1.8] They are specified as

with constitutive constants Πa,Πb and exponents α,β>0. [511.1.9] The angle 0≤ϑ≤π/2 is the angle between the direction of the column and the direction of gravity with ϑ=π/2 [page 512, §0]    corresponding to alignment. [512.0.1] In applications the parameters Πa,Πb and α,β are determined by measuring capillary pressure curves. [512.0.2] The capillary body forces in eqs. (10) reflect the wetting properties of the medium. [512.0.3] They oppose gravity and reduce buoyancy driven flows of the disconnected phases. [512.0.4] This is illustrated in the figures below.

[512.1.1] Finally, the momentum transfer terms are assumed to be given by linear viscous drag characterized by constitutive resistance coefficients Ri⁢j through the equations

where R12=0 and R34=0 was used because there is no common interface and hence no direct viscous interaction between these phase pairs. [512.1.2] Remember that the index i=5 represents the rock matrix. [512.1.3] For more details on these constitutive assumptions the reader is referred to the original papers [24, 25, 26].

[512.2.1] The balance laws (1b), (2) and (3) combined with the constitutive assumptions given above provide 9 equations for the 10 unknowns Si,vi,P1,P3 with i=1,2,3,4. [512.2.2] To close the system of equations the conditions v2=0 or v4=0 could be used. [512.2.3] These conditions apply when the nonpercolating phases are immobile as it is often observed in experiment. [512.2.4] It turns out, however, that there exists a less restrictive and, in our opinion, more natural selfconsistent closure.

[512.3.1] The selfconsistent closure condition used in this paper follows naturally from many limiting cases. [512.3.2] One such limit is the residual decoupling approximation close to hydrostatic equilibrium described in detail in [26, Section 5., p. 216ff]. [512.3.3] A second, more general limiting case is the limit of vanishing velocities, i.e. vi→0 for i=1,2,3,4. [512.3.4] Here we formulate the selfconsistent closure condition in its most general form as

[page 513, §1]    [513.2.1] This condition follows selfconsistently from the constitutive theory. [513.2.2] It expresses the experimental observation that the pressure difference P3-P1 depends more strongly on saturations than on velocities, and that it remains nonzero even for vanishing velocities. [513.2.3] Adding eqs. (3) for i=2 and i=3 and subtracting eqs. (3) with i=1 and i=4 from the result gives

[513.2.4] In this form the selfconsistent closure has been used in the numerical calculation below.

[513.3.1] The mathematical model defined above was introduced and discussed in [25, 24, 26]. [513.3.2] It was recently extended to include surface tension [27]. [513.3.3] Notation and model formulation in this paper follow [26]. [513.3.4] Approximations and analytical solutions for some special cases were given in [25, 24, 26]. [513.3.5] Here the system of equations will be solved by numerical methods. [513.3.6] To this end initial and boundary conditions are needed and will be discussed next.