Boris Andreianov
email: boris.andreianov@univ-fcomte.fr

Clément Konstantin Brenner
email: konstantin.brenner@unice.fr

Clément Cancès
email: cances@ann.jussieu.fr

Approximating the vanishing capillarity limit of two-phase flow in multi-dimensional heterogeneous porous medium

Keywords: two-phase flow, heterogeneous porous medium, discontinuous capillarity, vanishing capillarity, Buckley-Leverett model, oil trapping, admissible solution, finite volume approximation MSC2010 subject classification 76S05, 35L65, 65M08

come

Approximating the vanishing capillarity limit of two-phase flow in multi-dimensional heterogeneous porous medium

Introduction

Consider a flow in porous medium of a mixture of two immiscible phases (typically, oil and water). In this note, we start with a well-established parabolic-elliptic model that accounts for capillarity forces through a nonlinear function π(x, s) where x is a point location (the dependence on x accounts for physical properties of the medium) and s is the saturation of one of the phases. In this case, the flow is driven by convection (first-order) and diffusion/capillarity (second order) forces. Yet in large time and space scales which typically occur in geological applications, it is tempting to neglect the second-order effects. Formally, this results in a hyperbolic-elliptic system of equations. In one space dimension the elliptic equation is integrated explicitly, and the resulting Buckley-Leverett equation is a hyperbolic scalar conservation law. In the case of heterogeneous medium with sharp transitions between layers, both theory and numerical approximation of this equation were deeply investigated in recent years. In the present paper, our goal is to benefit both from the understanding of the coupling conditions at the interface between multidimensional porous media [START_REF] Cancès | An existence result for multidimensional immiscible two-phase flows with discontinuous capillary pressure field[END_REF] and the understanding of the relevant hyperbolic limit in the one-dimensional case [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF] in order to propose a consistent numerical strategy for the hyperbolic-elliptic model that appears in dimensions two and three. Starting from the the scheme exhibited in the work [START_REF] Brenner | Convergence of finite volume approximation for immiscible two-phase flows in porous media with discontinuous capillary pressure field in several dimensions[END_REF], we build and test numerical schemes that agree with the vanishing capillarity paradigm.

Let us first discuss the one-dimensional case following the references [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF][START_REF] Cancès | Asymptotic behavior of two-phase flows in heterogeneous porous media for capillarity depending only on space. I. Convergence to the optimal entropy solution[END_REF][START_REF] Cancès | Asymptotic behavior of two-phase flows in heterogeneous porous media for capillarity depending only on space. II. Nonclassical shocks to model oil-trapping[END_REF][START_REF] Kaasschieter | Solving the Buckley-Leverett equation with gravity in a heterogeneous porous medium[END_REF][START_REF] Kruzhkov | First order quasilinear equations with several independent variables[END_REF]. In these works, the Buckley-Leverett equation is seen as the singular limit of parabolic equation with capillarity or viscosity terms. The resulting solutions are characterized intrinsically as entropy solutions. In general, entropy solutions are weak (distributional) solutions satisfying additional integral inequalities that single out the physically relevant solution. In particular, in the fundamental case of a homogeneous medium, and also in the case of a heterogeneous medium with gradually changing physical properties, the classical Kruzhkov theory of entropy solutions [START_REF] Kruzhkov | First order quasilinear equations with several independent variables[END_REF] fully describes the vanishing capillarity limit (see, e.g., [START_REF] Kaasschieter | Solving the Buckley-Leverett equation with gravity in a heterogeneous porous medium[END_REF] and references therein).

The situation is very different when a jump heterogeneity is present in the medium, as this is the case at the interface between two different geological layers. In this case, we deal with a discontinuous-flux Buckley-Leverett equation (see [START_REF] Adimurthi | Godunov-type methods for conservation laws with a flux function discontinuous in space[END_REF][START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF][START_REF] Bürger | An Engquist-Osher-type scheme for conservation laws with discontinuous flux adapted to flux connections[END_REF][START_REF] Kaasschieter | Solving the Buckley-Leverett equation with gravity in a heterogeneous porous medium[END_REF] and references therein) that we see as the limit, as ε → 0, of the parabolic model φ(x) s t + f (x; s) x = ε(λ(x; s)π(x, s) x ) [START_REF] Adimurthi | Godunov-type methods for conservation laws with a flux function discontinuous in space[END_REF] with piecewise constant in x functions f , λ and π; we refer to Section 1 for details. The analysis of [START_REF] Cancès | On the effects of discontinuous capillarities for immiscible two-phase flows in porous media made of several rock-types[END_REF] (see [START_REF] Cancès | Asymptotic behavior of two-phase flows in heterogeneous porous media for capillarity depending only on space. I. Convergence to the optimal entropy solution[END_REF][START_REF] Cancès | Asymptotic behavior of two-phase flows in heterogeneous porous media for capillarity depending only on space. II. Nonclassical shocks to model oil-trapping[END_REF] for more details) shows that, according to the relation between the capillarity parameters on the two sides of the rock, two qualitatively different kinds of solution may appear. The optimal entropy conditions (see [START_REF] Adimurthi | Godunov-type methods for conservation laws with a flux function discontinuous in space[END_REF][START_REF] Kaasschieter | Solving the Buckley-Leverett equation with gravity in a heterogeneous porous medium[END_REF]) select the weak solution of the discontinuous-flux Buckley-Leverett equation that facilitates the migration of the phases across the interface between the two layers. The barrier entropy conditions of [START_REF] Cancès | Asymptotic behavior of two-phase flows in heterogeneous porous media for capillarity depending only on space. II. Nonclassical shocks to model oil-trapping[END_REF] select the weak solution of the same equation that minimizes the flux of the phases across the interface. Notice that this latter case can be seen as a rough mathematical model for the oil trapping phenomenon, responsible for the formation of natural oil reservoirs. Furthermore, it was shown in [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF] that infinitely many intermediate notions of entropy solution (accounting for a more or less stringent limitation of the flux conductivity of the interface) may appear at the vanishing capillarity limit. More exactly, given the convective fluxes f l,r (s) on the left (l) and on the right (r) from the interface, the appropriate entropy inequalities depend on the shape of capillary pressure functions π l,r (s). As a consequence, one may state that the physically relevant solution of the Buckley-Leverett equation takes into account the shape of (formally neglected) capillary pressure functions.

Roughly speaking, in the relevant solution the values of capillary pressure on each side from the interface should coincide or be "connected" by a suitable interface layer (see [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF] for details).

Due to this observation, the numerical approximation of the two-phase flow in heterogeneous media, at geological scales, becomes a particularly delicate issue. On the one hand, discretization of the full parabolic problem (as the one of [START_REF] Brenner | Convergence of finite volume approximation for immiscible two-phase flows in porous media with discontinuous capillary pressure field in several dimensions[END_REF]) may be too expensive with respect to the precision required. On the other hand, whenever the saturation is taken as the primal variable straightforward discretization of the singular limit must fail, because the formal equation at the limit does not contain any information on the capillary pressure functions; it becomes impossible to distinguish between "optimal", "barrier" and the intermediate kinds of solutions to the problem. Thus, in order to capture the physically relevant solution while discretizing the Buckley-Leverett equation, we have to complement a standard approximation strategy (such as fully implicit or IMplicit Pressure -Explicit Saturation (IMPES) scheme with time-explicit three-point finite volume approximation of the hyperbolic equation) with a specific procedure which goal is to connect the values of capillary pressure.

Such a scheme was developed in [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF] (see also [START_REF] Andreianov | The Godunov scheme for scalar conservation laws with discontinuous bell-shaped flux functions[END_REF]) in the one-dimensional case. Yet the multidimensional situation is different, because an elliptic equation for the pressure cannot be eliminated from the system. Therefore, we start with the multi-dimensional scheme of [START_REF] Brenner | Convergence of finite volume approximation for immiscible two-phase flows in porous media with discontinuous capillary pressure field in several dimensions[END_REF], set the capillarity parameter ε equal to zero, and introduce a series of simplifications aimed at reducing the computational complexity. The resulting schemes are compared numerically. Notice that, due to the low regularity of solutions of the equation giving fluid velocity, it does not seem possible to extend the theoretical study from the one-dimensional case to the multi-dimensional one. Thus, numerical evidence is the key ingredient of validation of the scheme, along with the convergence proof valid in the one-dimensional case only.

The paper is organized as follows. In Section 1 we describe the models. In Section 2 the information inferred from analysis of the one-dimensional case is outlined. Section 3 is the core of the paper: different versions of the finite volume scheme are presented. Section 4 provides numerical results.

Description of the models

In this paper, we consider a new model that is directly obtained from the one studied in [START_REF] Brenner | Convergence of finite volume approximation for immiscible two-phase flows in porous media with discontinuous capillary pressure field in several dimensions[END_REF][START_REF] Cancès | An existence result for multidimensional immiscible two-phase flows with discontinuous capillary pressure field[END_REF]; while we neglect the capillary diffusion within the rocks, we keep track of the singular effect linked to the capillary pressure discontinuity at the interface. In Section 1.1, we briefly present the case of non-zero capillary diffusion whose contribution away from inner interfaces is neglected in the model given in Section 1.2.

In the sequel we assume that the porous medium Ω ⊂ R d is made of two different rocks Ω 1 , Ω 2 separated by an interface Γ = ∂Ω 1 ∩ ∂Ω 2 . For the sake of simplicity, we assume that both Ω i are polygonal subdomains, and we restrict our study to the case of a finite time horizon T > 0.

The parabolic model

Our starting point is the model where immiscible incompressible two-phase flow in the homogeneous porous medium Ω i is governed by the coupling of the degenerate parabolic equation on the saturation

φ i ∂ t s + ∇ • (u t f i (s) + K i λ i (u)(-∇π i (u) + ρg)) = q o (s), (2)
with the uniformly elliptic equation

∇ • u t = q o (s) + q w (s), u t = -K i (M i (s)∇P -ζ i (s)g) . (3)
In equations ( 2) and (3), φ i ∈ (0, 1) denotes the porosity, s ∈ L ∞ (Ω × (0, T ); [0, 1]) is the saturation of the oil phase, u t is the total speed, the fractional flow f i (s) is an increasing function such that f i (0) = 0 and f i (1) = 1, K i is the intrinsic permeability tensor of the rock Ω i , λ i is a nonnegative function satisfying λ i (0) = λ i (1) = 0 and λ i (s) > 0 for s ∈ (0, 1), the capillary pressure function π i is increasing and belongs to L 1 (0, 1), ρ denotes the difference between the density of oil and the density of water, while g is the gravity. The oil and water source terms q o , q w satisfy q o (0) ≥ 0, q w (1) ≥ 0. The function M i (s) is such that M i (s) ≥ α > 0 for all s ∈ [0, 1], and ζ i (s) is a smooth function. The function P is the so called global pressure which is a kind of weighted mean of the phase pressure [START_REF] Chavent | Mathematical Models and Finite Elements for Reservoir Simulation: Single Phase, Multiphase, and Multicomponent Flows Through Porous Media[END_REF]. We stress here that since P is not a physical pressure, it can be discontinuous at the rock discontinuities.

Let us now focus on the transmission conditions across Γ. Due to conservation of mass, one has i=1,2

u t • n i = 0, i=1,2 (u t f i (s) + K i λ i (u)(-∇π i (u) + ρg)) • n i = 0. ( 4
)
It has been pointed out in [START_REF] Cancès | An existence result for multidimensional immiscible two-phase flows with discontinuous capillary pressure field[END_REF] that it is natural to extend the notion of phase pressure by considering π i as multivalued graphs, in order to give sense to the transmission conditions across rock discontinuities. Therefore, the capillary pressure has also to be considered as multivalued, yielding the introduction of the maximal monotone graph πi defined by

πi (s) = π i (s) if s ∈ (0, 1), πi (0) = [-∞, π i (0)] and πi (1) = [π i (1), +∞].
Requiring the "continuity" of the extended phase pressure turns to requiring the existence of a measurable function π : Γ × (0, T ) → R such that π ∈ π1 (s 1 ) ∩ π2 (s 2 ) and

P 1 -Z 1 (π) = P 2 -Z 2 (π), (5)
where

Z i (p) = p 0 f i • π-1 i (
a)da and s i , P i denote the traces of s, P on Γ×(0, T ) from the side of Ω i . The system (2)-( 5) is complemented by an initial condition on the saturation s(x, 0) = s 0 (x) ∈ [0, 1] for almost all x ∈ Ω, [START_REF] Andreianov | Finite volume schemes for locally constrained conservation laws[END_REF] and by boundary conditions. Here we consider the simplest choice, i.e. no-flux boundary conditions:

(u t f i (s) + K i λ i (u)(-∇π i (u) + ρg)) • n = u t • n = 0 on ∂Ω ∩ ∂Ω i × (0, T ). ( 7
)
We could deal with more general boundary conditions, but this is beyond the scope of this paper.

Let us point out that under suitable assumptions on the data, system (2)-( 7) admits at least one solution [START_REF] Cancès | An existence result for multidimensional immiscible two-phase flows with discontinuous capillary pressure field[END_REF], and it can be approximated by means of a fully implicit Finite Volume scheme of [START_REF] Brenner | Convergence of finite volume approximation for immiscible two-phase flows in porous media with discontinuous capillary pressure field in several dimensions[END_REF].

The singular limit as capillarity is neglected

As explained for example in [START_REF] Amaziane | Generalized nonequilibrium capillary relations for two-phase flow through heterogeneous media[END_REF][START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF][START_REF] Cancès | On the effects of discontinuous capillarities for immiscible two-phase flows in porous media made of several rock-types[END_REF], rock heterogeneities strongly influence the behavior of the flow, since a singular effect can be induced at the rock discontinuities. Therefore our purpose is to consider the capillarity only at the rock discontinuity, and to neglect it within the homogeneous domain where it plays a minor role. From a practical point of view, we slightly modify equation (2) to get

φ i ∂ t s + ∇ • (u t f i (s) + λ i (u)ρK i g) = q o (s), (8)
and we modify the flux continuity relations at the interface to get i=1,2

u t • n i = 0, i=1,2 (u t f i (s) + K i λ i (u)ρg) • n i = 0 on Γ × (0, T ). (9)
Finally, the boundary condition [START_REF] Andreianov | A theory of L 1 -dissipative solvers for scalar conservation laws with discontinuous flux[END_REF] takes the form

(u t f i (s) + λ i (u)ρK i g) • n = u t • n = 0 on ∂Ω ∩ ×(0, T ). ( 10
)
Notice that for a given total speed u t , equation ( 8) is a scalar conservation law. As noticed in [START_REF] Bürger | On the well-posedness of entropy solutions to conservation laws with a zero-flux boundary condition[END_REF], the no-flux boundary condition [START_REF] Bürger | An Engquist-Osher-type scheme for conservation laws with discontinuous flux adapted to flux connections[END_REF] is satisfied in a strong sense thanks to the fact that

u t • n = 0 and λ i (0) = λ i (1) = 0 for the flux function F defined by F(x, t; s) = u t (x, t)f i (s) + λ i (s)ρK i g.
In our simplified model, we should keep the formal jump conditions (5) so that the contribution of the capillarity is taken into account at the interface between the different rocks where it plays a crucial role. Let us stress that due to an interface layer that may develop as ε → 0, (5) may not hold for the solutions we construct. Yet, as shows the study of the one-dimensional case (Section 2), it is correct to prescribe (5) at the level of the numerical approximation; as the discretization step goes to zero, a numerical boundary layer may be observed. Except in Section 2, we do not investigate the precise meaning that should be given to [START_REF] Andreianov | The Godunov scheme for scalar conservation laws with discontinuous bell-shaped flux functions[END_REF]. Our goal is to solve numerically the problem described by equations ( 3),( 5),( 6),( 8)- [START_REF] Bürger | An Engquist-Osher-type scheme for conservation laws with discontinuous flux adapted to flux connections[END_REF].

Entropy solutions of the one-dimensional Buckley-Leverett model

In the one-dimensional case, equation ( 3) implies that u t does not depend on the space variable, so that the model reduces to a single scalar conservation law with discontinuous flux function for the sole unknown s (see [START_REF] Adimurthi | Godunov-type methods for conservation laws with a flux function discontinuous in space[END_REF][START_REF] Adimurthi | Optimal entropy solutions for conservation laws with discontinuous flux-functions[END_REF][START_REF] Andreianov | A theory of L 1 -dissipative solvers for scalar conservation laws with discontinuous flux[END_REF][START_REF] Bürger | An Engquist-Osher-type scheme for conservation laws with discontinuous flux adapted to flux connections[END_REF] and references therein for examples and general theory of such problems). Here we give a summary of the relevant analysis carried out in [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF], with emphasis on the numerical recipe of flux limitation and introduction of the pressure transmission parameter π.

In [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF], the Buckley-Leverett equation in piecewise homogeneous two-rocks' medium with interface situated at {x = 0} is seen as the singular limit, as ε → 0, of the problems the form (1) with

f (x, •) = f l (•), x < 0 f r (•), x > 0 , λ(x, •) = λ l (•), x < 0 λ r (•), x > 0 , π(x, •) = π l (•), x < 0 π r (•), x > 0 . (11)
Here, the phases are labeled in such way that f l,r are functions having at most one extremum within [0, 1] which is a local maximum. Recall that we have f l,r (0) = 0 and also f l (1), f r (1) coincide. For ε > 0, the weak traces of the total fluxes f l,r (s)ελ l,r (s)∂ x π l,r (s) at {x = 0} in (1) should match in order to ensure the conservation of the two phases. Moreover, the phase pressures should be continuous as long as the phase is present on both sides of the interface [START_REF] Cancès | An existence result for multidimensional immiscible two-phase flows with discontinuous capillary pressure field[END_REF][START_REF] Enchéry | Numerical approximation of a two-phase flow in a porous medium with discontinuous capillary forces[END_REF]; this leads to the trace relation πl (s(t, 0 -)) ∩ πr (s(t, ) + )) = ∅ with multivalued πl,r , as explained above. Existence of solutions for ε > 0 is shown in [START_REF] Cancès | Two-phase flows involving capillary barriers in heterogeneous porous media[END_REF] (see also [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF]) and compactness of the resulting family of solutions s ε comes under the additional assumption of strict nonlinearity of f l,r (the general case can be treated via a more specific BV loc argument, see [4, Appendix A.1]).

The limit s of (s ε ) ε satisfies the Kruzhkov entropy inequalities [START_REF] Kruzhkov | First order quasilinear equations with several independent variables[END_REF] away from the interface:

∂ t φ(x)|s -κ| + ∂ x Φ(x; s, κ) ≤ 0 in D ′ (R \ {0}) × (0, T ) , (12)
where

φ(x) = φ l 1 1 x<0 + φ r 1 1 x>0 and Φ(x; s, κ) = Φ l (s, κ)1 1 x<0 + Φ r (s, κ)1 1 x>0 , Φ l,r (•, κ) being the Kruzhkov entropy-fluxes sign (• -κ)(f l,r (•) -f l,r (κ))
. Further, due to the conservative weak formulation of the equation, the conservativity relation at x = 0 passes to the limit and yields the Rankine-Hugoniot condition f l (s(t, 0 -)) = f r (s(t, 0 + )) (for strictly nonlinear f l,r strong traces s(t, 0 ± ) exist due to Kruzhkov inequalities [START_REF] Cancès | Finite volume scheme for two-phase flow in heterogeneous porous media involving capillary pressure discontinuities[END_REF]; see [START_REF] Andreianov | A theory of L 1 -dissipative solvers for scalar conservation laws with discontinuous flux[END_REF] for more generality). While the relation connecting capillary pressures does not pass to the limit, the condition inherited as ε → 0 can be described implicitly by the specific global adapted entropy inequality:

∂ t φ(x)|s -c(x)| -∂ x Φ(x; s, c(x)) ≤ 0 in D ′ (R × (0, T )), (13)
with c(x

) := s π l , x < 0 s π r , x > 0. ( 14
)
Actually, the issue of interface coupling for the singular limit of (1) as ε → 0 reduces to the selection of the couple (s π l , s π r ) called connection, as described in [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF]. Two qualitatively different cases may occur:

(I) the fluxes are connected at the "optimal connection level":

f l (s π l ) = f r (s π r ) = F opt where F opt = min{max [0,1] f l , max [0,1] f r }, (15) with, moreover, f ′ l (s π l ) ≤ 0 and f ′ r (s π r ) ≥ 0. ( 16
)
In this case, the vanishing capillarity limits are the "optimal solutions" obtained by Kaasschieter [START_REF] Kaasschieter | Solving the Buckley-Leverett equation with gravity in a heterogeneous porous medium[END_REF]. One can observe that in this case, the pressures are not connected: π l (s π l ) = π r (s π r ). (II) both capillary pressures and fluxes are connected: [START_REF] Cancès | Two-phase flows involving capillary barriers in heterogeneous porous media[END_REF] holds, and πl (s π l ) ∩ πr (s π r ) = 0, and

f l (s π l ) = f r (s π r ) = F ≤ F opt . (17)
Notice that in both cases (I) and (II) we have [START_REF] Cancès | Two-phase flows involving capillary barriers in heterogeneous porous media[END_REF]. Except for very specific configurations of fluxes f l,r , the two inequalities are strict in the case (II), while only one of them is strict in the case (I). As a consequence, following the terminology of [START_REF] Adimurthi | Optimal entropy solutions for conservation laws with discontinuous flux-functions[END_REF][START_REF] Bürger | An Engquist-Osher-type scheme for conservation laws with discontinuous flux adapted to flux connections[END_REF] the connection (s π l , s π r ) consists in a strictly undercompressive discontinuity in case (II), while it is only marginally undercompressive in case (I).

Remark 2.1 The case (II) was judged exceptional and neglected in the pioneering work [START_REF] Kaasschieter | Solving the Buckley-Leverett equation with gravity in a heterogeneous porous medium[END_REF]; let us stress that both (I) and (II) can be relevant, depending on the shape of the nonlinearities f l,r and π l,r . We refer to [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF] for a detailed discussion, and to [START_REF] Cancès | Asymptotic behavior of two-phase flows in heterogeneous porous media for capillarity depending only on space. I. Convergence to the optimal entropy solution[END_REF][START_REF] Cancès | Asymptotic behavior of two-phase flows in heterogeneous porous media for capillarity depending only on space. II. Nonclassical shocks to model oil-trapping[END_REF] for examples.

To be precise, if (II) is possible then (s π l , s π r ) obeys (II), so that (I) only occurs when (II) cannot be realized. In practice, finding s π l,r corresponds to intersecting the graph of a strictly decreasing function with a non-decreasing graph. Hence the uniqueness of the intersection point in (II) follows, while (I) happens when existence fails. In addition, one can observe that among all undercompressive (in the above sense) shocks with connected fluxes, conditions (I),(II) minimize the gap |π l (s π l )π r (s π r )|, which is an intuitively appealing feature. We refer to [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF] for the study of elementary vanishing capillarity profiles that provides the justification of the selection rule (I),(II) .

As shown in [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF], the function c(x) in ( 14) chosen according to (I),(II) appears as an explicit vanishing capillarity limit for (1), i.e., it is an admissible solution for a very special initial datum. It is remarkable that the knowledge of (c π l , c π r ) and the rigid structure of the semigroups of solutions of (1) permit to characterize vanishing capillarity limits for all initial data. Indeed, we have Theorem 2.1 For every initial datum u 0 ∈ L ∞ (R) there exists a unique u ∈ C([0, T ]; L 1 loc (R)) weak solution of (1), [START_REF] Buzzi | Interface conditions for degenerate two-phase flow equations in one space dimension[END_REF] with ε = 0 that satisfies the Kruzhkov entropy inequalities (12) away from the interface and the global adapted entropy inequality (13) with s π l,r in [START_REF] Cancès | Asymptotic behavior of two-phase flows in heterogeneous porous media for capillarity depending only on space. II. Nonclassical shocks to model oil-trapping[END_REF] given by the selection rule (I),(II). Moreover, the limit, as ε → 0, of solutions u ε of (1) with datum u 0 exists and coincides with u.

This means that the property of being vanishing capillarity limit (which is very difficult to check directly) is fully characterized by entropy and well-chosen adapted entropy inequalities.

It is easily seen that extending such analysis to the multi-dimensional situation does not seem possible with respect to the existing state of the art. But some insight can be carried from the one-dimensional situation to the general case. Specifically, let us insist on the role of the parameter F which takes values in the interval [ F bar , F opt ]. It provides the connection level, which turns out to be the most relevant physical parameter at the vanishing capillarity limit. While the value F in (II) indicates the level of the interface flux for the particular solution c(x) given by ( 14), it is shown in [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF] that F is the largest possible value of the interface flux f l (s(t, 0 -)) = f r (s(t, 0 + )) for any solution of our singular limit problem. Thus the optimal entropy solution corresponding to F opt and case (I) optimizes (maximizes) the capacity of the interface to conduct the flow of the phases; the barrier entropy condition corresponding to F bar minimizes the flow, which can be seen as a model for the oil trapping phenomenon.

The value F depends on f l,r but also on π l,r that are formally absent from (1) when ε = 0. Let us denote it F π , in order to stress this dependence. The value F π is determined by an easy algorithm solving a monotone nonlinear equation on [0, 1]. It can be seen as the physical parameter that records the action of the capillary diffusion at the interface between the two rocks. But F π can also be used, in a strikingly simple way, in the numerical simulation in order to guarantee the convergence of approximations to the physically correct solution. Indeed, as explained in [4, Prop. 9] (see also [START_REF] Andreianov | The Godunov scheme for scalar conservation laws with discontinuous bell-shaped flux functions[END_REF]), it is enough to take the Godunov numerical flux G suitable for the approximation of the optimal solution (see [START_REF] Adimurthi | Godunov-type methods for conservation laws with a flux function discontinuous in space[END_REF]) and consider F int : (z l , z r ) → min F π , G(z l , z r ) , with the optimal flux G : (z l , z r ) → min f l min{z l , sl } , f r max{s r , z r } [START_REF] Chavent | A Finite-Element simulator for incompressible two-phase flow[END_REF] where sl,r := argmax [0,1] f l,r . Finally, let us stress that the flux limitation in the Godunov scheme can be achieved by a less explicit method that nevertheless has a clear numerical interpretation; we will use it thoughout Section 3. Denote by G l,r the Godunov fluxes (see, e.g., [START_REF] Adimurthi | Godunov-type methods for conservation laws with a flux function discontinuous in space[END_REF][START_REF] Andreianov | The Godunov scheme for scalar conservation laws with discontinuous bell-shaped flux functions[END_REF][START_REF] Eymard | Finite Volume Methods. Handbook of Numerical Analysis[END_REF]) associated with f l,r , respectively; they can be computed using sl,r .

Lemma 2.2 (see [4, Lem. 8 and Prop. 9]) Let (z l , z r ) ∈ [0, 1] 2 , then the system

π ∈ πl ∩ πr G l (z l , π -1 l (π)) = G r (π -1 r (π), z r ) ( 19
)
admits at least one solution π ∈ R. The value F int (z l , z r ) := G l (z l , π -1 l (π)) = G r (π -1 r (π), z r
) is defined uniquely by [START_REF] Chavent | Mathematical Models and Finite Elements for Reservoir Simulation: Single Phase, Multiphase, and Multicomponent Flows Through Porous Media[END_REF], furthermore, it coincides with the value computed by the explicit formula [START_REF] Chavent | A Finite-Element simulator for incompressible two-phase flow[END_REF].

The lemma introduces an additional variable π on the interface that serves to connect the pressures in the vanishing capillarity limit; while it is rather useless for the one-dimensional problem, the value π will be instrumental in the multi-dimensional schemes of Section 3. Following closely the proof of [4, Prop. 9] one finds the algorithm to determine the value π without solving [START_REF] Chavent | Mathematical Models and Finite Elements for Reservoir Simulation: Single Phase, Multiphase, and Multicomponent Flows Through Porous Media[END_REF]; but this is not our goal. Heading towards a simpler practical implementation, let us generalize [START_REF] Chavent | Mathematical Models and Finite Elements for Reservoir Simulation: Single Phase, Multiphase, and Multicomponent Flows Through Porous Media[END_REF] as follows.

Proposition 2.3 If the Godunov fluxes G l,r are replaced in the scheme of [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF] (including equation [START_REF] Chavent | Mathematical Models and Finite Elements for Reservoir Simulation: Single Phase, Multiphase, and Multicomponent Flows Through Porous Media[END_REF] for F int ) by consistent with f l,r monotone Lipschitz continuous numerical fluxes R l,r (see, e.g., [START_REF] Eymard | Finite Volume Methods. Handbook of Numerical Analysis[END_REF]), then the scheme still converges to the solution given by Theorem 2.1.

Proof: It is easily checked that F int (a, b) remains increasing in a and decreasing in b. Thus the scheme is monotone, and due the the Crandall-Tartar lemma, discrete L 1 contraction inequalities follow. Therefore the convergence arguments as in [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF]Prop. 6] apply provided we show that the solution (14) of equation ( 1) avec ε = 0 is a limit of some solutions of the scheme. This is obvious in case (II) because ( 14) gives a stationary solution to the scheme, due to the consistency of fluxes R l,r . In case (I), we use a numerical counterpart of the construction of [4, Prop. 3(ii)]: using the monotonicity of the numerical flux, one constructs an appropriate stationary numerical profile. E.g., the situation where s π l = sl we take the value sπ l ∈ [0, sl ] such that π l (s π l ) = π r (s π r ). There exists a solution to the scheme

u = (u i ) 0 i=-∞ such that u 0 = sπ l , R l (u -1 , u 0 ) = f r (s π r ) and ∀i ≤ -1 R l (u i-1 , u i ) = R l (u i , u i+1 ),
that we complement by the constant values u i := s π r for i > 0. The profile is monotone in i thanks to the monotonicity of R l , hence there exists ℓ := lim i→-∞ u i . Due to the continuity of R l , from the above scheme we obtain that R

l (ℓ, ℓ) = lim i→-∞ R l (u i-1 , u i ) = f r (s π r ), which yields ℓ = s opt l . By a simple scaling argument, u h (x) = i∈Z u i 1 1 [(i-1)h,ih] converges, as h → 0, to the function (14).
Remark 2.2 A particular case of construction [START_REF] Chavent | Mathematical Models and Finite Elements for Reservoir Simulation: Single Phase, Multiphase, and Multicomponent Flows Through Porous Media[END_REF] of the interface fluxes appears, under the form

F int (z l , z r ) := G l (z l , s) = G r (s, z r ) for some s ∈ [0, 1], ( 20
)
in the early work of Chavent, Cohen and Jaffré [START_REF] Chavent | A Finite-Element simulator for incompressible two-phase flow[END_REF]. The choice of numerical fluxes [START_REF] Diehl | A uniqueness condition for nonlinear convection-diffusion equations with discontinuous coefficients[END_REF] means that connection of the saturation values is expected at the interface; that is, π l ≡ π r should be assumed in our model. Therefore the so obtained solutions correspond to the mathematical vanishing viscosity method. Notice that (20) was also used by Diehl [START_REF] Diehl | A uniqueness condition for nonlinear convection-diffusion equations with discontinuous coefficients[END_REF] to characterize vanishing viscosity limits for general conservation laws with discontinuous flux (see [START_REF] Diehl | A uniqueness condition for nonlinear convection-diffusion equations with discontinuous coefficients[END_REF]Section 3]). Let us stress that if π l ≡ π r , the vanishing viscosity regularization of the Buckley-Leverett equation leads, in general, to non-physical solutions.

Finite volume approximation at vanishing capillarity limit

We propose three finite volume schemes designed for approximating problem (3), ( 5), ( 6), ( 8)- [START_REF] Bürger | An Engquist-Osher-type scheme for conservation laws with discontinuous flux adapted to flux connections[END_REF]. The first one is described in Section 3.2. It is obtained by neglecting the inner capillary diffusion in the scheme proposed in [START_REF] Brenner | Convergence of finite volume approximation for immiscible two-phase flows in porous media with discontinuous capillary pressure field in several dimensions[END_REF], which yields a fully implicit (and thus fully coupled) scheme whose unknowns are the saturations (s n+1 K ) and the pressures (P n+1 K ) in the control volumes as well as the capillary pressure transmission parameters (π n+1 σ ) on the edges σ lying in the inner interface Γ. The main advantage of this scheme is that it preserves energy bounds, ensuring some stability.

Its main drawback is of course its computational cost, that we aim to reduce in the variants of the scheme. The first alternative we propose in Section 3.3 consists in the IMPES scheme decoupling the saturations (s n+1 K ) from the pressure and capillary pressure unknowns (P n+1 K ) and (π n+1 σ

). Further, in Section 3.4 we propose an algorithm that fully decouples the saturations (s n+1 K ), the pressures (P n+1 K ) and the interface capillary pressures (π n+1 σ

). For the sake of simplicity, we suppose in the sequel that the porous medium is isotropic, i.e. K i = k i I d , so that the elliptic equation ( 3) can be approximated by means of a two-point flux FV-scheme. The consistency of such a method requires the so-called orthogonality condition on the mesh that is recalled in Section 3.1. The natural extensions to the case of anisotropic porous media and to more general meshes are not discussed in this paper.

Discretization of Ω × (0, T )

Since the scheme used for discretizing the elliptic equation is based on a two-point flux approximation, the space discretization has to fulfill some admissibility condition. The notion of admissible discretization of Ω is inspired from the one given in [START_REF] Eymard | Finite Volume Methods. Handbook of Numerical Analysis[END_REF], and additionally requires that each cell K is either in Ω 1 or in Ω 2 . We note T i the set of all the cells included in Ω i , so that Ω i = K∈Ti K i . As a consequence, the interface Γ is made of a union of a subset of edges denoted by E Γ , i.e. Γ = σ∈EΓ σ. For each cell K ∈ T , there is a point x K ∈ K called center of K such that, denoting by K|L the edge between K and its neighboring cell L, the straight line (x K x L ) is orthogonal to K|L. We refer to [8, Definition 2.1] for a complete definition of the notion of admissible discretization of Ω that we use here.

The size of the mesh is h = max K∈T diam(K); we ask the usual regularity assumption for a family of meshes, namely, that there exist α > 0 such that |K| ≥ αh d , |∂K| = σ∈EK |σ| ≤ 1 α h d-1 for K ∈ T . The time discretization consists of a sequence of time steps ∆t n ) 0≤n≤N which have to fulfill a CFL condition when explicit schemes are considered. We set t 0 = 0, t n+1 = t n + ∆t n and assume t N = T .

The coupled implicit scheme

The scheme proposed now consists in fully neglecting the capillary diffusion (that is, merely setting ε = 0) in the scheme proposed and validated in [START_REF] Brenner | Convergence of finite volume approximation for immiscible two-phase flows in porous media with discontinuous capillary pressure field in several dimensions[END_REF] for ε > 0. Let us stress that for the onedimensional case, the scheme of [START_REF] Brenner | Convergence of finite volume approximation for immiscible two-phase flows in porous media with discontinuous capillary pressure field in several dimensions[END_REF] is asymptotic-preserving, i.e., as ε → 0 it converges a numerical scheme that does approximate the entropy solution described in Section 2 (see [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF]Cor. 10]).

The first equation of ( 3) is approximated by

σ∈EK |σ|u n+1 K,σ = |K|q o (s n+1 K ) + |K|q w (s n+1 K ) ∀K ∈ T , (21)
where

u n K,σ is an approximation of 1 ∆t n |σ| t n+1 t n σ u t • n K,σ dxdt,
where n K,σ denotes the normal to σ outward w.r.t. K. The source term q o is supposed to be decreasing and such that q o (0) = 0, while the source term q w is supposed to be increasing and such that q w (1) ≥ 0. The quantity u n+1 K,σ is given by

u n+1 K,σ = 0 if σ ⊂ ∂Ω, u n+1 K,σ = k i M n+1 K|L P n+1 K -P n+1 L d(x K , x L ) + ζ n+1 K|L if σ = K|L ∈ E i , (22)
where

M n+1 K|L = MK (s n+1 K )ML(s n+1 L )d(xK ,xL) MK (s n+1 K )d(xK,σ)+ML(s n+1 L )d(xL,σ) and ζ n+1
K,L denotes the Rusanov numerical flux for the flux function ζ i (s)g • n K,σ and for the Riemann problem with initial data s n+1 K 1 1 x<0 + s n+1 L 1 1 x>0 . The case σ = K|L ∈ E Γ is treated by means of the introduction of two new unknowns P n+1 K,σ and P n+1 L,σ to enforce the continuity of the flux (4) and the pressure condition ( 5) at the discrete level. Denoting by g K,σ = g • n K,σ , we prescribe at the interface that for all σ = K|L ∈ E Γ ,

u n+1 K,σ + u n+1 L,σ = 0 with u n+1 J,σ = k J M J (s n+1 J ) P n+1 J -P n+1 J,σ d(x J , σ) + ζ J (s n+1 J )g J,σ , J = K, L. (23)
Moreover, the interface pressures P n+1 K,σ and P n+1 L,σ have to fulfill

P n+1 K,σ -Z K (π n+1 σ ) = P n+1 L,σ -Z L (π n+1 σ ), (24)
where the unknowns (π n+1 σ ) σ∈EΓ are used to ensure the pressure transmission. Notice that for every σ ∈ E Γ , unknowns P n+1 K,σ and P n+1 L,σ are readily eliminated thanks to the linear relations ( 23)-( 24). The initial datum s 0 is discretized by setting

s 0 K = 1 |K| K s 0 (x)dx. Equation (8) is discretized by s n+1 K -s n K ∆t n |K| + σ∈EK |σ|F n+1 K,σ = |K|q o (s n+1 K ), (25)
F n+1 K,σ = 0 if σ ⊂ ∂Ω, F n+1 K,σ = R F n+1 K,σ ; s n+1 K , s n+1 L if σ = K|L ∈ E i , (26)
where

F n+1 K,σ denotes the function s → f K (s)u n+1 K,σ + k K ρλ K (s)g K,σ
and where R(F ; s L , s R ) denotes the Rusanov flux corresponding to the flux function F and to left and right states given by s L , s R . To handle the fluxes across the interfaces σ ∈ E Γ , additional unknowns π n+1 σ introduced in (24) are exploited. Namely, for σ = K|L ∈ E Γ , the flux F n+1 K,σ is defined by

F n+1 K,σ + F n+1 L,σ = 0 where for J = K and J = L, F n+1 J,σ = R F n+1 J,σ ; s n+1 J , π-1 J (π n+1 σ ) . (27)
Relations (27) yield one equation per unknown π n+1 σ , σ ∈ E Γ , thus the system is closed. Notice that [START_REF] Kaasschieter | Solving the Buckley-Leverett equation with gravity in a heterogeneous porous medium[END_REF],( 27) are directly inspired by the approach (19) of the one-dimensional case, but where the Godunov flux has been replaced by the Rusanov flux by virtue of Proposition 2.3.

The main advantage of this approximation is that it is unconditionally stable with respect to the time discretization. The following proposition can be proved by following the method proposed in [START_REF] Brenner | Convergence of finite volume approximation for immiscible two-phase flows in porous media with discontinuous capillary pressure field in several dimensions[END_REF]. Proposition 3.1 For all admissible discretization of Ω, and for all n ∈ {0, . . . , N -1}, there exists

(s n+1 K ) K∈T , (P n+1 K ) K∈T , (π n+1 σ ) σ∈EΓ ∈ [0, 1] #T × R #T × R #EΓ solution to the scheme (21)-(27).
It is possible, following the sketch of the proof of Proposition 3.5 in [START_REF] Brenner | Convergence of finite volume approximation for immiscible two-phase flows in porous media with discontinuous capillary pressure field in several dimensions[END_REF], to derive the following estimate on the discrete L 2 ((0, T ); H 1 (Ω i )) semi-norm (cf. [START_REF] Enchéry | Numerical approximation of a two-phase flow in a porous medium with discontinuous capillary forces[END_REF][START_REF] Eymard | Finite Volume Methods. Handbook of Numerical Analysis[END_REF]):

i=1,2 N -1 n=0 ∆t n σ=K|L∈Ei |σ| (P n+1 K -P n+1 L ) 2 d(x K , x L ) ≤ C. (28)
From this estimate, we can prove that the piecewise constant approximate pressure P D converges weakly in L 2 (Ω i × (0, T )) towards some pressure P ∈ L 2 ((0, T ); H 1 (Ω i )). We also deduce from the fact that 0 ≤ s n+1 K ≤ 1 for all K ∈ T that the piecewise constant approximate saturation s D converges in the L ∞ (Ω × (0, T )) weak-⋆ sense towards a function s in L ∞ (Ω × (0, T )), with 0 ≤ s ≤ 1 a.e. in Ω × (0, T ). This is obviously not sufficient to pass to the limit in the nonlinearities if the problem. Nevertheless, a complete convergence study of the scheme presented above has been performed in [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF] in the one-dimensional framework, where it was possible to establish rigorous convergence results.

The (P, π)-s partially decoupled explicit scheme

As it will appear in the Section 4 where the numerical results are exhibited, the scheme proposed in Section 3.2 is quite expensive due to the fact that it is fully coupled and implicit. Our goal is then to uncouple some discrete equations, introducing an IMPES (IMplicit Pressure -Explicit Saturation) discretization of equations ( 3), ( 5), ( 6), ( 8)- [START_REF] Bürger | An Engquist-Osher-type scheme for conservation laws with discontinuous flux adapted to flux connections[END_REF]. To do so, given s n K K∈T , we compute first the pressure variables (P n+1 K ) K∈T , (π n+1 σ ) σ∈EΓ , and then we use these new values of the pressures to compute s n+1 K K∈T by means of an explicit finite volume scheme. a) Pressure step. We assume that s n K K∈T is given in [0, 1] #T . The system we intend to solve is then

σ∈EK |σ|u n+1 K,σ = |K|q o (s n K ) + |K|q w (s n K ) ∀K ∈ T , (29)
u n+1 K,σ = 0 if σ ⊂ ∂Ω, u n+1 K,σ = k i M n K|L P n+1 K -P n+1 L d(x K , x L ) + ζ n K|L if σ = K|L ∈ E i , (30)
u n+1 K,σ + u n+1 L,σ = 0, u n+1 J,σ = k J M J (s n J ) P n+1 J -P n+1 J,σ d(x J , σ) + ζ J (s n J )g J,σ , J = K,L, if σ ∈ E Γ . (31)
while equation [START_REF] Kaasschieter | Solving the Buckley-Leverett equation with gravity in a heterogeneous porous medium[END_REF] still has to hold. After eliminating P n+1 K,σ and P n+1 L,σ thanks to [START_REF] Kaasschieter | Solving the Buckley-Leverett equation with gravity in a heterogeneous porous medium[END_REF] and (31), we obtain a system of #T equations (29) for the #T + #E Γ unknowns P n+1 K K∈T , π n+1 σ K∈T . In order to close the system we impose that, for all σ ∈ E Γ , Equation ( 8) is now discretized by

F n K,σ + F n L,σ = 0 where for J = K and J = L, F n J,σ = R F n+1 J,σ ; s n J , π-1 J (π n+1 σ ) . (32
s n+1 K -s n K ∆t n |K| + σ∈EK |σ|F n K,σ = |K| -L qo s n+1 K + q o (s n K ) + L qo s n K , (33)
b) Capillary pressure step. The previous step also yields the values u n+1 K,σ , thus the expressions of the fluxes F n+1 K,σ are made available. Now we can compute the values π n+1 σ for each σ ∈ E Γ by solving

F n K,σ + F n L,σ = 0 with F n J,σ = R F n+1 J,σ ; s n J , π-1 J (π n+1 σ ) for J = K, L (38)

Numerical experiments

All the three schemes of Section 3 preserve the physical bounds on the saturation s. It seems much more delicate to achieve stability estimates like (28), especially for the fully decoupled scheme, and convergence analysis is beyond our reach. Therefore we now switch to numerical tests. Below, the fully decoupled scheme is used in all the test cases; the fully implicit and the IMPES schemes actually yield quite close numerical solutions but a longer computational time.

4.1 Comparison with 1D reference solution obtained by a convergent scheme.

Let Ω = (-1, 1) 2 be 2-dimensional domain composed of two horizontal layers Ω 1 = {(x, y) ∈ Ω|y < 0} and Ω 2 = {(x, y) ∈ Ω|y > 0}, which corresponds to different rock-types. We prescribe a zero total flux on the whole boundary ∂Ω. The boundary conditions for the saturation equation are given by s = 0.5 on Γ 1 = {(x, y) ∈ ∂Ω s.t. |y| = 1}, (u t f i (s) + λ i (u)K i g) • n = 0 on Γ 2 = {(x, y) ∈ ∂Ω s.t. |x| = 1}, the initial saturation is set to s 0 = 0.5. We set g = -9.81e y , µ o = 1, µ w = 3, ρ o = 0.87, ρ w = 1, and

f = s s -µ o /µ w (1 -s) , λ = (ρ o -ρ w ) s(1 -s) µ w s -µ o (1 -s) , ζ = ρ o µ o s + ρ w µ o (1 -s)
The porosity φ is constant and set to φ = 1; the absolute permeability K is given by K 1 = 10, K 2 = 5. We take capillary pressure curves π i (s) := π i,0ln(1-s) that differ only by the value of the entry pressure: π 1,0 is set to 0, while π 2,0 takes values 0.5 (Test 1a) or 2 (Test 1b). The flow is driven by buoyancy and by the spatial discontinuity: roughly speaking, the second rock has a smaller pores, so that it is more attractive for water. The total flow u t being zero, the oil phase is displaced from the rock Ω 2 . This is actually the 1D test borrowed from [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF]; there exists an x-independent solution. However, using triangular meshes we ensure that the scheme we test is truly 2D. Indeed, discrete solution is computed by the fully decoupled IMPES scheme on structured triangular meshes with different number of volumes. It is compared to the reference 1D numerical solution evaluated using 1000 grid points (this is indeed a reference solution, due to Th. 2.1 and Prop. 2.3). We present the results in Fig. 1 (Test 1a), Fig. 2 (Test 1b). In each case, 2D schemes for two different choices of numerical fluxes (namely the Godunov and the Rusanov ones) are compared. For the final time T = 0.2, the error in the norm L 2 ((0, T ) × Ω) is plotted at the right of Fig. 1,2. We observe convergence rates close to 1 (for the Godunov flux) and 0.5 (the Rusanov flux). Notice that in Test 1a the two discretizations exhibit quite different behavior. In fact, for both choices of fluxes the connecting capillary pressure π is such that π -1 2 (π) > 0. Thus, a small amount of oil propagating in the domain Ω 2 is due only to the extra diffusion introduced by Rusanov flux. This is not the case for the Godunov flux, which provides the exact solution to the Riemann problem.

A 2D flow constrained by capillary barriers

We now turn to a truly 2D test case inspired from [23, test case 4.3]. The 2D domain Ω (see Fig. 3, top left), which mostly consists of rock Ω 1 , is initially saturated in water: s 0 = 0. The flow is constrained by a presence of two barriers (rock Ω 2 ) having a higher entry pressure (π 0,1 = 0, π 0,2 = 1.5). The vertical boundaries are assumed to be impermeable, so on Γ N zero flux is prescribed for both phases. At the bottom and the top of Ω we prescribe a constant rate of a total flux, i.e

u t • n = -1 on Γ in , u t • n = 1 on Γ out .
The constant saturation value s = 0.5 is imposed on Γ in . We use the same set of data as previously, except for the absolute permeability K 1,2 now assumed to be equal to 1 in both sub-domains. We use here the Rusanov scheme as an approximate Riemann solver. Computational domain and the oil saturation field at time t = 0.075, 0.2, 0.4, 0.6, 0.8. Even if both rocks Ω i have the same intrinsic permeability (so that the flux function is continuous across the interface), it appears that some oil remains trapped under the capillary barriers. Indeed, due to singular capillary effects at rock discontinuities, the relevant solution is not the usual entropy solution [START_REF] Andreianov | Vanishing capillarity solutions of Buckley-Leverett equation with gravity in two-rocks' medium[END_REF]. A flux limitation [START_REF] Andreianov | Finite volume schemes for locally constrained conservation laws[END_REF] is active at the lower boundary of the barrier, hampering the progression of the oil phase. This singular effect depends on the capillary pressure functions π i , despite the capillary diffusion has been neglected within the rocks Ω i .

Figure 1 :

 1 Figure 1: Test 1a (left to right): Oil saturation profile at t = 0.1 and t = 0.2 (Godunov flux); Godunov/Rusanov comparison; Error in L 2 , red for Godunov (the bottom line) and black for Rusanov; dashed lines: slopes 1 and 1/2.

Figure 2 :

 2 Figure 2: Test 1b (left to right): Oil saturation profile at t = 0.1 and t = 0.2 (Godunov flux); Godunov/Rusanov comparison; Error in L 2 , red for Godunov (bottom line) and in black for Rusanov; dashed lines: slopes 1 and 1/2.

Figure 3 :

 3 Figure 3: Test 2.Computational domain and the oil saturation field at time t = 0.075, 0.2, 0.4, 0.6, 0.8.

  (we have one nonlinear equation per edge σ ∈ E Γ , solved separately for each σ).

		Existence of a
	solution π n+1 σ	and uniqueness of the resulting fluxes F n K,σ follow by a monotonicity argument, as
	for Lemma 2.2.
	c) Saturation step. The saturations (s n+1 K ) K∈T are computed explicitly by equations (33) and
	the fluxes readily computed with (38) and (34).

Acknowledgement. The work of the first author was partially supported by the French ANR project CoToCoLa. The work of the second and third authors was partially supported by GNR MoMaS, CNRS-2439 (PACEN/CNRS, ANDRA, BRGM, CEA, EDF, IRSN).

where L qo denotes the Lispchitz constant of q o , and where

(in addition, recall that the fluxes F n K,σ for σ ∈ E Γ have already been computed during the pressure step in (32)). A CFL condition is mandatory to ensure stability of the scheme at this step. Denoting by L fi (resp. L λi ) the Lispchitz constant of f i (resp. λ i ), then we require that

, and given

The above result follows from the fact that equation (33) can be rewritten as

the right-hand side of the above expression being a nondecreasing function with respect to (w.r.t.) each of the variables s n K , K ∈ T . Notice that the interface fluxes F n K,σ given by (32) depend in a monotone and Lipschitz continuous way on s n K and s n L because, thanks to Lemma 2.2, they can be rewritten under the form [START_REF] Chavent | A Finite-Element simulator for incompressible two-phase flow[END_REF] on which monotonicity and Lipschitz continuity is obvious. Remark also that the discretization of the source term in (33) is such that it is linear w.r.t. s n+1 K (thus the scheme remains explicit), and so that it is increasing w.r.t. s n+1 K and nonincreasing w.r.t. s n K . This yields a stabilization of the scheme in the sense that the CFL condition (35) does not depend on the source term q o .

The fully decoupled explicit scheme

Notice that in Pressure Step of the above IMPES scheme, relations [START_REF] Kaasschieter | Solving the Buckley-Leverett equation with gravity in a heterogeneous porous medium[END_REF] and (32) involving (π n+1 σ ) σ∈EG are fully coupled. In particular, the expressions of the flux functions F n+1 K,σ needed to compute F n K,σ depend on the whole set of values of (π n+1 σ ) σ∈EG via the values u n+1 K,σ . Here, we propose a variant of the previous scheme where the computation of the capillary pressure transmission parameters π n+1 σ σ∈EΓ is decoupled from the computation of the global pressures P n+1 K K∈T . For this scheme, we need to initialize the capillary pressure. To do so, we build π 0 σ as the linear interpolation between the neighboring cell values of the capillary pressure:

a) Global pressure step. Values s n K K∈T and (π n σ ) σ∈EΓ being given, we obtain P n+1 K K∈T as solution to (29)-(31), closed by the interface relation

This system is linear and one can check that it admits a unique solution P n+1 K K∈T ∈ R #T .


Keywords: coloration détectable, NP-complétude, graphe biparti, graphe cubique

Indice de détection des graphes bipartis et des graphes cubiques

Résumé : Pour un graphe connexe G d'ordre |V (G)| ≥ 3 et un k-étiquetage c : E(G) → {1, 2, . . . , k} des arêtes de G, le code d'un sommet v de G est le k-uplet (ℓ 1 , ℓ 2 , . . . , ℓ k ), où ℓ i est le nombre d'arêtes incidentes à v qui sont étiquetées i. Le k-étiquetage c est détectable si, quels que soient deux sommets adjacents de G, leurs codes sont distincts. Le plus petit entier strictement positif k pour lequel G a un k-étiquetage détectable est l'indice de détection det(G) de G. Dans ce rapport, nous montrons qu'il est NP-complet de décider si l'indice de détection d'une graphe cubique vaut 2. Nous montrons également que l'indice de détection de tout graphe biparti de degré minimum au moins 3 est au plus 2. Enfin, nous donnons des conditions suffisantes pour qu'un graphe cubiquesoit d'indice de détection 3.

Introduction

For graph-theoretical terminology and notation, we in general follow [START_REF] Balakrishnan | A Textbook of Graph Theory[END_REF]. In this paper, we assume that the graphs G in discussion are finite, connected, undirected and simple with order |V (G)| ≥ 3. Let c : E(G) → {1, 2, . . . , k} be a labelling of the edges of G, where k is a positive integer. The color code of a vertex v of G is the ordered k-tuple code c (v) = (ℓ 1 , ℓ 2 , . . . , ℓ k ), where ℓ i is the number of edges incident with v that are labelled i for i ∈ {1, 2, . . . , k}. Therefore, (v), the degree of v in G. The labelling c is called a detectable coloring of G if any pair of adjacent vertices of G have distinct color codes. The detection number or detectable chromatic number of G, denoted det(G), is the minimum positive integer k for which G has a detectable k-coloring. We call G k-detectable if G has a detectable k-coloring.

ℓ 1 + ℓ 2 + • • • + ℓ k = d G
The concept of detection number was introduced by Karoński, Luczak and Thomason in [START_REF] Karoński | Edge weights and vertex colours[END_REF], inspired by the basic problem in graph theory that concerns finding means to distinguish the vertices of a connected graph and to distinguish adjacent vertices of a graph, respectively, with the minimum number of colors. For a survey on vertex-distinguishing colorings of graphs, see [START_REF] Escuadro | Vertex-distinguishing colorings of graphs: A survey of recent developments[END_REF].

In [START_REF] Karoński | Edge weights and vertex colours[END_REF], Karoński, Luczak and Thomason conjectured that det(G) ≤ 3. In [START_REF] Addario-Berry | Vertex colouring edge partitions[END_REF], Addario-Berry, Aldred, Dalal and Reed proved that: (i) det(G) ≤ 4 and (ii) if χ(G) ≤ 3, then det(G) ≤ 3. However, as observed by Khatirinejad et al. [START_REF] Khatirinejad | Vertex-colouring edge-weightings with two edge weights[END_REF], it seems NP-complete to decide if a graph is 2-detectable. Conjecture 1.1 (Khatirinejad et al. [START_REF] Khatirinejad | Vertex-colouring edge-weightings with two edge weights[END_REF]). It is NP-complete to decide whether a given graph is 2-detectable.

As an evidence, Dudek and Wajc [START_REF] Dudek | On the complexity of vertex-colouring edge-weightings[END_REF] showed that closely related problems are NP-complete. In Section 2, we settle this conjecture by showing that deciding if a cubic graph is 2-detectable is an NP-complete problem.

On the other hand, Khatirinejad et al. [START_REF] Khatirinejad | Vertex-colouring edge-weightings with two edge weights[END_REF] believed that for a given bipartite graph, deciding if it is 2-detectable should be easy. For

m 1 ≤ m 2 ≤ • • • ≤ m d , let Θ(m 1 , .
. . , m d ) be the graph constructed from d internally disjoint paths between distinct vertices x and y, in which the i-th path has length m i . Such a graph is called a Theta and the two vertices {x, y} are its poles. It is bad if m 1 = 1 and m i ≡ 1 mod 4 for all 2 ≤ i ≤ d. Khatirinejad et al. [START_REF] Khatirinejad | Vertex-colouring edge-weightings with two edge weights[END_REF] proved that a Theta is 2-detectable if and only if it is not bad, and asked whether all bipartite graphs except the bad Thetas were 2-detectable. This was answered in the negative by Davoodi and Omoomi [START_REF] Davoodi | On the 1-2-3-conjecture[END_REF] who gave a new family of non-2-detectable bipartite graphs, the Theta trees. A Theta tree is a graph obtained from a tree T by replacing each vertex t of V (T ) by a bad Theta with poles u t and v t and every edge st of E(T ) by a path P st of length p st between u t and u s and a path Q st of length q st between v t and v s such that p st and q st are odd and p st + q st ≡ 0 mod 4. Hence, they raised the following question. Problem 1.2. Except from bad Thetas and Theta trees, is there any bipartite graph which is not 2-detectable?

We partially answer to this question by showing (Theorem 3.1) that every bipartite graph with minimum degree at least 3 is 2-detectable. In particular, every cubic bipartite graph is 2-detectable.

We then restrict our attention to cubic graphs. For such graphs, by Brooks' theorem, if G = K 4 , then χ(G) ≤ 3, and hence by the result of Addario-Berry, Aldred, Dalal and Reed that det(G) ≤ 3. In [START_REF] Escuadro | A three-color problem in graph theory[END_REF], Escuadro, Okamoto and Zhang observed for some cubic graphs that: det(K 4 ) = 3; det(K 3,3 ) = 2, where K r,s is the complete bipartite graph with partite sizes r and s; det

(C 3 K 2 ) = 3, det(C 4 K 2 ) = 2, det(C 5 K 2 ) = 3 and if n ≥ 6 is an integer, then det(C n K 2 ) = 2,
where denotes the Cartesian product, and C n denotes the cycle of length n. We then exhibit some infinite families of cubic graphs with detection number 3. This allow us to characterize all cubic graphs up to ten vertices according to their detection number.

NP-completeness for cubic graphs

The aim of this section is to prove the following theorem. The proof of this theorem is a reduction from NOT-ALL-EQUAL 3SAT, which is defined as follows: Input: A set of clauses each having three literals. Question: Does there exists a suitable truth assignment, that is such that each clause has at least one true and at least one false literal?

This problem was shown NP-complete by Schaefer [START_REF] Schaefer | The complexity of satisfiability problems[END_REF].

In order to construct gadgets and proceed with the reduction, we need some preliminaries.

The halter is the graph depicted Figure 1. The vertices a and b are the ends of the halter, and the edges aa ′ and bb ′ are its reins. Proof. Let c be a detectable 2-coloring of G. Without loss of generality assume that c(uv

) = 1. If c(ua ′ ) = c(va ′ ) = c(ub ′ ) = c(vb ′ ), then code(u) = code(v), a contradiction.
Out of the four edges ua ′ , va ′ , ub ′ and vb ′ , assume that exactly three are of same color. By symmetry, assume that c(va

′ ) = c(ub ′ ) = c(vb ′ ). Suppose c(ua ′ ) = 1. If c(aa ′ ) = 1, then code(a ′ ) = code(u), a contradiction; if c(aa ′ ) = 2, then code(a ′ ) = code(v), a contradiction. Hence, c(ua ′ ) = 2. If c(bb ′ ) = 1, then code(b ′ ) = code(v), a contradiction; if c(bb ′ ) = 2, then code(b ′ ) = code(u), a contradiction.
Consequently, among the four edges ua ′ , va ′ , ub ′ , and vb ′ , two are of color 1 and the remaining two are of color 2.

If c(ua

′ ) = c(ub ′ ) and c(va ′ ) = c(vb ′ ), then code(u) = code(v), a contradiction.
By symmetry, assume that c(ua

′ ) = c(ub ′ ). So c(va ′ ) = c(vb ′ ) and c(ua ′ ) = c(va ′ ). Assume without loss of generality that c(ua ′ ) = 1. Since code(v) = (1, 2), c(aa ′ ) = c(bb ′ ) = 1. Consequently, we have c(aa ′ ) = c(ua ′ ) = c(uv) = c(ub ′ ) = c(bb ′ ) = 1 and c(va ′ ) = c(vb ′ ) = 2. See Figure 2 (a). Similarly, if c(uv) = 2, then we have c(aa ′ ) = c(ua ′ ) = c(uv) = c(ub ′ ) = c(bb ′ ) = 2 and c(va ′ ) = c(vb ′ ) = 1. See Figure 2 (b). INRIA Lemma 2.3.
Let G be a cubic graph. If a vertex x is the end of two halters in G, then in any detectable 2-coloring of G, x has code (3, 0) or (0, 3).

Proof. Assume for a contradiction, that the code of x is neither (3, 0) nor (0, 3). By symmetry, we may assume that x has code (2, 1). Therefore x is incident to two edges colored 1 and thus at least one of the rein e incident to it is colored 1. Therefore by Lemma 2.2, the neighbor of x through e has code (2, 1), a contradiction.

Proof of Theorem 2.1. Let C be a collection of clauses of size three over a set U of variables. We construct a cubic graph G = G(C ,U) as follows.

For every clause C ∈ C , we create a vertex v(C). For every variable u ∈ U, let C u be the set of clauses in which one of the two literals u and ū appears. We Clearly, the resulting graph G is cubic. Let us now prove that G is 2-detectable if and only if C admits a suitable assignment.

Suppose first that G admits a detectable 2-coloring. Let us establish few claims. The first one follows directly from Lemmas 2.2 and 2.3.

Claim 1.

In the variable gadget of every variable u, all the p(u,C) have the same code, which is either (3, 0) or (0, 3). Proof. By Lemma 2.3, p(u,C) has code (3, 0) or (0, 3). Without loss of generality, we may assume that code(p(u,C)) = (3, 0). Then by successive applications of Lemmas 2.2 and 2.3, all the p i (u,C), 0 ≤ i ≤ 5, and all the p ′ i (u,C), 1 ≤ i ≤ 5, have code (3, 0). Hence the edges p 4 (u,C)q(u,C) and p 5 (u,C)q(u,C) are both colored 1. Hence q(u,C) has code (2, 1) and so the reins of the halter with ends q(u,C) and p( ū,C) are colored 2, by Lemma 2.2. Similarly, one shows that the reins of the halter with ends q ′ (u,C) and p( ū,C) are colored 2. Hence by Lemma 2.3, the code of p( ū,C) is (0, 3).

Claim 3. For every clause C, the three neighbours of v(C) do not have the same code.

Proof. By construction, the neighbours of v(C) are all ends of two halters, and so have code in {(3, 0), (0, 3} by Lemma 2.3. Assume for a contradiction, that they all have the same code, say (3, 0), then the three edges incident to v(C) are colored 1 and the code of v(C) is also (3, 0), a contradiction. With these claims in hand, we can now prove that C admits a suitable assignment. Let φ be the truth assignment defined by φ(u) = true if all the p(u,C) of its variable gadget have code (3, 0), and φ(u) = f alse if all the p(u,C) of its variable gadget have code (0, 3). This assignment is well-defined because of Claim 1. Now by Claim 2, for any negated literal ū in a clause C, p( ū,C) has code (3, 0) if φ( ū) = true and p( ū,C) has code (0, 3) if φ( ū) = f alse. Now, by Claim 3 the three neighbours of v(C), which corresponds to the three literals of C do not have the same code. This implies that the corresponding literals do not have the same value. Therefore the truth assignment φ is suitable.

Conversely, suppose that C admits a suitable truth assignment φ. For each variable u, color the edges incident to each p(u,C) with 1 if φ(u) = true and with 2 if φ(u) = f alse. Similarly, color the edges incident to p( ū,C) with 1 if φ( ū) = true and with 2 if φ( ū) = f alse. It can easily be seen that such a coloring extends using the colorings of halter shown in Figure 2 to variable and negation gadgets, so that no two adjacent vertices in these gadget have the same code. It remains to show that every vertex v(C) has a code distinct from its neighbours. But since φ was p(u,C) suitable, at least one literal is false so the edge between v(C) and the vertex corresponding to this literal is colored 2, and at least one literal is true so the edge between v(C) and the vertex corresponding to this literal is colored 1. This implies that the code of v(C) is in {(2, 1), (1, 2)}. But in our coloring the code of the neighbours of v(C) are either (3, 0) or (0, 3). Hence we have a detectable 2-coloring.

p 1 (u,C) p 2 (u,C) p 3 (u,C) v(C) p 4 (u,C) p( ū,C) q(u,C) q ′ (u,C) p 0 (u,C) p ′ 2 (u,C) p ′ 1 (u,C) p ′ 3 (u,C) p ′ 5 (u,C) p 5 (u,C) p ′ 4 (u,C)

Bipartite graphs

In this section, our aim is to prove the following theorem.

Theorem 3.1. Every bipartite graph with minimum degree at least 3 is 2-detectable.

If any one of the parts of the bipartite graph have even number of vertices, Theorem 3.1 is an immediate consequence of Theorem 3.3 of [START_REF] Khatirinejad | Vertex-colouring edge-weightings with two edge weights[END_REF]. For sake of completeness, we give its proof here.

Theorem 3.2 (Khatirinejad et al. [8]). If G = ((A, B), E) is a connected bipartite graph with |B| even, then G admits an edge labelling c : E(G) → {1, 2} such that every vertex in A is incident to an even number of edges labelled 1 and every vertex in B is incident to an odd number of edges labelled 1. In particular, det(

G) = 2. Proof. Set B = {b 1 , b 2 , . . . , b 2p }. For every i, 1 ≤ i ≤ p, let P i be a path joining b 2i-1 to b 2i .
We start with all edges labelled 2. Then, for each i, 1 ≤ i ≤ p, one after antoher, we exchange the labels along P i . Hence at the end of this process, every vertex of A is incident to an even number of edges labelled 1 and every vertex in B is incident to an odd number of edges labelled 1.

To complete the proof of Theorem 3.1, we need some preliminaries. Let G be a graph. The closed neighborhood of vertex v is the set

N[v] := N(v) ∪ {v}.

Lemma 3.3. Let G = ((A, B), E) be a connected bipartite graph. Then there exists a nonempty set of twins S such that G -N[S] is connected.

Proof. If for some vertex v of G, G -N[v] is connected, then we have the result with S = {v}. So assume that for every vertex v of G, G -N[v] is not connected. Let us choose a vertex v of G such that G -N[v] has a component of largest possible size. Moreover, we choose v with largest possible degree among such vertices. Without loss of generality assume that v ∈ A.

By assumption, G-N[v] is not connected and let C be the vertex set of a component of In both cases, H is bipartite and V (H) ∩ B is B \ X so has even size. Therefore, by Theorem 3.2, H admits an edge labelling c : E(H) → {1, 2} such that every vertex in A \ S is incident to an even number of edges labelled 1 and every vertex in B \ X is incident to an odd number of edges labelled 1. Observe moreover that, when k is even, the edge ut is necessarily labelled 1.

G-N[v] of largest size. Then every vertex u in N(v) is adjacent to a vertex in C, for otherwise N[C] is included in a component of G -N[u], which contradicts our choice of v since |E(G[N[C]])| > |E(G[C])| as G is connected. Hence N[C] = C ∪ N(v). Set S = (V (G) \ N[C]) ∩ A. Let w be a vertex in S \ {v}. Then N(v) ⊆ N(w), for otherwise N(v) \ N(w)
Pick a vertex v ∈ S and extend c by labelling 1 to all the edges from v to X and all remaining edges incident to a vertex in N[S] with 2. Then, for every vertex b in

B \ N[S], code c (b) = (α, d G (b) -α), α ≡ 1 mod 2, for every vertex a in A\ {v}, code c (a) = (β, d G (a)-β), β ≡ 0 mod 2, for every vertex x in N(S), code c (x) = (1, d G (x)-1), and code c (v) equals (k, 0) if k is odd and equals (k -1, 1) if k is even. Hence c is a detectable 2-coloring because k ≥ 3.

Cubic graphs with detection number 3

In this section, our aim is to exhibit some infinite families of cubic graphs with detection number 3.

First, suppose there is a K 3 , say I, in G. Let V (I) = {x, y, z}. So E(I) = {xy, yz, zx}. Let x ′ , y ′ , z ′ be, respectively, the neighbors of x, y, z not belonging to I in G. Assume that x ′ = y ′ , y ′ = z ′ and z ′ = x ′ . Define subgraph I ′ by V (I ′ ) = V (I) ∪ {x, y, z} and E(I ′ ) = E(I) ∪ {xx ′ , yy ′ , zz ′ }. See Figure 4.

Lemma 4.1. If I ′ is a subgraph of a cubic graph G and if G has a detectable 2-coloring, then the edges of I receive both the colors and {code

(x), code(y), code(z)} is either {(3, 0), (2, 1), (1, 2)} or {(0, 3), (2, 1), (1, 2)}. Proof. Let c be a detectable 2-coloring of G. Suppose c(xy) = c(yz) = c(zx) = 1.
(The other possibility is similar.) Out of the three edges xx ′ , yy ′ , zz ′ , at least two are of same color. Without loss of generality assume that c(xx ′ ) = c(yy ′ ). Then code(x) = code(y), a contradiction.

Among the three edges xy, yz and zx, if color 1 appears twice, then {code(x), code(y), code(z)} = {(3, 0), (2, 1), (1, 2)}; if color 2 appears twice, then {code(x), code(y), code(z)} = {(0, 3), (2, 1), (1, 2)}; see Figure 5 Let M be a subgraph of a cubic graph G with edges Without loss of generality assume, by Lemma 2.2, that code(v 3 ) = code(v 6 ) = (2, 1). Then c(v Proof. Suppose G has a detectable 2-coloring c.

v 1 v 2 , v 2 v 3 , v 3 v 4 , v 4 v 5 , v 5 v 6 , v 6 v 7 , v 7 v 8 , v 2 v 4 , v 5 v
(2, 1) (1, 2) (b) (0, 3) (1, 2) (2, 1) Figure 5: Possible colorings of I ′ v 8 v 1 v 2 v 3 v 4 v 5 v 6 v 7 Figure 6: The graph M v 8 v 1 v 2 v 3 v 4 v 5 v 6 v 7 v 8 v 1 v 2 v 3 v 4 v 5 v 6 v 7 v 8 v 1 v 2 v 3 v 4 v 5 v 6 v 7 v 8 v 1 v 2 v 3 v 4 v 5 v 6 v 7 Figure 7: Possible colorings of M RR n°8115 v 1 v 2 v 3 v 4 v 5 v 6 v 7 v 8 v 9 v 10 v 11 v 12 v 13 v 14 N 3 v 1 v 2 v 3 v 4 v 5 v 6 v 7 v 8 v 9 v 10 v 11 v 12 v 1 v 2 v 3 v 4 v 5 v 6 v 7 v 8 v 9 v 10 N 1 N 2 Figure 8: The graphs N 1 , N 2 and N 3 v 2 v 3 v 4 v 5 v 6 v 7 v 1 v 8
2 v 3 ) = c(v 6 v 7 ) = 1. If c(v 2 v 7 ) = 1, then c(v 1 v 2 ) = c(v 7 v 8 ) = 1, and hence code(v 2 ) = code(v 7 ), a contradiction. If c(v 2 v 7 ) = 2, then c(v 1 v 2 ) = c(v 7 v 8 ) = 2,
If c(w 1 z 1 ) = c(w 1 z 2 ) = c(w 2 z 3 ) = c(w 2 z 4 ), then code(w 1 ) = code(w 2 ), a contradiction. If c(w 1 z 1 ) = c(w 1 z 2 ), then code(w 1 ) / ∈ {(1, 2), (2, 1)
). Similarly, if c(w 2 z 3 ) = c(w 2 z 4 ), then code(w 2 ) / ∈ {(1, 2), (2, 1)
}, again a contradiction to c(w 1 w 2 ). Hence, c(w 1 z 1 ) = c(w 1 z 2 ), say, 1 and c(w 2 z 3 ) = c(w 2 z 4 ) = 2. Note that code(w 1 ) = (2, 1) and code(w 2 ) = (1, 2), again a contradiction to c(w 1 w 2 ).

Hence, det(G) = 2. Now we construct a family of cubic graphs L n , n ≥ 2, with det(L n ) = 3 as follows: Begin with C 5n , the cycle of length 5n, say, v 0 v 1 v 2 . . . v 5n-1 v 0 ; add chords of distance 2, v 5r+1 v 5r+3 and v 5r+2 v 5r+4 for r ∈ {0, 1, 2, . . . , n -1}. If n is even, pair the vertices in {v 0 , v 5 , v 10 , . . . , v 5n-5 } in any order and join these pairs as edges; if n is odd, except three vertices in {v 0 , v 5 , v 10 , . . . , v 5n-5 }, pair the remaining vertices in any order and join these pairs as edges and add a new vertex v and join v to the omitted three vertices.

By Theorem 4.3 with i = 4 and Theorem 4.4, for n ≥ 4, det(L n ) = 3. We have to consider the cases n = 3 and n = 2. For n = 3, suppose L 3 has a detectable 2-coloring c. Consider the claw with center v and ends v 0 , v 5 , v 10 .

For ℓ ∈ {0, 1, 2}, if c(v 5ℓ v 5ℓ-1 ) = c(v 5ℓ v 5ℓ+1 ), then code(v 5ℓ ) / ∈ {(1, 2), (2, 1)}, a contradiction to c(v 5ℓ v). Hence, for every ℓ ∈ {0, 1, 2}, c(v 5ℓ v 5ℓ-1 ) = c(v 5ℓ v 5ℓ+1 ). There are 5 nonisomorphic cubic graphs on 8 vertices [START_REF] Meringer | Fast generation of regular graphs and construction of cages[END_REF]. It is known that det(C 4 K 2 ) = 2, see [START_REF] Escuadro | A three-color problem in graph theory[END_REF]. In the remaining four graphs, exactly two have detection number 3, and they are shown in Figure 10. Similar to the proof of Theorem 4.3 with i = 1, the graph in Figure 10 (a) has detection number 3. By ad hoc arguments one can check that the graph in Figure 10 (b) has detection number 3.

There are 19 nonisomorphic cubic graphs on 10 vertices [START_REF] Meringer | Fast generation of regular graphs and construction of cages[END_REF]. Out of these, exactly 6 have detection number 3. It is known that det(C 5 K 2 ) = 3, see [START_REF] Escuadro | A three-color problem in graph theory[END_REF]

Theorem 2 . 1 .

 21 The following problem is NP-complete. Input: A cubic graph G. Question: Is G 2-detectable?

Figure 1 : 2 . 2 .Figure 2 :

 1222 Figure 1: The halter

  construct a variable gadget associated to u, by considering a cycle on the |C u | vertices {p(u,C) | C ∈ C u } and replacing each edge ab of this cycle by a halter with ends a and b. Now for each variable u and clause C ∈ C u , we connect v(C) and p(u,C) with an edge if the literal u appears in C, and with the negation gadget depicted Figure 3 if the literal ū appears in C.

Claim 2 .

 2 In every negation gadget for (u,C), we have {code(p(u,C)), code(p( ū,C))} = {(3, 0), (0, 3)}.

Figure 3 :

 3 Figure 3: The negation gadget

  For a set S of vertices, we set N(S) := s∈S N(s) and N[S] := s∈S N[s]; and G[S] denotes the subgraph induced by S. Two vertices x and y are twins if N(x) = N(y).

Figure 4 :

 4 Figure 4: The graph I ′

7 and v 3 v 6 .Figure 6 . 4 . 2 .

 6642 SeeLemma If M is a subgraph of a cubic graph G and if G has a detectable 2-coloring, then the edges of M receive colors shown in any one of the Figure7.

INRIA

Figure 9 :

 9 Figure 9: The graph N 4

Theorem 4 . 5 .

 45 Let ℓ ∈ {0, 1, 2}. If c(v 5ℓ v 5ℓ-1 ) = c(v 5ℓ v 5ℓ+1 ) = 1, then code(v 5ℓ ) = (2, 1) implies that c(v 5ℓ v) = 1; if c(v 5ℓ v 5ℓ-1 ) = c(v 5ℓ v 5ℓ+1 ) = 2, then code(v 5ℓ ) = (1, 2) implies that c(v 5ℓ v) = 2. Since c(v 0 v 1 ) = c(v 4 v 5 ), c(v 5 v 6 ) = c(v 9 v 10 ), c(v 10 v 11 ) = c(v 14 v 0 ), we have c(v 0 v 1 ) = c(v 4 v 5 ) = c(v 5 v 6 ) = c(v 9 v 10 ) = c(v 10 v 11 ) = c(v 14 v 0 ) = c(v 0 v) = c(v 5 v) = c(v 10 v). Consequently, code(v 0 ) = code(v 5 ) = code(v 10 ) = code(v), a contradiction.Hence, det(L 3 ) = 3. Similarly, one can verify that det(L 2 ) = 3. Thus we have For each n, there exists a cubic graph of order 5n satisfying det(G) = 3.

Figure 10 :

 10 Figure 10: Cubic graphs on eight vertices with detection number 3.

  and hence code(v 2 ) = code(v 7 ), again a contradiction. Let w 1 w 2 be an edge of a connected cubic graph G. Suppose G -{w 1 , w 2 } contains four disjoint subgraphs J 1 , J 2 , J 3 , J 4 , where J i ∈ {K 4e, M -{v 1 , v 8 }} for i ∈ {1, 2, 3, 4}, and if w 1 is adjacent to a 2 degree vertex z 1 of J 1 and a 2 degree vertex z 2 of J 2 , and w 2 is adjacent to a 2 degree vertex z 3 of J 3 and a 2 degree vertex z 4 of J 4 , in G, then det(G) = 3.

	Theorem 4.4.

  } by Lemmas 2.2 and 4.2, we have a contradiction to c(w 1 w 2

RR n°8115

Proof. Let c be a detectable 2-coloring of G. Without loss of generality assume that c(v 4 v 5 ) = 1. If c(v 2 v 4 ) = c(v 3 v 4 ) and c(v 5 v 6 ) = c(v 5 v 7 ), then code(v 4 ) = code(v 5 ), a contradiction. Hence, either c(v 2 v 4 ) = c(v 3 v 4 ) or c(v 5 v 6 ) = c(v 5 v 7 ). Assume by symmetry that c(v 2 v 4 ) = c(v 3 v 4 ).

Suppose c(v

Hence the claim is true. Then c(v

In conclusion, we have only two possibilities for c(v

(Observe that code(v 2 ) = (2, 1) and code(v 7 ) = (1, 2).) See Figure 7 (b).

If c(v 4 v 5 ) = 2, then we have:

2) and code(v 7 ) = (2, 1).) See Figure 7 (d).

In all the four possibilities: c(v