Unsaturated subsurface flow with surface water and nonlinear in- and outflow conditions

0.98 MB·English

by  Heiko Berninger

Unsaturated subsurface flow with surface water
and nonlinear in- and outflow conditions∗
Heiko Berninger†, Mario Ohlberger‡, Oliver Sander§, Kathrin Smetana¶
January 14, 2013
3
1
0
2
Abstract

n
We analytically and numerically analyze groundwater flow in a ho-
a
J mogeneous soil described by the Richards equation, coupled to surface
water represented by a set of ordinary differential equations (ODE’s) on
1
1 parts of the domain boundary, and with nonlinear outflow conditions of
Signorini’s type. The coupling of the partial differential equation (PDE)

] and the ODE’s is given by nonlinear Robin boundary conditions. This
A
article provides two major new contributions regarding these infiltration
N conditions. First, an existence result for the continuous coupled problem
isestablishedwiththehelpofaregularizationtechnique. Second,wean-
.
h alyze and validate a solver-friendly discretization of the coupled problem
t
a basedonanimplicit–explicittimediscretizationandonfiniteelementsin
m space. The discretized PDE leads to convex spatial minimization prob-
lems which can be solved efficiently by monotone multigrid. Numerical
[
experiments are provided using the Dune numerics framework.

1
v
8 Keywords: saturated-unsaturated porous media flow, Kirchhoff transforma-
8 tion,convexminimization,finiteelements,monotonemultigrid,nonlineartrans-
4 mission problem
2
.
1 AMS Subject Classification: 35K61, 65N30, 65N55, 76S05
0
3
1
:
v
i
X
r
a
∗This work was supported by the BMBF under contract numbers 03MOPAF1 and
03KOPAF4
†SectiondeMathématiques,UniversitédeGenève,2-4rueduLièvre,CP64,1211Genève
4,Switzerland,[email protected]
‡InstituteforComputationalandAppliedMathematics,UniversityofMuenster,Einstein-
str. 62,48149Muenster,Germany,[email protected]
§InstitutfürGeometrieundPraktischeMathematik,RWTHAachenUniversity,Templer-
graben55,52062Aachen,Germany,[email protected]
¶InstituteforComputationalandAppliedMathematics,UniversityofMuenster,Einstein-
str. 62,48149Muenster,Germany,[email protected]
1
Figure 1: Typical shapes of the coefficient functions k(s) and p (s)
c
1 Introduction
This article is concerned with existence of solutions, and efficient numerical
approximation of unsaturated groundwater flow with surface water and nonlin-
ear in- and outflow conditions. To start with, let us introduce the considered
model. Let Ω⊂Rd, d=1,2,3 denote a bounded domain occupied by a homo-
geneous soil with boundary ∂Ω of class C1 and outer normal ν. Let Σ ,Σ
in out
and Σ ⊂∂Ω denote three relatively open, pairwise disjoint d−1-dimensional
N
C1-manifolds, representing the infiltration, outflow and Neumann boundaries,
such that ∂Ω=Σ ∪Σ ∪Σ and Σ ∩Σ =∅. The unsaturated ground-
in out N in out
water flow is modelled in the time interval [0,T] using Richards equation (see
e.g. [9]) for the water saturation s : Ω×[0,T] → [0,1] and the water pressure
p:Ω×[0,T]→R
n∂ s−div(cid:0)k(s)µ−1(∇p+e)(cid:1)=f in Ω×[0,T]. (1)
t
Here n denotes the porosity, µ the viscosity, k the permeability, and e the
gravity vector. For simplicity we set n=µ=1 and e=(0,0,1)T in this paper.
Due to our assumption of a homogeneous soil, the permeability k depends only
onthesaturations. IntheRichardsmodel,theairpressureintheporespaceis
assumed to be constant. We consider it normalized to p =0 and replace the
gas
waterpressurein(1)bythecapillarypressurefunctionp˜ (s). Thus,thecapillary
c
pressure p has to be seen as a multivalued graph with p (s) = {p˜ (s)} for
c c c
s<1 and p (1)=[p˜ (1),∞). Furthermore, we assume the residual saturation,
c c
meaning the saturation level below which no flow of water occurs, to be zero,
such that p is defined on [0,1]. This is necessary for the definition of the
c
infiltration boundary conditions. Additionally, we suppose p˜ (1) = 0, which is
c
essential for the limiting process on the outflow boundary. Typical shapes of
k(s)andp (s)aredepictedinFig.1. Itcanbeseenthatthecoefficientfunctions
c
degenerate in the sense that k(0) = 0 and that k and p may have slopes that
c
are unbounded in the neighbourhood of s = 1 and s = 0, respectively. These
degeneracies of the coefficient functions are one of the major challenges when
treating (1) analytically.
Weimposethreedifferenttypesofboundaryconditions. OnΣ weprescribe
N
homogeneousNeumannboundaryconditions. OntheoutflowboundaryΣ we
out
2
choose the Signorini-type boundary conditions
q·ν ≥0, p≤0, (q·ν)p=0 on Σ ×[0,T] (2)
out
with q = −k(s)(∇p+e) the Darcy velocity of the fluid. This models a hy-
drophilic porous medium that is in contact with air (cf. [31]). The inequalities
(2) express that (i) water may not enter the domain, (ii) the pressure cannot
be positive as p has been normalized to 0, and (iii) water can exit only if
air
p = 0. Concerning the analysis given below, a major problem is to give sense
to the product of traces in the last equality of (2). On the infiltration bound-
ary Σ , we only consider inflow of water into the soil and ponding of water on
in
the surface (without runoff). We assume that there is an external water source
r(x,t),x∈Σ . The water at x either infiltrates directly into the soil, or accu-
in
mulates in x (ponding). Hence the state of the surface water can be described
by the water table w :Σ ×[0,T]→R. At the absence of surface water runoff,
in
w is modelled as an ensemble of ordinary differential equations
∂ w =q·ν+r on Σ ×[0,T], (3)
t in
which is coupled to the subsurface boundary flux q (to be understood in a
suitable weak sense). Note that although Σ is the only part of the boundary
in
through which inflow can occur, we may also have outflow through Σ . A
in
possiblemodelforthefluxqonΣ istoconsideritasaverticalleakagethrough
in
averysmallsemiperviouslayerofthicknessbandhydraulicconductivityK [24],
h
which can be modelled (cf. [9, Chapt. 6.4]) as
p−w
q·ν = . (4)
c
Herec=b/K istheresistanceofthesemiperviouslayer. Itsinversec−1 isalso
h
known as leakage coefficient. Unfortunately, (4) is only valid for w > 0, as in
the case of w =0 and a negative pressure p, we would obtain an inflow into the
ground without water present on the infiltration boundary. Following [24] we
modify (4) by a w-dependent factor to obtain a new law valid for all w
p +p min(cid:8)1,(cid:0)w(cid:1) (cid:9)−w
q·ν = + − σ + , (5)
c
where p =max{p,0} and p =min{p,0}. Now, w =0 in (5) implies q·ν =0
+ −
ifp≤0. Theparameterσ isaregularizationthreshold,with(5)beingthesame
as (4) for w ≥ σ. As we will see later, solutions of this model fulfill w ≥ 0
automatically.
Existence results for unsaturated flow in porous media go back to the pio-
neering work of Alt, Luckhaus, and Visintin [2]. The basis of their approach -
which we will also use in this article - is the so called Kirchhoff transformation.
[37]picksupthemethodsusedin[2]andprovestheexistenceofsolutionsfor(1)
by regularization techniques, with a stronger solution concept than in [2]. The
vitalimprovementisthat[37]allowsdryregions,meaningregionswherethesat-
uration lies below the (positive) residual saturation, on the outflow boundary.
Hence,onehastomodifytheboundaryconditions,asthepressureisnotdefined
in dry regions. In [37] Robin boundary conditions are used in the regularized
problem and convergence to the limit on the outflow boundary is obtained via
3
defect measures. Without considering outflow conditions, [18] also studies dry
regions and proves the continuity of the free boundaries between the saturated
and unsaturated regions and between the unsaturated and dry regions. In [33]
uniqueness of the solution is proved for the problem introduced in [2] by adapt-
ing the methods developed in [32]. In order to show L1-contraction and with
this the uniqueness of the solution, doubling of variables techniques are used as
introduced in [29]. With respect to the condition on the infiltration boundary,
in [24] condition (5) is used to couple the Richards equation with a hyperbolic
PDE modelling ponding of water and surface runoff on the infiltration bound-
ary. Under the restriction ∂ s ∈ L1(Ω) existence and uniqueness is shown for
t
this coupled system. We remark that we do not need this assumption in our
approach.
Concerning the numerical treatment of the Richards equation, a rich litera-
ture can be found. For an overview we refer to [10, Sec.2.2] and [13] and the
literaturecitedtherein. TheRichardsequationhasbeendiscretizedusingfinite
volume methods [22, 27], mixed finite element methods [3, 35], finite elements
[25, 26], and discontinuous Galerkin schemes [5, 20]. In most cases, the re-
sulting algebraic systems were solved using Newton’s method, which, however,
suffers from ill-conditioning problems due to the degeneracies in the parameter
functions (Fig. 1). In contrast, Berninger et al. [13] discretize the Kirchhoff-
transformed Richards equation in a way that the resulting spatial problems
can be solved efficiently using a monotone multigrid method. Since this ap-
proach is based on convexity rather than on smoothness, it is also well-suited
for non-smooth Signorini-type boundary conditions. Numerical studies in [13]
demonstrate the efficiency of the solver, as well as robustness for extreme soil
parameters which correspond to parameter functions that degenerate into step
functions. An extension of this approach to heterogeneous soil using domain
decomposition techniques can be found in [12]. Concerning the numerical cou-
pling of Richards equation with surface water, we mention [20, 21, 39], where
shallowwaterequationsareusedasthesurfacewatermodel. Althoughdifferent
discretizations are applied, all these approaches enforce continuity of pressure
and normal flux, i.e., mass conservation, across the interface. These coupling
conditions have also been considered in an approach based on the Kirchhoff
transformed subsurface flow in [8] and [14], where an implicit time discretiza-
tion of the coupled subsurface and surface water models has been solved with
a heterogeneous domain decomposition method. In [14], however, the pressure
continuitywasreplacedbyaleakageconditionsuchas(4)thatrepresentedclog-
ging of a nearly impermeable river bed. Accordingly, a Robin–Neumann type
iteration was applied in [14], whereas in [8] a Dirichlet–Neumann type itera-
tion was used to solve the coupled system. In the present article, we choose an
explicit time discretization for the surface water as in [10, Sec.4.3].
Let us now give an outline of the article. In the first part of this article
we prove a new existence result for the Richards equation (1) endowed with
nonlinearoutflowconditionsofSignorini’stype(2)andcoupledtosurfacewater
by nonlinear Robin conditions (3, 5). The proof is done without requiring the
assumption∂ s∈L1(Ω). Toachievethisgoal,inSection2weintroduceaglobal
t
pressure u with the help of the Kirchhoff transformation and reformulate the
systemaccordingly. InSection3,westatethemainexistencetheoremandderive
its proof by combining two central approaches for existence and uniqueness
results for nonlinear, degenerated parabolic PDEs, namely the regularization
4
technique (cf. [31, 37]) and the method of L1-contraction (cf. [24, 33]).
Inthesecondpartofthearticle, startingwithSection4, wepresentourdis-
cretizationofthecouplednonlinearproblemthatisimplicit–explicitintimeand
uses finite elements in space. Via variational inequalities, the implicit–explicit
time discretization leads to spatial convex minimization problems which con-
tain nonlinear outflow and nonlinear physical Robin boundary conditions but
no Dirichlet conditions. We prove unique solvability of these problems if the
relative permeability is non-degenerate (e.g., regularized) and the water table
is non-zero. In Section 5, we give a numerical example combining the coupled
surface–subsurface problem and Signorini outflow conditions. We observe that
the proposed time discretization leads to a stable method for reasonable time
step sizes. We also observe that the surface water height remains nonnega-
tive as predicted by the theory, with the exception of certain singular points.
Our numerical approach proves to be useful even beyond the confines of the
assumptions needed for the existence result.
2 Global pressure formulation of the fully cou-
pled system and assumptions on the data
For a better assessment of the unsaturated flow equation, we use Kirchoff’s
transformation to obtain an equivalaent formulation for a global pressure u.
Wethensummarizetheresultingfullycoupledsystemwithnonlinearboundary
conditions and detail all assumptions on the data that we will need for our
analysis.
Let us define the Kirchoff transformation
p(cid:90)c(s) (cid:40)Φ˜(s) if s<1,
Φ˜(s):= k(p−1(q))dq, Φ(s)= (6)
c [Φ˜(1),∞) if s=1,
0
which yields ∇Φ(s)=k(s)∇p (s). We then define the global pressure function
c
(cid:40)
Φ˜(s(x,t)) if s(x,t)<1,
u(x,t):=
Φ˜(1)+k(1)(p(x,t)−p˜ (1)) if s(x,t)=1.
c
We thus have ∇u = k(s)∇p for a sufficiently smooth solution (s,p) of (1).
Therefore, we obtain as equivalent versions of the Richards equation
n∂ s=µ−1div(∇u+k(s)·γe)+f(s), u∈Φ(s),
t
(7)
or n∂ Φ−1(u)=µ−1div(∇u+k(Φ−1(u))·γe)+f(Φ−1(u)).
t
Conversely, the inverse Kirchhoff transformation κ−1 :u(cid:55)→p can be defined as
follows (cf. [13])
(cid:40) 0
−h(u) if u≤0, (cid:90) 1
p=κ−1(u):= , h(u) := dτ. (8)
u+ else k(Φ−1(τ))
k(1)
u−
Withtheglobalpressureformulation(7)wearenowpreparedtostateourfully
coupled system.
5
2.1 Fully coupled system for unsaturated flow in porous
mediawithinfiltrationandSignorini-typeoutflowcon-
ditions
Let us define Ω =Ω×[0,T], Ω =Ω×{0} and Σ =Σ ×[0,t] for t∈[0,T]
T 0 i,t i
with i = in,out,N. With given initial values s and w , we get the following
0 0
system of differential equations for the unknown global pressure u and surface
water height w
u∈Φ(s), and ∂ s=div(∇u+k(s)·e)+f(s) in Ω ,
t T
−(∇u+k(s)·e)·ν ≥0 on Σ ∪Σ ,
out,T N,T
u≤Φ˜(1) and k2(s)s−k2(1)1≤0 on Σ ,
out,T
0≤−((∇u+k(s)·e)·ν)·(k2(s)s−k2(1)1) on Σ , (9)
out,T
(∇u+k(s)·e)·ν =0 on Σ ,
N,T
(∇u+k(s)·e)·ν =g(u,w) on Σ ,
in,T
s=s in Ω ,
0 0
∂ w =(r−g(u,w)) in Σ ,
t in,T
w =w on Σ ,
0 in,0
where
u w h(u)ψ(w) (cid:26) (cid:16)w(cid:17) (cid:27)
g(u,w)=− + + + , and ψ(w)=min 1, .
k(1)c c c σ +
As the mappings s (cid:55)→ k2(s)s and s (cid:55)→ p (s) are both monotone functions, and
c
as p =0 and p˜ (1)=0, the outflow conditions stated in (2) and the outflow
gas c
conditions prescribed in (9) are formally equivalent.
It is also possible to allow p˜ (1)<0 as it is the case, e.g., for Brooks–Corey
c
parameterfunctions[16]. Butthentheoutflowconditionsin(2)andtheoutflow
conditions in (9) are no longer equivalent. We would have to deal either with
the conditions in (2) or with similar conditions formed by replacing p by u in
(2)asdone,e.g,in[10]. Insuchcaseitwouldbeverydifficulttogiveameaning
to the traces in the last equality in (2), as it demands an a priori estimate for
(cid:107)div(q)(cid:107) , which in turn necessitates an estimate for (cid:107)∂ s(cid:107) . In order
L2(Ω) t L2(Ω)
to solve this problem, one can either assume (cid:107)∂ s(cid:107) ≤C or Φ(cid:48)(s)≥c >0.
t L2(Ω) 1
If one makes one of these assumptions, it is also possible to show an existence
result for the system with original outflow conditions (2), based on the results
in[34]. Otherwise,itisonlyfeasibletoderiveaweightedL2-estimatefordiv(q)
[37] and the original outflow conditions have to be replaced by some surrogate
ones as suggested in [37] and (9).
2.2 Assumptions on the data
In order to guarantee among others the well-definedness of the boundary con-
ditions, the coefficient functions have to fulfill several assumptions. Additional
assumptions are needed for the proof of our existence result (see Theorem 3.1).
The relative permeability and capillary pressure functions are supposed to
fulfill the following assumption.
6
soil type k(1)γ/µ [cm/day] α [cm−1] l [−]
Hygiene sandstone 108.0 0.0079 10.4
Touchet Silt Loam G.E.3 303.0 0.005 7.09
Silt Loam G.E.3 4.96 0.00423 2.06
Guelph Loam (drying) 31.6 0.0115 2.03
Guelph Loam (wetting) - 0.02 2.76
Beit Netofa Clay 0.082 0.00152 1.17
Table 1: Typical values for k(1)γ/µ, α and l used for the van Genuchten
parametrization [40] for five different soil types.
Assumption 2.1. The permeability
k ∈C1([0,1],[0,∞)) is monotonically non-decreasing
and k(s) = 0 for s = 0. The capillary pressure p : (0,1] → {0,1}R is a
c
monotone graph given by a function
p˜ ∈C1((0,1),R), monotonically increasing.
c
We set p (s) = {p˜ (s)} for s ∈ (0,1) and p (1) = [p˜ (1),∞). Furthermore, we
c c c c
suppose p˜ (1) = 0 (hydrophilic case). Finally, we assume that with a constant
c
c >0 the following conditions hold:
0
∂ p (s)≥1/c ∀s∈(0,1),
s c 0
|∂ k(s)|2 ≤c k(s) ∀s∈(0,1),
s 0
(10)
(cid:112)
k(s)|p (s)|+ k(s)∂ p (s)≤c ∀s∈(0,1/2),
c s c 0
(cid:112)
(1−s) ∂ p (s)≤c ∀s∈(1/2,1).
s c 0
In order to guarantee that the condition on the infiltration boundary is well-
defined, we additionally presume that the coefficient functions are chosen in
such a way that
u ∈L∞(Ω ) =⇒p ∈L∞(Ω ). (11)
− T − T
Taking Assumption 2.1 for granted, we can prove u ∈ L∞(Ω ). Together
− T
with (11) this yields the well-definedness of the conditions on the infiltration
boundary.
Tojustifytheaboveassumptionsfromaphysicalviewpoint,wecheckifthey
can be fulfilled for the van Genuchten model [40] of the coefficient functions
k(s)=s1/2·(cid:0)1−(1−s1/m)m(cid:1)2, p (s)=−1(s−1/m−1)1/l,
c α
with m=1−1/l. In Table 1 some typical values for k(1)γ/µ — the hydraulic
conductivity at full saturation — and two variable parameters α and l for five
different soil types are listed. The values are extracted from [40]. Starting with
theimplication(11),weseethatthisholdstrueforl<3. Concerningestimates
(10),weobservethattheyareallfulfilledexceptforthesecondone. Thisisdue
to the fact that the description of the relative permeability k(s) proposed by
7
van Genuchten is bounded for s=1, but on the other hand its derivative blows
up. However, since unbounded derivatives of k are hydrologically unrealistic,
we can modify k(s) for s near 1 so that all estimates (10) are satisfied. We
emphasizetheimportanceofthethirdinequalityin(10),whichguaranteesthat
k declines fast enough for small s to absorb the pressure p , which goes to −∞
c
in this case.
Anotherparametrizationofk(s)andp (s)isgivenbyBrooksandCorey[16],
c
see also [17] and [40]. One can check that all estimates (10) are satisfied for
this parametrization, however, with regard to Assumption 2.1, we always have
p˜ (1) < 0. Moreover, it is a characteristic property of the Brooks–Corey func-
c
tionsthatthecrucialimplication(11)isneverfulfilled. Onthecontrary,thein-
verseKirchhofftransformation(8)isalwaysill-posedaroundtheminimalglobal
pressure u = Φ(0) with κ−1(u ) = −∞, cf. [10, Sec.1.3]. With regard to the
c c
analysis in physical variables, we can consider regularizations of Brooks–Corey
functions which exhibit the non-degenerate case k(s)>c >0, for which impli-
0
cation (11) is satisfied, cf. [10, Sec.1.4.3] and Section 4. However, even though
thenumericsfortheRichardsequationischallengingforthedegenerateBrooks–
Corey functions, we use them in our numerical examples in Section 5 analogue
to [13].
Concerningthesourceterm, weassumethatf :Ω ×[0,1]→Risbounded,
T
continuously differentiable in all entries, and monotonically decreasing in s.
Moreover, f shall satisfy f(x,t,1) ≤ 0 and f(x,t,0) ≥ 0 for all (x,t) ∈ Ω .
T
Eventually,wesupposethats(cid:55)→f(x,t,s)isaffineon(0,λ)forsmallλ. Further-
more,weassumethatr :Σ →R+ isboundedandLipschitzcontinuous. The
in,T
initial conditions shall be given by functions s :Ω→[0,1] and w :Σ →R+,
0 0 in
with s ∈L∞(Ω) and w ∈L∞(Σ ) and
0 0 in
s(x,0)=s (x) a.e. in Ω and w(x,0)=w (x) a.e. on Σ .
0 0 in
We impose the following compatibility condition for s :
0
∃p :Ω→R, p ∈p (s ) a.e. on {k(s )>0},
0 0 c 0 0
p ∈H1,2(Ω)∩L∞(Ω), s | >0.
0 0 Σin
3 Existence result for the degenerate coupled sys-
tem
The goal of this section is to prove an existence result for system (9). Unfortu-
nately, (9) exhibits some unpleasant properties. First, the system is degenerate
in the sense that k(s) → 0 for s → 0 and that the derivative p(cid:48)(s) may be
c
unbounded in the neighbourhood of s = 0 and s = 1. As a result, the solution
oftheRichardsequationlacksregularity,makingthetreatmentoftheboundary
conditions challenging, as it is difficult to give sense to the appearing traces of
the functions. Moreover, we deal with a coupled system due to the boundary
conditionsattheinfiltrationboundary. Toachieveourexistenceresult,wethus
firstregularizethecoefficientfunctionsk(s)andp (s)withasmallregularization
c
parameter δ in order to obtain a non-degenerate system. Afterwards, we fix δ
anddecouplethesystem. Existenceofsolutionsforthedecoupledsub-problems
8
thenfollowsbystandardarguments. ByapplyingBanach’sfixedpointtheorem,
we get a unique solution (u ,w ) for the regularized coupled system system, for
δ δ
fixed δ. The crucial part when proving a priori estimates independent of δ is
the derivation of a maximum principle for u and w . This can be achieved
δ δ
by considering classical solutions of a regularized parabolic system. Finally, we
pass to the limit in δ and get a solution for our original problem (9). Our main
result summarizes as follows.
Theorem 3.1 (Existence of a solution pair (u,w) of system (9)). Let the as-
sumptions from Section 2 be fulfilled and let the pair (u ,w ) be a sequence
δ δ
of regularized solutions in L∞(Ω )∩L2((0,T),H1,2(Ω))×C0([0,T],L1(Σ ))∩
T in
L∞((0,T),L2(Σ )). Then there exists a subsequence (u ,w ) and a solution
in δ δ
pair (u,w) ∈L2((0,T),H1,2(Ω)) ×C0([0,T],L1(Σ )) such that
in
u (cid:42)u in L2((0,T),H1,2(Ω)),
δ
w →w in C0([0,T],L1(Σ )),
δ in
and (u,w) are solutions of (9) in the distributional sense. Additionally, w ≥0
holds a.e. on Σ for all t∈[0,T].
in
In the following subsections we will prove Theorem 3.1 as indicated above.
3.1 Regularization of the coupled problem
At first, we regularize the coefficient functions.
Assumption3.2(Regularizedcoefficientfunctions). Theregularizedcoefficient
functions should fulfill the following conditions:
k ∈C1([0,1],(0,∞)) and ρ ∈C0([0,1],R) piecewise C1
δ δ
are monotonically increasing. For δ → 0 we have k (cid:38) k uniformly on [0,1],
δ
ρ →p uniformly on compact subsets of (0,1) and k (0)=δ2. We suppose that
δ c δ
(cid:83) ρ ([0,1]) = R and, for simplicity, that p (1/2) = ρ (1/2) holds. Finally, k
δ δ c δ δ
and ρ shall fulfill the inequalities of Assumption 2.1.
δ
The Assumptions 2.1 and 3.2 allow the definition of a regularized global
pressure via the Kirchhoff transformation
ρ(cid:90)δ(s)
u (x,t)=Φ (s)= k (ρ−1(q))dq,
δ δ δ δ
0
and guarantee that Φ is bounded from below and that Φ → Φ uniformly on
δ δ
compact subsets of [0,1). We set u =Φ(0) as the minimal global pressure.
c
Example 3.3. Let Assumption 2.1 be fulfilled and suppose that c s2 ≤ k(s) ≤
1
c s2 on(0,1)forconstants0<c ≤c . Thenthefollowingregularizationfulfills
2 1 2
Assumption 3.2.

p (δ)+ s−δ ∀s∈[0,δ],
 c δ
k (s):=δ2+k(s) ∀s∈[0,1], ρ (s):= p (s) ∀s∈(δ,1−δ],
δ δ c
p (1−δ)+ s−(1−δ) ∀s∈(1−δ,1].
c δ2
9
Note that all inequalities (10), which have to be satisfied by the regularized
coefficient functions, except for inequality (10) , are needed for the a priori
3
estimates. If (10) is not fulfilled, we would still be able to prove convergence
3
of the approximate solution pair to a limit (u,w), but this limit would not
necessarily solve our problem.
Using the regularized coefficient functions, we get the regularized Richards
equation
∂ s =div(∇u +k (s )·e)+f , u =Φ (s ),
t δ δ δ δ δ δ δ δ
with the pressure p = ρ (s ) and f = f(x,t,s (x,t)). In a straightforward
δ δ δ δ δ
mannerweobtaintheconditionsontheNeumannandtheinfiltrationboundary:
(∇u +k (s )·e)·ν =g (u ,w ) on Σ
δ δ δ δ δ δ in,T
(∇u +k (s )·e)·ν =0 on Σ .
δ δ δ N,T
Inspired by [37], we impose the following nonlinear Robin condition on the
outflow boundary:
1
−(∇u +k (s )·e)·ν = k (s )(p ) . (12)
δ δ δ δ δ δ δ +
In this way, we keep the idea behind the original ouflow conditions, as water
can flow out only if p ≥ 0 and if p > 0 we have a very large outflow due to
the scaling 1/δ. On the other hand, boundary condition (12) is more handy
than the original one, as we do not have to deal with a variational inequality.
Rewriting (12) with q =−(∇u +k (s )·e) one gets
δ δ δ δ
(p ) −0
q ·ν =k (s ) δ + ,
δ δ δ δ
which, similar to our condition on the infiltration boundary, can be interpreted
asapressuredifferencebetweenthepressureinsidetheporousmediumandthe
pressure outside p =0. The initial conditions are replaced by
gas
(cid:40)
ρ−1(p (x)) if s >1/2,
s (x,0)=sδ(x):= δ 0 0
δ 0 s (x) if s ≤1/2,
0 0
for all x∈Ω, guaranteeing Φ (sδ)∈H1,2(Ω).
δ 0
We get the following regularized system for the unknowns u and w
δ δ
u =Φ (s ), ∂ s =div(∇u +k (s )·e)+f(s ) in Ω ,
δ δ δ t δ δ δ δ δ T
1
−(∇u +k (s )·e)·ν = k (s )(p ) , on Σ ,
δ δ δ δ δ δ δ + out,T
(∇u +k (s )·e)·ν =0 on Σ ,
δ δ δ N,T
(13)
(∇u +k (s )·e)·ν =g (u ,w ) on Σ ,
δ δ δ δ δ δ in,T
(cid:40)
ρ−1(p (x)) if s (x)> 1,
s (x,0)=sδ := δ 0 0 2 in Ω ,
δ 0 s (x) if s (x)≤ 1, 0
0 0 2
∂ w =(r−g (u ,w )) in Σ ,
t δ δ δ δ in,T
w =w on Σ ,
δ 0 in,0
10