Asymptotic Analysis of the Laminar Viscous Flow Over a Porous Bed.
We consider the laminar viscous channel flow over a porous  surface.  The size of the pores is much smaller than the size of  the channel, and it is important to determine the effective  boundary conditions at the porous surface.  We study the  corresponding boundary layers, and,  by a rigorous asymptotic  expansion, we obtain Saffman's modification of the interface condition  observed by Beavers and Joseph. The effective coefficient in the  law is determined through an auxiliary boundary-layer type  problem, whose computational and modeling aspects are discussed in  detail. Furthermore, the approximation errors for the velocity and  for the effective mass flow are given as powers of the  characteristic pore size $\ep$. Finally, we give the interface condition  linking the effective pressure fields in the porous medium and in  the channel, and we determine the jump of the effective pressures  explicitly.
Introduction
Finding eective boundary conditions at the surface which separates a channel
ow and a porous medium is a classical problem.
Supposing a laminar incompressible and viscous
ow, we nd out immediately
that the eective
ow in a porous solid is described by Darcy's law. In the free
uid we obviously keep the Navier-Stokes system. Hence we have two completely
dierent systems of partial dierential equations. First, Darcy's law combined
with the incompressibility gives a second order equation for the pressure and a
rst order system for the velocity. In the Navier-Stokes system, the orders of the
corresponding dierential operators are dierent, and it is not clear what kind
of conditions one should impose at the interface between the free
uid and the
porous part. Clearly, due to the incompressibility, the normal mass
ux should
be continuous. Other classically used conditions are the continuity of the pressure
and, for a free
uid, the vanishing of the tangential velocity at the interface.
Let us discuss the mathematical background of the interface conditions. It
is well-known that Darcy's law is a statistical result giving the average of the
momentum equation (the Navier-Stokes equations) over the pore structure. Its
rigorous derivation involves the weak convergence in
div) (respectively the
two-scale convergence) of velocities, and only the continuity of the normal velocities
is preserved. Other continuity conditions at the interface are generally lost, such
that further analysis is required.
Concerning other interface conditions used in engineering literature, the vanishing
of the tangential velocity is found to be an unsatisfactory approximation,
and in [2] a new condition is proposed. The condition reads
where ~u e is the eective velocity in the channel, ~v f is the mean ltration velocity
given by Darcy's law, ~ is a tangent vector to the interface, ~ is the normal into the
Asymptotic analysis of the laminar viscous
ow over a porous bed 3
uid, K is the permeability of the porous medium, and the scalar  is a function
of the geometry of the porous medium.
In [2], this law is derived by heuristic arguments and justied experimentally.
A theoretical attempt to derive (1.1) is undertaken in [15] and, using a statistical
approach, a Brinkman type approximation in the transition layer is derived. A
matching argument then allows to obtain the formula
The interested reader can also consult the lecture note [4].
Dierent considerations can be found in [5] and [11]. They distinguish two
cases:
(a) The pressure gradient on the side of the porous solid of the interface is normal
to the interface. Consequently, we have a balanced
ow on both sides of the
interface. Then, using an asymptotic point of view, the following laws are
obtained in [11]: ( ~u e
on the interface. This case describes the
ows in cavities. The mathematical
justication is in [9]. We shall not consider it in this paper.
(b) The pressure gradient on the side of the porous solid at the interface is
not normal. This case is considered in the fundamental paper [5]. After
discussing the orders of magnitude of the unknowns it is found out that
on the interface the velocity of the free
uid is zero, and the pressure is
continuous.
All results cited above are not mathematically rigorous. Furthermore, dierent
approaches give dierent results and two natural questions arise immediately:
(Q1) What are the correct matching conditions (i.e. conditions at the interface)
between those two
ow equations?
(Q2) What are the eective constants entering the matching conditions?
We are going to answer those questions in the following. In Section 2, we
dene our problem and discuss some simple approximations. In Section 3, we
introduce an important auxiliary problem of boundary layer type which we need
to construct a better approximation, and Section 4 gives additional results for an
analogous problem on a nite strip. Then, in Section 7, the eective equations
with the Beavers-Joseph type boundary condition are presented together with
improved error estimates. Finally, in Section 6, we show how the computation of
the constants involved in the interface conditions works in practice. Especially, it
turns out that the dierence between Darcy's pressure in the porous part and the
eective channel pressure is equal to a constant multiple of the eective channel
shear stress, thus contradicting [11].
Asymptotic analysis of the laminar viscous
ow over a porous bed 5
e
e
Y
Y *
Z *

Figure

1: The
ow
region
2 Setting of the problem
This section deals with the equations describing
uid motion over a porous bed,
under a uniform pressure gradient. We assume the condition of the experiment by
Beavers and Joseph [2], i.e. a stationary laminar incompressible viscous
ow.
For simplicity reasons, we consider a
ow over a periodic porous medium with
a characteristic pore size ". The
ow
region
of two parts, see Figure 1. The upper
part
and the lower
part
2 is the porous medium which is obtained by putting solid
obstacles of size "  b into the
domain
the (permeable) interface
betweenand
2 . More precisely, let
made of two complementary parts, the solid part Z  and the
uid part Y  . It is
assumed that Z  is a smooth closed subset of Y , strictly included in Y , Y  \Z
and Y  [Z Now assume that both b and L are integer multiples of ". Then
the
domain
can be covered by a regular mesh of N(") cells of
size ". Each cell Y "
divided into a
uid part "(Y  + k) and a solid
6 Willi Jager, Andro Mikelic and Nicolas Neu
k). The
uid
part
2 of the porous medium is therefore
"=
and the whole
ow region
is
After specifying the geometry, we consider the equations determining the velocity
eld and the pressure eld p " in the Beavers{Joseph experiment:
in
in
(@
For small ", this problem is extremely dicult to solve, so that one has to look
for approximations. Classically, the system (2.3)-(2.8) used to be approximated
by a Poiseuille
ow
in
for
and, indeed, it could be shown in [10] that (2.9) is an approximation in the following
sense:
Z
Z
Z
Asymptotic analysis of the laminar viscous
ow over a porous bed 7
The above estimates indicate that in the L 2 -sense
in
in
in
in
in
on  and the eective mass
ow behaves as H 3
As shown in [2] and [15], this O(") approximation often is not good enough.
Therefore, we would like to continue with the asymptotic expansion for ~u " and p " .
Energy estimates from [10] imply that in the interior
of
1 , but
globally in
Further, they show that there is an oscillatory
boundary layer conned to the neighborhood of  which can be represented in the
~  bl;"
bl ( x
where ~
are solutions to a Stokes problem on an innite strip Z bl , which we
discuss in the following Section 3.
Z -2
Z bl

Figure

2: The boundary layer cell Z bl
3 The auxiliary boundary layer problem
Figure 2. For later use, we also
introduce the following notation
Z <k := Z bl \ (0; 1)  (1;
l := Z bl \ (0; 1)  (l;
Now f ~
are the solutions to the problem
~
div y
~
[(r y
~
~
Asymptotic analysis of the laminar viscous
ow over a porous bed 9
Further, in order to dene ! bl uniquely, we require
A variational formulation of this problem is: nd a locally square integrable
vector eld ~  bl 2 W satisfying
Z
Z bl
r ~  bl
Z
where W is the function space which contains all y 1 -periodic, divergence-free, locally
square integrable vector elds ~z dened on Z bl , having nite viscous energy
(i.e.
R
satisfying the no-slip boundary condition at
the solid boundaries [ 1
In [9] (Proposition 3.22, pages 462-463) the following result is proved:
Theorem 1 The problem (3.11) has a unique solution ~
bl 2 W . It is locally
innitely dierentiable outside of the interface S. Furthermore, a pressure eld ! bl
exists, which is unique up to an additive constant and locally innitely dierentiable
outside of the interface S, such that (3.5) holds. In a neighborhood of S we have
~
bl
The following lemma states some simple properties of f ~  bl
l
Z
Z bl
jr ~
Proof. (3.12) follows immediately from (3.6) by
Z 1@ bl@y 2
Z 1@ bl@y 1
Then, by integrating the second component of the momentum equation (3.5) over
the rectangle (0; 1)  (a; b) we obtain
and (3.13) follows by applying (3.12).
Next, integrate the rst component of the momentum equation (3.5) over (0; 1)
with respect to y 1 . Using (3.13) we nd out that
linear
function of y 2 . Since it is bounded in the limit y it has to be a constant
for y 2  0, which proves (3.14).
Finally, after taking
bl as the test function in (3.11), we obtain (3.15).
We expect that the problem (3.5){(3.10) represents a boundary layer. This
means that changes of the velocity and the pressure elds are concentrated around
the interface S and vanish very rapidly with increasing distance from S. In linear
elasticity, results of this type are called Saint-Venant's principle. Saint-Venant's
principle is also valid in our case, and we want to discuss this in more detail.
We start our considerations with Z (the part of Z bl lled by the
uid). Here,
one can obtain sharp decay estimates by using results for the decay of solutions of
general elliptic equations (see Theorem 10.1 from [12], which is an application of
Tartar's lemma). However, in our particular situation we can give a direct proof:
Theorem 3 Let
Then, for every a > 0, we have
Asymptotic analysis of the laminar viscous
ow over a porous bed 11
Proof. Applying the curl -operator to equation (3.5), we see that
ag : (3.19)
Furthermore, by (3.14) and periodicity, we conclude that
Consequently, the solution  bl of (3.19) may be written as
We see, that the leading decay exponent is proportional to 2, such that (3.17)
and (3.18) follow.
With the help of this estimate for
bl , we can prove exponentially
fast stabilization of ~
Corollary 4 Let
Then for every a > 0 and every  < 2 we have
~
and   D  ~
Proof. We follow [9]: by (3.16) and (3.6) the functions  bl
bl
bl
@
Using the explicit form of the right hand side given by (3.20), a variation of constants
yields also an explicit representation for  bl
bl
(D 1
(D 1
bl
(D 2
1;n
(D 2
e 2ny2 (3.27)
with the additional relations
4n
4n
The representation (3.26),(3.27) allows us to conclude (3.22) and (3.23).
For the pressure eld ! bl , we have:
Corollary 5 For
we have   ! bl (y
Proof. Taking the divergence of the momentum equation (3.5), we nd that
bl is square integrable and the averages over
sections fy ag of ! bl C bl
are zero, ! bl can be written as
from which the assertion follows.
Now we turn our attention to the porous part Z . Due to the presence of the
solid obstacles, our estimates will be much less precise in this case.
Lemma 6 Let the distance between the solid obstacles and the boundary of the
unit cell Y be bigger than or equal to 2d  . Let ~
bl be the solution of (3.5){(3.10),
and let l 2 Z; l < 1. We introduce a function ~
lower by
~
lower
bl for y 2  l (3.32)
Asymptotic analysis of the laminar viscous
ow over a porous bed 13
~
lower
lower
bl
d
Z y1 bl
Z l
bl
Z 2l y 2
l
lower

bl
Z 2l y 2
l
bl
for
Then ~  lower 2 W \ H 1 (Z bl

r ~
lower
C lower

r ~  bl
where the constant C lower satises
C lower d 3
with C P denoting the Poincare constant appearing in

bl

r bl
If we further assume that the ball B  (1=2; 1=2) with radius  and center (1=2; 1=2)
is contained in Z  , an easy calculation in polar coordinates yields the estimate
which (together with d   1) results in the estimate

r ~  lower

r ~
bl
Proof. By a straightforward calculation, one veries that ~  lower
This calculation also yields the estimates
Z
dy d 2
bl   2
dy
Z
bl   2
dy
14 Willi Jager, Andro Mikelic and Nicolas Neu
Z
@ lower@y 2
dy d 4
Z
bl   2
dy
Z
bl   2
dy
+d 4
Z
@ lower@y 2
dy
Z l
l d
bl
d
Z
dy  4
Z
where Z l;l we use the simple trace estimates
Z
bl   2
dy
Z
dy (3.44)
Z l
l 1
bl
Z
bl   2
dy
Z
We insert (3.44),(3.45) and (3.38) into the sum of (3.41){(3.43) and use d   1 to
obtain (3.36) and (3.37).
With the help of this lemma, we are able to deduce exponential decay of r ~
Proposition 7 Assume that l 2 Z with l < 0. Dene Z <l as in (3.2). Then the
solution ~
bl of problem (3.11) fullls

r ~
bl


r ~
bl

e
lower
where
lower =2 ln

C lower
Proof. If we test (3.11) with some function ~  bl ~  lower , where ~
lower is the
function constructed in the previous lemma, this results in the estimate

r ~
bl
C lower

r ~
bl
Using the hole-lling technique as in[9], Lemma 2.4 (we add C lower

r ~
bl
to both sides of (3.48) and evaluate the resulting recursion), we obtain exponential
decay with rate
lower given by (3.47).
Applying local regularity results, one immediately obtains:
Asymptotic analysis of the laminar viscous
ow over a porous bed 15
Corollary 8 For any a < 0,  2 IN , the solution f ~
exponentially
for
C(a; )e
lower
and   D  ! bl (y
C(a; )e
lower
for all y 2 < a < 0.
Proof. From Proposition 3.7 of [9] we know that

Z
Z l

C

r ~
bl

Z
Z l
Z
Z l+1
C

r ~
bl

which implies exponentially fast stabilization to a constant (which has to be zero
because of (3.10)):


Ce
lower jlj for l <
Pointwise estimates for ~
bl , ! bl and their derivatives can then be obtained as usual
by dierentiating the equations obtaining estimates for higher derivatives, which
can then be used with the Sobolev embedding theorem.
The above results imply that ~
bl is a boundary layer type velocity eld and
bl is a boundary layer type pressure. Only the constants C bl
1 and C bl
! from (3.21)
and (3.29) will enter the eective
ow equations in the channel. They contain the
information about the geometry of the porous bed.
Remark 9 As we shall see from the numerical examples, in general C bl
However, if the geometry of Z is axisymmetric with respect to re
ections around
the axis y
This result is
obtained by the following simple argument:
Let Z be axisymmetric around the axis y
bl be a solution
for (3.5){(3.10). Then ( bl
a solution. By uniqueness, it must be equal to ( ~
we conclude that
Asymptotic analysis of the laminar viscous
ow over a porous bed 17
4 Approximation of the boundary layer problem
on a nite domain
In this section we propose a scheme for computing the actual values of C bl
1 and
cases where the geometry of the porous medium is known. Since problem
resp. (3.11) is dened on an innite domain Z bl , the rst step is
to approximate its solution with solutions of problems dened on nite domains
of the form Z k
l := Z bl \ (0; 1)  (l;
bl
k;l g be
the solution of
~
bl
l [ Z k
div y
~
bl
l [ Z k
bl
k;l
[fr y
~
bl
~
bl
bl
k;l g is y
with the following additional boundary conditions motivated by Corollary 4 and
Proposition 7:
~
bl
Further, to dene the pressure eld ! bl
k;l uniquely, we require
Z
Analogous to the estimates on ( ~
of the previous section, one can prove
the following properties for the solution ( ~  bl
k;l ) of (4.1){(4.8):
0 , the solution f ~
bl
k;l g of problem (4.1){(4.8) can be
represented as
(D 1


e 2nk
where
1;k;l :=
Z 1( bl
and
k;l (y) =X
with
!;k;l :=
Z 1! bl
Z 1! bl
This representation immediately yields, that for every 0 < a < k, 0 <  < 2, and
are constants C(a; ); C(a; ; ) such that for all y 2 Z k
a we have
~
bl
D  ~
bl
Asymptotic analysis of the laminar viscous
ow over a porous bed 19
and   ! bl
!;k;l
For

~
bl
k;l

l
Ce
lower jmj (4.18)
with
lower from (3.47). Furthermore, for y 2 Z 0
l
k;l (y)
Ce
lower
~  bl
k;l (y) (C bl
Ce
lower
and   ! bl
Ce
lower jy
From these estimates we obtain:
Proposition 11 Let f ~
bl
k;l g be the solution of problem (4.1){(4.8). Then, for
every  < 2, a constant C exists such that

r ~
bl
k;l r ~
bl

l
lower
Proof. Let  := ~
bl ~
bl
k;l . Then (; !) is y 1 -periodic, vanishes
on
l
k=1 (@Z  (0; k)) and solves
in Z k
l . Set
lower
Jager, Andro Mikelic and Nicolas Neu
where  lower is derived from  in the same way as ~
lower was derived from ~
bl in
Lemma 6. Testing (4.23) with ^
, we obtain
Z
Z l+1
Z 1rr lower =
Now note that ! may be replaced by ! c for an arbitrary c 2 IR due to
R
=const  Then we can apply the exponential stabilization results for
f ~
k;l g from (3.22), (3.23), (3.30), (4.15), (4.15), and (4.17) to
obtain (4.22).
The approximation error between C bl
1 and C bl
1;k;l can then be estimated as follows

Corollary 12 For every  < 2 there is a constant C such that
1;k;l
lower
Proof. Note that
1;k;l
Z
bl

bl

which can be estimated as desired by using Poincare's inequality on Z 0 [ Z 1
together with (4.22).
In order to obtain estimates for the pressure dierence ! bl
we need the
following result:
Lemma 13 For each F 2 L 2 (Z k
l
R
l
there is a function ~
l ) vanishing on the boundaries y
l
satisfying div ~
together with the stability estimate
l
l
More generally if F 2 H r (Z k
l ) for some integer r  0, then ~
' can be chosen such
that
l
Asymptotic analysis of the laminar viscous
ow over a porous bed 21
Proof. The proof is similar to the proof of Lemma 3.4 in [9]. We search for ~
in the form
with div correcting the non-zero boundary values of r.
More precisely, let  be the solution to
@
@
n=l
(@Z  (0; n)) [ fy
Since
R
l
the testing of (4.32) with  yields
Z
l
Z
l
Z
l
l
Z
l
l
l
and therefore the estimate
l
Note that we have used the Poincare inequality

l
Z
l

l
l
with a constant C being independent of k and l. This estimate can easily be proved
by extending  from Z k
l to the rectangle (0; 1)  (l; k) and by using the Poincare
estimate there. Furthermore, by localizing we also get estimates for higher order
derivatives of  in the form
l
l
The function # from (4.31) has to correct the non-zero boundary values of
and should therefore fulll
22 Willi Jager, Andro Mikelic and Nicolas Neu
@#
@
@
on
n=l (@Z  (0; n)) [ fy lg. Since @Z  is smooth, such a function
can be constructed by a local H r -lift such that
l
l
with C being independent of k and l.
The combination of (4.37) and (4.41) then yields the desired regularity estimate
for ~
', and the lemma is proved.
With the help of this lemma we obtain:
Proposition 14 Let ! bl and ! bl
k;l be the pressure elds determined by (3.5){(3.10),
resp. (4.1){(4.8). Then, for every  < 2, there is a constant C such that

k;l

l
lower
For the dierence C bl
!;k;l we have the better estimate
!;k;l
C
lower
Proof. Set
l
Z
l
and let F :=
k;l
!. Since
R
l
yields a function ~
which
we can use as test function in the dierence of momentum equations
bl ~
bl
Inserting
! in the pressure term, testing with ~
and doing a partial integration
we obtain Z
l
r( ~
bl ~
bl
k;l )r~'
Z
l
Asymptotic analysis of the laminar viscous
ow over a porous bed 23
If we now use the stability estimate (4.29) together with (4.22), this yields


l
lower
Finally, to estimate
!, let
Z
Then, obviously, we have
By Theorem 3.7 of [9], we also have

r( ~
bl ~
bl

which can be used together with (4.21) to estimate
lower

r( ~
bl ~  bl

l
Using the triangle inequality in (4.47) together with (4.51) and (4.22), we obtain
(4.42).
In order to get (4.43), we note that
!;k;l
such that the estimate follows directly from (4.51) and (4.22).
We complete this section with a regularity result, which we will need in the
following.
Proposition
l f; g be the y 1 -periodic solution of
in Z k
l with boundary conditions
l
and pressure normalization Z
Then we have the estimate
l
l


l


l
l
C
l
with a constant C independent of k and l.
Proof. We rst note that we can bound k 2 k L 2 (Z k
l ) and


l
in
terms of krk L 2 (Z k
l ) . In the lower region, this follows from Poincare's inequality
applied on every cell. In the upper region, we rst have
l

@

l
because
setting
we also obtain
l

@

Additionally, the Hardy inequality



@

together with an estimate of S 1 (0) by trace inequality and Poincare inequality on


l
Asymptotic analysis of the laminar viscous
ow over a porous bed 25
Next, for deriving a bound for krk L 2 (Z k
l ) , we test the momentum equation (4.53)
with  to obtain
l
Z
l
On Z 0
l , we may again use Poincare's inequality to obtain
Z
l
f
l
l
and because of
also the estimate
Z
can be obtained by applying Poincare's inequality. In order to get an estimate for
R
and write
Z
Z
Z
Here, the rst term can be estimated as
Z
while for the second term we have
Z


using again (4.62). Combining these estimates, we obtain the desired bound for
26 Willi Jager, Andro Mikelic and Nicolas Neu
Next, we localize the problem by multiplying with smooth cut-o functions
are identically 1 on Z i and vanish for y 2
Denote with S the region Z i+1
l . By shifting the pressure by a suitable constant
we may assume that
R
(otherwise, the boundary conditions on
@Z  would be violated), we may also shift  by constant multiples of e 1 such that
we may also assume that
R
The resulting function f ~
is a
solution of
~
div ~
in S with Dirichlet and periodic boundary conditions in the case i < k 1 and a
combination of Dirichlet, periodic and slip boundary conditions when
We examine the right hand side of (4.70). Since
Z
Z
Z
r ~
we have kr~k H 1 (S)  C

r ~

. Since
R
0, we can apply Proposition 1.2,
Ch. I of [16] to obtain
C

r ~

Then, however, the right hand side of (4.70) is in L 2 while the right hand side
of (4.71) is in H 1 , and we can apply Proposition 2.2 of Ch. I of [16] to get that


and kr~k L 2 (S) can be bounded in terms of krk L 2 (S) (with a constant
depending on the geometry of the porous inclusion). We note that the slip boundary
fy easily eliminated by making an even extension for  1 ;  and f 1
and an odd extension for  2 and f 2 (the zeroth order re
ection). Summing up
these local estimates we get that kD 2 k L 2 (Z k
l ) and krk L 2 (Z k
l ) are bounded by a
multiple of krk L 2 (Z k
l ) with a constant that does not depend on k and l. Thus,
Proposition 15 is proved.
Asymptotic analysis of the laminar viscous
ow over a porous bed 27
Corollary
bl
k;l g be the solution to (4.1){(4.8), and let
~
:= ~
bl
k;l ~
y 22
Then

l

r! bl
k;l

l
where C can be chosen independent from k and l.
28 Willi Jager, Andro Mikelic and Nicolas Neu
5 Discretization
Now, we turn our attention to the discretization of problem (4.1){(4.8). Essen-
tially, we use a stabilized nite element discretization in the sense of [8],[3]. Unfor-
tunately, in both of these papers only polygonal domains were considered, so that
the direct application of these results to our domain Z k
l is unsatisfying because
of the curved boundaries
l
We resolve this diculty by using
generalized domain partitions given by nonlinear mappings of the elements while
essentially keeping the approximation results known for usual domain partitions
with linear mappings, see [19], [20], [1], [14]. While this is very convenient theoret-
ically, the practical implementation will usually be too complicated. Therefore, it
is important that the use of a simpler polygonal approximation of the domain can
be interpreted as a perturbation of this approach, see the discussion at the end of
this section.
Let
be a Lipschitz domain, and let T h be a partition
of
in subsets e, called
elements, where each element e is the image of a reference element ^
e under a
mapping
e is either the reference triangle ^
or the reference quadrangle ^
We require the following properties of the partition
1.
2. For two elements
3. A side of an element e (which is dened as the image of a side of ^ e) is either
a subset of
@
or the side of exactly one other element e 0 6= e. In the second
case, we require that the mapping  1
restricted to the corresponding
side of ^
e is linear.
4. The
are bi-Lipschitzian mappings with
Asymptotic analysis of the laminar viscous
ow over a porous bed 29
e
e (y)j
for some constant C 1 > 0.
5.  e 2 W 2;1 , with kD 2  e k 1  C 2 diam(e).
6. For simplicity reasons, we consider only quasiuniform partitions, i.e. a constant
exists such that for h := max e2T h
diam(e) we have diam(e)
Arbitrary ne triangulations of this kind exist, which can be shown by modifying
triangulations of polygonal approximations of the
domain
see [19]. The quality
of the domain partition T h is determined by the constants C
renement of such a partition results in a new partition which still has the above
properties (especially property 5).
Next, let
where P 1 is the space of linear polynomials, and Q 1 is the span of P 1 and the
polynomial x 1 x 2 . With these nite elements, we can dene the space
The following approximation result is then the substitute for approximation
results on triangulations with linear element mappings.
Theorem 17
Let
h be as above. Then, for
, an u h 2 S h
exists such that
Second, for
, an u h 2 S h exists such that


Additionally, if u has zero boundary values on some component of
@
, then also u h
can be chosen such that it has zero boundary values on that boundary component.
In all cases, the constants only depend on the smoothness of the domain and
the quality of the domain partition T h .
Proof. See [19], [20] for the case of (deformed) triangular elements, and [14]
where also the case of (deformed) quadrilaterals is handled by using a generalized
interpolation operator of Clement type.
Remark In contrast to standard interpolation results for triangulations with
linear element mappings, the term kruk L
appears on the right hand side of (5.6).
This is due to the fact that linear functions are no more interpolated exactly (in
contrast to constant functions). Note that this estimate is only possible because the
e approximate linear mappings in the limit diam(e) ! 0 (see property 5 required
for a partition T h ).
let
l , and let T h be a domain partition of Z k
l (tting across the
lateral boundaries, where periodic boundary conditions are prescribed for problem
(3.5){(3.10)). Then dene S h analogously to (5.4) as
l
the ansatz space for the velocity eld as
~
l
and the ansatz space for the pressure eld as
Z
We now search for ( ~
bl
being the solution to
Z
l
r ~  bl
Z
l
r! bl
k;l;h ~
Z
~
Z
l
div ~
bl
Z
r! bl
Asymptotic analysis of the laminar viscous
ow over a porous bed 31
for all (~' k;l;h ;  k;l;h
. The second term in (5.11) must be included for the
above pair of ansatz spaces to stabilize the discretization, see [8],[3]. The constant
can be any positive number. The following error estimates then hold:
Proposition 19 Let T h ; ~
be dened as above. Assume that the interior
boundary aligned with the sides of the elements of T h , let ( ~
bl
be the solution of problem (4.1){(4.8), and let ( ~
bl
H h  L h be the
solution of (5.10){(5.11). Furthermore, let
bl

r ~
bl
k;l

l

bl
k;l


bl
k;l

l

r! bl
k;l

l
Then we have an error estimate in the viscous energy norm of the form

r( ~
bl
k;l;h
~
bl

l
ChR( ~  bl
together with a stability estimate for the pressure gradient of the form@ X

r(! bl
ChR( ~
bl
We further have the following L 2 -error estimate for the pressure

k;l;h (! bl
k;l  Z k
l
Z
k;l dy)

l
bl
and the velocity

~
bl
k;l;h
~  bl
k;l

l
bl
with a constant C which is independent of k and l.
Proof. The proof can be done as in [8], [3], or [7] with the following modi-
cations: rst, on our more general domain partitions, the approximation results
from theorem 17 replace the standard ones. Second, [8], [3], [7] handle only the
Dirichlet case. However, one can easily check, that their proofs can be transfered
with almost no changes because the slip boundary condition still allows for partial
integration with vanishing boundary terms. For (5.14), the proof uses the
assertion of Lemma 13 such that the factor k needs to be included in the
estimate. For (5.15), one needs the H 2 -regularity estimate from proposition 15
which introduces the factor k + 1.
We now describe the denition of the discrete approximations to C bl
1;k;l and
!;k;l . First, we set
Z
Second, in order to obtain good error estimates for the numerical approximation of
!;k;l , we dene C bl
!;k;l;h as a smoothly weighted average of the pressure eld in the
following be a Lipschitz-continuous function satisfying
~
elsewhere
It is obvious that Z
~
Z
Z l
such that we may approximate C bl
!;k;l by
Z
l
~
Then we can prove:
Proposition 20 For C bl
!;k;l;h given by (5.16), (5.19) we have the estimates:
1;k;l;h C bl
1;k;l
!;k;l;h C bl
!;k;l
lower
Asymptotic analysis of the laminar viscous
ow over a porous bed 33
Proof. In order to show (5.20), we use the denitions (4.11) and (5.16) in the
momentum equations (4.1) and (5.10) to obtain
1;k;l;h C bl
Z
k;l;h
~
bl
ds
Z
l
r ~
bl
k;l r( ~
bl
k;l
~  bl
Z
l
r! bl
bl
k;l
~  bl
Z
l
r( ~
bl
k;l
~
bl
Z
l
r ~
bl
k;l;h r( ~  bl
k;l
~
bl
Z
l
r! bl
bl
k;l
~  bl
Z
l
r( ~
bl
k;l
~
bl
dx
Z
l
r(! bl
bl
Z
l
r! bl
bl
k;l
~
bl
Here, the rst term is of order O(h 2 ) by (5.12), the second term can be estimated
as
Z
l
r(! bl
bl
Z
l
r(! bl
k;l;h
~
bl
k;l
by (5.15), and, if one applies (5.11), the third term can be estimated as
Z
l
r! bl
bl
k;l
~
bl
Z
l
r(! bl
bl
k;l
~
bl
k;l;h )+
Z
r! bl
k;l;h
Next, if we use (5.15) and (5.13), we can see that the resulting terms are of order
respectively. Thus, we have shown (5.20).
In order to show the estimate (5.21), we rst note that
!;k;l C bl
Z
l
(! bl
k;l;h (y)) ~ (y) dy +O(e
lower jlj ) (5.24)
by (4.21), (4.14), and (5.19). Next, let ~
l ) be the vector eld given by
Lemma 13 which satises div ~
l , and


l
34 Willi Jager, Andro Mikelic and Nicolas Neu
Furthermore, theorem 17 yields the existence of an interpolating vector eld ~
~
H h which fullls
l


l
l


l
By using this, we can estimate the rst term from the right hand side of (5.24) as
follows:
Z
l
(! bl
Z
l
(! bl
dy
Z
l
r(! bl
dy
Z
l
r(! bl
Z
l
r(! bl
Z
l
r(! bl
Z
l
r( ~
bl
k;l
~
bl
Here, the rst term is of order C(k by (5.13), (5.26), and (5.25). The
second term can be estimated by replacing rst ~
' h by ~
' which can be done up
to an error of order C(k due to (5.12), (5.27), and (5.25). For the rest, a
partial integration yields
Z
l
r( ~
bl
k;l
~  bl
Z
l
k;l
~
bl
k;l;h )~' dy  C(k
if one uses (5.15) and (5.25). Thus, also (5.21) is proved.
As mentioned at the beginning of this section, the practical implementation of
domain partitions with nonlinear element mappings is rather complicated. It is
usually easier to approximate a
domain
with curved boundaries by polygonal
domains
h which then can be partitioned in triangles and/or quadrangles (see
for example [17]). However, it can be shown that, on a discrete level, this simpler
approach is equivalent to our theoretically more convenient formulation up to a
quadrature error of optimal order, see [18].
Asymptotic analysis of the laminar viscous
ow over a porous bed 35

Figure

3: Symmetric cell and its coarsest grid.
6 Numerical Results
We are now ready to demonstrate our method for computing the constants C bland C bl
on specic examples. First, we consider the symmetric geometry shown
on the left part of Figure 3. The boundary @Z  is a circle with radius 0:25 and
center at (0:5; 0:5). On the right-hand side of Figure 3, the initial polygonal
approximation and the initial grid T h0 with quadrangle elements (see section 5)
is depicted. T h0 is then uniformly rened, yielding further grids T
. By
using the discretization described in the previous section, we obtain for every
a system of linear equations which must be solved to obtain the
discrete solution f ~
bl
k;l;h g.
Since the arising linear systems are very large, we have applied the multigrid
method which is known to be of optimal complexity for a large range of problems,
see [6]. However, due to the pressure stabilization and the polygonal approximation
of the smooth boundary
l
k=1 (@Z  (0; k)), we are not in a Galerkin setting.
In this case, it is known that the simple multigrid V-cycle (one coarse-grid correction
between pre- and post-smoothing) does not have to converge independently
of the number of levels. We therefore used the W-cycle (two coarse-grid corrections
between smoothing). Both pre- and post-smoothing are done by two steps
of a block incomplete decomposition where each block contains the unknowns of
one grid node (the corners of the elements). And indeed, our numerical observa-
36 Willi Jager, Andro Mikelic and Nicolas Neu
(j 1)

Table

1: Results for the symmetric cell
tions conrm that this method is robust with respect to the number of levels and
variations of the parameters k; l, see also Table 5 below.
In

Table

1, the results for a computation with are shown. Starting
from the coarsest level which contains 20 elements, we rene 5 times, which yields
a grid with 20480 elements (61440 unknowns). On each level
the discretized equation and compute the approximations C bl
1;k;l;h and C bl
!;k;l;h given
by (5.16) and (5.19) (with the choice (5.17)). The value given for
computed by polynomial extrapolation. We see that both C bl
1;k;l;h and
!;k;l;h converge to limit values with rate O(h 2 ) which we expected from Proposition
19. As shown in Remark 9, C bl
must be zero. Since our grid is symmetric,
the solution of the discrete problem also has the same symmetry property, such
that the approximations C bl
!;k;l;h are zero up to machine accuracy.
Let us now assume that the porous part is generated by the unsymmetric cell
from

Figure

4, where the boundary curve is given by the ellipse
0:5 0:25
For this domain we obtain the results shown in Table 2. Both C bl
1;k;l;h and C bl
!;k;l;h
converge again with order O(h 2 ), which is in accordance with Proposition 19. The
Asymptotic analysis of the laminar viscous
ow over a porous bed 37

Figure

4: Unsymmetric cell and coarsest grid.
(j 1)
(j 1)

Table

2: Results for the unsymmetric cell
l

Table

3: C bl
1;k;l (extrapolated) for varying k; l.
l

Table

4: C bl
!;k;l (extrapolated) for varying k; l.
error arising from the cutting of the domain is not noticeable any more already
2. This is shown in Table 3, where the results for varying k and l are
given. Note that even the values for are accurate up to the extrapolation
error.

Figure

5 shows the three solution components  bl
From here and
from the values of C bl
in

Table

2 it is obvious that a pressure jump occurs inside
the boundary layer.
Finally, Table 5 shows the convergence rates of our multigrid iteration. As
we expected from the discussion above, it is perfectly robust with respect to the
number of levels. So far, we did not observe a signicant dependency on k or
l, even if one might expect some deterioration, since the error estimates from
Proposition 19 are needed in the usual W-cycle convergence theory. The reason
might be that sharper error estimates are possible in suitably weighted norms.
Asymptotic analysis of the laminar viscous
ow over a porous bed 39

Table

5: Multigrid convergence rates.

Figure

5: Detail from a picture showing  bl
Asymptotic analysis of the laminar viscous
ow over a porous bed 41
7 The eective equations
It turned out that ~  bl stabilizes exponentially fast to a constant vector (C bl
(and to 0 for y 2 !1), which translates into ~  bl;" tending to "(C bl
for This forces us to consider the corresponding counter-
ow in the
channel, which is described by the following 2D Oseen-Couette system
in
d
in
div ~
in
~
If we assume that Re :=
is not too big, the problem (7.1){(7.4)
has a unique solution in the form of 2D Couette
ow ~
Following the ideas from [9], we write down the correction of order " for the
velocity. Essentially, it corresponds to the elimination of the tangential component
of the normal stress at , caused by the approximation by Poiseuille's
ow and its
contribution to the energy estimate. The correction reads
bl ( x
The correction of the velocity eld by the oscillatory boundary layer velocity ~
bl;"
involves the introduction of the boundary layer pressure eld ! bl;" . Hence, as usual
for
ow problems, it is necessary to correct the pressure eld simultaneously. The
corresponding pressure correction reads
Here, p 1;" is an appropriate regularization of the eective pressure p in the porous
42 Willi Jager, Andro Mikelic and Nicolas Neu
bed, dened by
in
where the permeability tensor K is dened as
Z
Y
Here, the w j
are solutions of the auxiliary problems
div y ~
~
Z
Y
The pressure eld p is a C 1 -function outside the corners. However, due to
the discontinuities of the traces at (0; 0) and at (b; 0), its gradient is not square
integrable
in
2 . We must regularize the values at the upper corners and p 1;" is
such a regularization, satisfying
").
Let us now introduce the dierence U "
0 between the velocity eld ~u " and its
expansion up to order O("), i.e. we set
~
bl;" @v 0@x 2
Then, after constructing the corresponding outer boundary layer, it was proved
in [10] that
jr ~
Z bj ~
Z bZ Hj ~
Asymptotic analysis of the laminar viscous
ow over a porous bed 43
It should be noted that the presence of the logarithmic term is a consequence of
the corner singularities in the eective pressure. It was proved in [9] that in the
absence of the boundary singularities, the above estimates hold without the log "
term. Therefore, in the interior of the
domain
the expansion (7.5) is of order
O(" 3=2 ). Globally, it is of order O(" 3=2 j log "j).
The estimates (7.16)-(7.18) are sucient for calculating the next order correction
at the interface . The estimate (7.18) gives us a possibility of approximating
the velocity values at  by an oscillatory velocity eld. Computationally, it is not
very useful.
In view of the problem setting in [9], the Beavers{Joseph law corresponds to
taking into the account the next order corrections for the velocity. In fact, we
formally get on the interface :
"j log
"C bl
x
x
Integrals of the absolute value of the right hand side are not small, since our
estimates do not give any pointwise control of ru on . Nevertheless, the right
hand side of (7.21) is small in the appropriate Sobolev norm, of a negative order.
The precise estimates are in [10]. Hence, we get the eective law
which is exactly Saman's modication of the Beavers and Joseph law from (1.2)
with 1
K
44 Willi Jager, Andro Mikelic and Nicolas Neu
Beavers-Joseph profile
bl
bl4
Poiseuille
profile1
bl
bl

Figure

Beavers-Joseph versus Poiseuille prole.
Let us introduce the eective
ow equations
in
1 through the following boundary
value problem: nd a velocity eld ~u e and a pressure eld p e such that
in
div
in
Under the same assumption of laminarity as for problem (7.1){(7.4), the problem
(7.23)-(7.28) has a unique solution
"C bl
Asymptotic analysis of the laminar viscous
ow over a porous bed 45
(see

Figure

and the eective mass
ow rate through the channel is
Z Hu e
where C bl
By using the theory of the very weak solutions for the Stokes system, the
following approximation properties of f~u e ; are proved in [10]:
The estimates (7.32)-(7.34) justify Saman's modication of the law by Beavers
and Joseph. Furthermore, we are able to calculate their proportionality
it is equal toC blK 1=2 . We note that (7.28) is not the only possible interface
law. If one replaces Beavers and Joseph's law by u e
(0) the estimates
remain valid. However, such a condition involves the knowledge of
the zeroth order approximation v 0
1 , and is not really an interface condition.
At this stage, we can consider the approximation to the channel
ow as satis-
but we still must determine the ltration velocity and the pressure in the
porous medium. We have already mentioned the in
uence of corner singularities
on the solution of the problem (7.8)-(7.10). In order to avoid the discussion of such
eects, we limit ourselves to the behavior in the interior of the porous medium.
The inertia eects are negligible and we can use the theory from [9], Theorem 3.
In the interior, all boundary layer eects are exponentially small and we have
~
where ~
are dened by (7.12). Hence, the ltration velocity is given
through the Darcy law
in
Far away from the corners, the pressure eld p approximates p " at order O(
").
It can only be determined after we have found the eective pressure eld in the
channel and the stabilization constant C bl
giving the pressure dierence at innities
for the auxiliary problem. As shown in Section 6, this stabilization constant is
generally dierent from zero and we must use the interface law
This shows that, contrary to the intuition, the eective pressure in the system
channel
ow/porous medium is not always continuous. Thus, the continuity assumption
for the eective pressure from [11] is not correct in general.
Finally, let us note that we have some freedom in choosing the position of the
interface. It could be set at any distance of order O(") from the solid obstacles.
Intuitively, the law of Beavers and Joseph should be independent of such a choice.
This invariance result can be established rigorously, and we refer for details to [13]
where it was proved that a perturbation of order O(") in the position of the interface
implies a perturbation in the solution of order O(" 2 ). Consequently, the
eective law doesn't change and the physical intuition is conrmed.
Asymptotic analysis of the laminar viscous
ow over a porous bed 47



--R

Finite element methods for problems with rough coe
Boundary conditions at a naturally permeable wall

Lectures on the Mathematical Theory of Multiphase Flow
Equations et ph
of Computational Mathematics
Eine robuste und e



On boundary conditions for uid ow in porous media
Some Methods in the Mathematical Analysis of Systems and Their Control

PhD thesis
On the boundary condition at the interface of a porous medium
Theory and Numerical Analysis


Curved elements in the

--TR
finite elements;unbounded domain;homogenization;multigrid;boundary layer;periodic structures;interface law;beavers-joseph;navier-stokes