A phase-by-phase upstream scheme that converges to the vanishing capillarity solution for countercurrent two-phase flow in two-rock media

We discuss the convergence of the upstream phase-by-phase scheme (or upstream mobility scheme) towards the vanishing capillarity solution for immiscible incompressible two-phase flows in porous media made of several rock types. Troubles in the convergence were recently pointed out by Mishra and Jaffré (Comput. Geosci. 14, 105–124, 2010) and Tveit and Aavatsmark (Comput. Geosci. 16, 809–825, 2012). In this paper, we clarify the notion of vanishing capillarity solution, stressing the fact that the physically relevant notion of solution differs from the one inferred from the results of Kaasschieter (Comput. Geosci. 3, 23–48, 1999). In particular, we point out that the vanishing capillarity solution depends on the formally neglected capillary pressure curves, as it was recently proven in by Andreianov and Cancès (Comput. Geosci. 17, 551–572, 2013). Then, we propose a numerical procedure based on the hybridization of the interfaces that converges towards the vanishing capillarity solution. Numerical illustrations are provided.     

Model reduction and discretization using hybrid finite volumes for flow in porous media containing faults

In this paper, we study two different model reduction strategies for solving problems involving single phase flow in a porous medium containing faults or fractures whose location and properties are known. These faults are represented as interfaces of dimension N ? 1 immersed in an N dimensional domain. Both approaches can handle various configurations of position and permeability of the faults, and one can handle different fracture permeabilities on the two inner sides of the fracture. For the numerical discretization, we use the hybrid finite volume scheme as it is known to be well suited to simulating subsurface flow. Some results, which may be of use in the implementation of the proposed methods in industrial codes, are demonstrated.     

Vanishing capillarity solutions of Buckley–Leverett equation with gravity in two-rocks’ medium

For the hyperbolic conservation laws with discontinuous-flux function, there may exist several consistent notions of entropy solutions; the difference between them lies in the choice of the coupling across the flux discontinuity interface. In the context of Buckley–Leverett equations, each notion of solution is uniquely determined by the choice of a “connection,” which is the unique stationary solution that takes the form of an under-compressive shock at the interface. To select the appropriate connection, following Kaasschieter (Comput Geosci 3(1):23–48, 1999), we use the parabolic model with small parameter that accounts for capillary effects. While it has been recognized in Cancès (Networks Het Media 5(3):635–647, 2010) that the “optimal” connection and the “barrier” connection may appear at the vanishing capillarity limit, we show that the intermediate connections can be relevant and the right notion of solution depends on the physical configuration. In particular, we stress the fact that the “optimal” entropy condition is not always the appropriate one (contrarily to the erroneous interpretation of Kaasschieter’s results which is sometimes encountered in the literature). We give a simple procedure that permits to determine the appropriate connection in terms of the flux profiles and capillary pressure profiles present in the model. This information is used to construct a finite volume numerical method for the Buckley–Leverett equation with interface coupling that retains information from the vanishing capillarity model. We support the theoretical result with numerical examples that illustrate the high efficiency of the algorithm.     

Sequential implicit vertex approximate gradient discretization of incompressible two-phase Darcy flows with discontinuous capillary pressure

This work deals with sequential implicit schemes for incompressible and immiscible two-phase Darcy flows which are commonly used and well understood in the case of spatially homogeneous capillary pressure functions. To our knowledge, the stability of this type of splitting schemes solving sequentially a pressure equation followed by the saturation equation has not been investigated so far in the case of discontinuous capillary pressure curves at different rock type interfaces. It will be shown here to raise severe stability issues for which stabilization strategies are investigated in this work. To fix ideas, the spatial discretization is based on the Vertex Approximate Gradient (VAG) scheme accounting for unstructured polyhedral meshes combined with an Hybrid Upwinding (HU) of the transport term and an upwind positive approximation of the capillary and gravity fluxes. The sequential implicit schemes are built from the total velocity formulation of the two-phase flow model and only differ in the way the conservative VAG total velocity fluxes are approximated. The stability, accuracy and computational cost of the sequential implicit schemes studied in this work are tested on oil migration test cases in 1D, 2D and 3D basins with a large range of capillary pressure parameters for the drain and barrier rock types. It will be shown that usual splitting strategies fail to capture the right solutions for highly contrasted rock types and that it can be fixed by maintaining locally the pressure saturation coupling at different rock type interfaces in the definition of the conservative total velocity fluxes. The numerical investigation of the sequential schemes is also extended to the widely used finite volume Two-Point Flux Approximation spatial discretization.

     

Effect of 2-D Random Field Discretization on Failure Probability and Failure Mechanism in Probabilistic Slope Stability

This paper investigates, using the random field theory and Monte Carlo simulation, the effects of random field discretization on failure probability, p _f, and failure mechanism of cohesive soil slope stability. The spatial sizes of the discretized elements in random field Δx, Δy in horizontal and vertical directions, respectively, are assigned a series of combinational values in order to model the discretization accuracy. The p _f of deterministic critical slip surface (DCSS) and that of the slope system both are analyzed. The numerical simulation results have demonstrated that both the ratios of Δy/λ _y (λ _y  = scale of fluctuation in vertical direction) and Δx/λ _x (λ _x  = scale of fluctuation in horizontal direction) contribute in a similar manner to the accuracy of p _f of DCSS. The effect of random field discretization on the p _f can be negligible if both the ratios of Δx/λ _x and Δy/λ _y are no greater than 0.1. The normalized discrepancy tends to increase at a linear rate with Δy/λ _y when Δx/λ _x is larger than 0.1, and vice versa for p _f of DCSS. The random field discretization tends to have more considerable influence on the p _f of DCSS than on that of the slope system. The variation of p _f versus λ _x and λ _y may exhibit opposite trends for the cases where the limit state functions of slope failure are defined on DCSS and on the slope system as well. Finally, the p _f of slope system converges in a more rapid manner to that of DCSS than the failure mechanism does to DCSS as the spatial variability of soil property grows from significant to negligible.     

Higher resolution total velocity <Emphasis Type="Italic">Vt</Emphasis> and <Emphasis Type="Italic">Va</Emphasis> finite-volume formulations on cell-centred structured and unstructured grids

Novel cell-centred finite-volume formulations are presented for incompressible and immiscible two-phase flow with both gravity and capillary pressure effects on structured and unstructured grids. The Darcy-flux is approximated by a control-volume distributed multipoint flux approximation (CVD-MPFA) coupled with a higher resolution approximation for convective transport. The CVD-MPFA method is used for Darcy-flux approximation involving pressure, gravity, and capillary pressure flux operators. Two IMPES formulations for coupling the pressure equation with fluid transport are presented. The first is based on the classical total velocity Vt fractional flow (Buckley Leverett) formulation, and the second is based on a more recent Va formulation. The CVD-MPFA method is employed for both Vt and Va formulations. The advantages of both coupled formulations are contrasted. The methods are tested on a range of structured and unstructured quadrilateral and triangular grids. The tests show that the resulting methods are found to be comparable for a number of classical cases, including channel flow problems. However, when gravity is present, flow regimes are identified where the Va formulation becomes locally unstable, in contrast to the total velocity formulation. The test cases also show the advantages of the higher resolution method compared to standard first-order single-point upstream weighting.     

A study on iterative methods for solving Richards’ equation

This work concerns linearization methods for efficiently solving the Richards equation, a degenerate elliptic-parabolic equation which models flow in saturated/unsaturated porous media. The discretization of Richards’ equation is based on backward Euler in time and Galerkin finite elements in space. The most valuable linearization schemes for Richards’ equation, i.e. the Newton method, the Picard method, the Picard/Newton method and the L-scheme are presented and their performance is comparatively studied. The convergence, the computational time and the condition numbers for the underlying linear systems are recorded. The convergence of the L-scheme is theoretically proved and the convergence of the other methods is discussed. A new scheme is proposed, the L-scheme/Newton method which is more robust and quadratically convergent. The linearization methods are tested on illustrative numerical examples.     

Gravitational darkening coefficients of the stars in the semidetached binary V Puppis

We have analyzed high-precision vby light curvesfor the semi-detached binary V Pup in a Roche model. They are consistent with the standard gravitational darkening coefficient for hot stars, β = 0.25, rather than the value β = 1.36 ± 0.04 derived by Kitamura and Nakamura [1] using a simpler model. We rigorously estimate the confidence intervals for the allowed gravitational darkening coefficients for a star filling its Roche lobe to be β = (?0.24, +1.29) for the 99% confidence level and β = (?0.21, +1.26) for the 67% confidence level.     

Method of negative saturations for flow with variable number of phases in porous media: extension to three-phase multi-component case

Various EOR methods lead to the appearance of specific macroscopic surfaces called interfaces of phase transition (IPT) such that the number of phases on two sides of an IPT is different, and fluids separated by an IPT are in non-equilibrium. Therefore, the flow equations are also different on two sides of an IPT and cannot be deduced from each other by a continuous degeneration, which imposes difficulties in numerical modelling. To describe such systems, we developed a new conceptual mathematical method based on the replacement of the real single-phase fluid by an imaginary multi-phase multi-component continuum having fictitious properties. As the result, the fluid over all zones becomes multi-phase and can be described by uniform multi-phase hydro- and thermodynamic equations, which allows applying the direct numerical simulation. The equivalence principle determines the physical properties of the fictitious multi-phase fluid, as well as the structure of the uniform multi-phase equations. It also proves that the saturation of each phase becomes an extended function negative or higher than unity in non-equilibrium zones, which becomes the efficient method of tracking the interfaces, the number of phases at any point, and their degree of disequilibrium. The method was developed in [1, 2] for the two-phase case. In the present paper, the new version of the method is developed for the three-phase case with gravity, diffusion, and capillarity. We have obtained the new equivalent uniform multi-phase equations which contain additional non-classical terms responsible for the diffusion and gravity across an IPT. The comparison with classical method is presented. The presentation is illustrated by several examples of simulation by means of the code developed by the research group; their concern: EOR by miscible methods and CO ₂ bubble raising in aquifer.     

A numerical method for the solution of the nonlinear turbulent one-dimensional free surface flow equations

Free surface flow of an incompressible fluid over a shallow plane/undulating horizontal bed is characteristically turbulent due to disturbances generated by the bed resistance and other causes. The governing equations of such flows in one dimension, for finite amplitude of surface elevation over the bed, are the Continuity Equation and a highly nonlinear Momentum Equation of order three. The method developed in this paper introduces the “discharge” variable q = η U, where η = elevation of the free surface above the bed level, and U = average stream-wise forward velocity. By this substitution, the continuity equation becomes a linear first-order PDE and the momentum equation is transformed after introduction of a small approximation in the fifth term. Next, it is shown by an invertibility argument that q can be a function of η: q = F(η), rendering the momentum equation as a first order, second degree ODE for F(η), that can be be integrated by the Runge-Kutta method. The continuity equation then takes the form of a first order evolutionary PDE that can be integrated by a Lax-Wendroff type of scheme for the temporal evolution of the surface elevation η. The method is implemented for two particular cases: when the initial elevation is triangular with vertical angle of 120 ^° and when it has a sinusoidal form. The computations exhibit the physically interesting feature that the frontal portion of the propagating wave undergoes a sharp jump followed by tumbling over as a breaker. Compared to other discretization methods, the application of the Runge-Kutta and an extended version of the Lax-Wendroff scheme is much easier.     

Analytical expressions for three-phase generalized relative permeabilities in water- and oil-wet capillary tubes

We analyze three-phase flow of immiscible fluids taking place within an elementary capillary tube with circular cross-section under water- and oil-wet conditions. We account explicitly for momentum transfer between the moving phases, which leads to the phenomenon of viscous coupling, by imposing continuity of velocity and shear stress at fluid-fluid interfaces. The macroscopic flow model which describes the system at the Darcy scale includes three-phase effective relative permeabilities, K _{i j,r} , accounting for the flux of the ith phase due to the presence of the jth phase. These effective parameters strongly depend on phase saturations, fluid viscosities, and wettability of the solid matrix. In the considered flow setting, K _{i j,r} reduce to a set of nine scalar quantities, K _{i j,r} . Our results show that K _{i j,r} of the wetting phase is a function only of the fluid phase own saturation. Otherwise, K _{i j,r} of the non-wetting phase depends on the saturation of all fluids in the system and on oil and water viscosities. Viscous coupling effects (encapsulated in K _{i j,r} with i ≠ j) can be significantly relevant in both water- and oil-wet systems. Wettability conditions influence oil flow at a rate that increases linearly with viscosity ratio between oil and water phases.     

Finite volume approximation for an immiscible two-phase flow in porous media with discontinuous capillary pressure

We consider an immiscible incompressible two-phase flow in a porous medium composed of two different rocks so that the capillary pressure field is discontinuous at the interface between the rocks. This leads us to apply a concept of multivalued phase pressures and a notion of weak solution for the flow which have been introduced in Cancès and Pierre (SIAM J Math Anal 44(2):966–992, 2012). We discretize the problem by means of a numerical algorithm which reduces to a standard finite volume scheme in each rock and prove the convergence of the approximate solution to a weak solution of the two-phase flow problem. The numerical experiments show in particular that this scheme permits to reproduce the oil-trapping phenomenon.     

Minimum velocity dispersion in stable stellar disks. Numerical simulations

N-body dynamical simulations are used to analyze the conditions for the gravitational stability of a three-dimensional stellar disk in the gravitational field of two rigid spherical components—a bulge and halo whose central concentrations and relative masses vary over wide ranges. The number of point masses N in the simulations varies from 40 to 500 000 and the evolution of the simulated systems is followed over 10–20 rotation periods of the outer edge of the disk. The initially unstable disks are heated and, as a rule, reach a quasi-stationary equilibrium with a steady-state radial-velocity dispersion c_r over five to eight turns. The radial behavior of the Toomre stability parameter Q_T(r) for the final state of the disk is estimated. Simple models are used to analyze the dependence of the gravitational stability of the disk on the relative masses of the spherical components, disk thickness, degree of differential rotation, and initial state of the disk. Formal application of existing, analytical, local criteria for marginal stability of the disk can lead to errors in c_r of more than a factor of 1.5. It is suggested that the approximate constancy of Q_T?1.2–1.5 for r?(1–2)×L (where L is the radial scale of disk surface density), valid for a wide range of models, can be used to estimate upper limits for the mass and density of a disk based on the observed distributions of the rotational velocity of the gaseous component and of the stellar velocity dispersion.     

Ordering-based nonlinear solver for fully-implicit simulation of three-phase flow

The Fully Implicit method (FIM) is often the method of choice for the temporal discretization of the partial differential equations governing multiphase flow in porous media. The FIM involves solving large coupled systems of nonlinear algebraic equations. Newton-based methods, which are employed to solve the nonlinear systems, can suffer from convergence problems—this is especially true for large time steps in the presence of highly nonlinear flow physics. To overcome such convergence problems, the time step is usually reduced, and the Newton steps are restarted from the solution of the previous (converged) time step. Recently, potential ordering and the reduced-Newton method were used to solve immiscible three-phase flow in the presence of buoyancy and capillary effects (e.g., Kwok and Tchelepi, J. Comput. Phys. 227(1), 706–727 2007). Here, we improve the robustness of the potential-based ordering method in the presence of gravity. Furthermore, we also extend this nonlinear approach to interphase mass transfer. Our algorithm deals effectively with mass transfer between the liquid and gas phases, including phase disappearance (e.g., gas going back in solution) and reappearance (e.g., gas coming out of solution and forming a separate phase), as a function of pressure and composition. Detailed comparisons of the robustness and efficiency of the potential-based solver with state-of-the-art nonlinear/linear solvers are presented for immiscible two-phase (Dead-Oil), Black-Oil, and compositional problems using heterogeneous models. The results show that for large time steps, our nonlinear ordering-based solver reduces the number of nonlinear iterations significantly, which leads to gains in the overall computational cost.     

Ordering-based nonlinear solver for fully-implicit simulation of three-phase flow

The Fully Implicit Method (FIM) is often the method of choice for the temporal discretization of the partial differential equations governing multiphase flow in porous media. The FIM involves solving large coupled systems of nonlinear algebraic equations. Newton-based methods, which are employed to solve the nonlinear systems, can suffer from convergence problems—this is especially true for large time steps in the presence of highly nonlinear flow physics. To overcome such convergence problems, the time step is usually reduced, and the Newton steps are restarted from the solution of the previous (converged) time step. Recently, potential ordering and the reduced-Newton method were used to solve immiscible three-phase flow in the presence of buoyancy and capillary effects (e.g., Kwok and Tchelepi, J. Comput. Phys. 227(1), 706–727 9). Here, we improve the robustness of the potential-based ordering method in the presence of gravity. Furthermore, we also extend this nonlinear approach to interphase mass transfer. Our algorithm deals effectively with mass transfer between the liquid and gas phases, including phase disappearance (e.g., gas going back in solution) and reappearance (e.g., gas coming out of solution and forming a separate phase), as a function of pressure and composition. Detailed comparisons of the robustness and efficiency of the potential-based solver with state-of-the-art nonlinear/linear solvers are presented for immiscible two-phase (Dead-Oil), Black-Oil, and compositional problems using heterogeneous models. The results show that for large time steps, our nonlinear ordering-based solver reduces the number of nonlinear iterations significantly, which leads to gains in the overall computational cost.     

Abundances of carbon,nitrogen, oxygen,and other elements in the atmospheres of the giants 25 Mon and HR 7389

An analysis of high-resolution CCD spectra of the giant 25 Mon, which shows signs of metallicity, and the normal giant HR 7389 is presented. The derived effective temperatures, gravitational accelerations, and microturbulence velocities are T_eff = 6700 K, log g = 3.24, and ξ _t = 3.1 km/s for 25 Mon and T_eff = 6630 K, log g = 3.71, and ξ _t = 2.6 km/s for HR 7389. The abundances (log ε) of nine elements are determined: carbon, nitrogen, oxygen, sodium, silicon, calcium, iron, nickel, and barium. The derived excess carbon abundances are 0.23 dex for 25 Mon and 0.16 dex for HR 7389. 25 Mon displays a modest (0.08 dex) oxygen excess, with the oxygen excess for HR 7389 being somewhat higher (0.15 dex). The nitrogen abundance is probably no lower than the solar value for both stars. The abundances of iron, sodium, calcium (for HR 7389), barium, and nickel exceed the solar values by 0.22–0.40 dex for both stars. The highest excess (0.62 dex) is exhibited by the calcium abundance for 25 Mon. Silicon displays a nearly solar abundance in both stars—small deficits of ?0.03 dex and ?0.07 dex for 25 Mon and HR 7389, respectively. No fundamental differences in the elemental abundances were found in the atmospheres of 25 Mon and HR 7389. Based on their T_eff and log g values, as well as theoretical calculations, A. Claret estimated the masses, radii, luminosities, and ages of 25 Mon (M/M_⊙ = 2.45, log(R/R_⊙) = 0.79, log(L/L_⊙) = 1.85, t = 5.3 × 10⁸ yr) and HR 7389 (M/M_⊙ = 2.36, log(R/R_⊙) = 0.50, log(L/L_⊙) = 1.24, t = 4.6 × 10⁸ yr), and also of the stars 20 Peg (M/M_⊙ = 2.36, log(R/R_⊙) = 0.73, log(L/L_⊙) = 1.79, t = 4.9 × 10⁸ yr) and 30 LMi (M/M_⊙ = 2.47, log(R/R_⊙) = 0.73, log(L/L_⊙) = 1.88, t = 4.8 × 10⁸ yr) studied by the author earlier.     

Fractal dimension of pore space in carbonate samples from Tushka area (Egypt)

To investigate inhomogeneous and porous structures in nature, the concept of fractal dimension was established. This paper briefly introduces the definition and measurement methods of fractal dimension. Three different methods including mercury injection capillary pressure (MICP), nuclear magnetic resonance (NMR), and nitrogen adsorption (BET) were applied to determine the fractal dimensions of the pore space of eight carbonate rock samples taken from West Tushka area, Egypt. In the case of fractal behavior, the capillary pressure P _c and cumulative fraction V _c resulting from MICP are linearly related with a slope of D-3 in a double logarithmic plot with D being the value of fractal dimension. For NMR, the cumulative intensity fraction V _c and relaxation time T ₂ show a linear relation with a slope of 3-D in a double logarithmic plot. Fractal dimension can also be determined by the specific surface area S _por derived from nitrogen adsorption measurements and the effective hydraulic radius. The fractal dimension D shows a linear relation with the logarithm of S _por . The fractal dimension is also used in models of permeability prediction. To consider a more comprehensive data set, another 34 carbonate samples taken from the same study area were integrated in the discussion on BET method and permeability prediction. Most of the 42 rock samples show a good agreement between measured permeability and predicted permeability if the mean surface fractal dimension for each facies is used.     

Effect of pore connectivity on reflection amplitudes of an inhomogeneous wave in a composite porous solid saturated by two immiscible fluids

Present paper aims to study the phenomenon of reflection and transmission when an inhomogeneous wave strikes some discontinuity in a composite porous medium saturated by two immiscible viscous fluids. The incident wave splits into six reflected and six transmitted waves at the interface. All reflected and transmitted waves are inhomogeneous in nature with different directions of propagation vector and attenuation vector. A dimensionless parameter \(\varsigma \in [0, 1]\) is introduced to represent the extent of connection among the pores at the interface. Expression of Umov–Poynting vector is derived to obtain energy flux vector. Continuity of energy flux vector at the interface gives the required boundary conditions for the system. Connecting parameter \(\varsigma \) is also employed in boundary conditions to model the partial connection of pores at the interstices of two media. For numerical discussion we consider a porous medium composed of sandstone and ice, saturated with oil and water. The effect of parameter \(\varsigma \) and angle of incidence is determined numerically on the amplitude and the energy ratios of reflected and transmitted waves.