Numer. Math. 63, 195 211 (1992)

Numerische Mathematik

9 Springer-Verlag1992

The L.P.D.E.M., a mixed finite difference method for elliptic problems M. Jai 1, G. Bayada 1'2, and M. Chambat 3 1Centre de math. U.R.A. 740 C.N.R.S., I.N.S.A. de LYON Bfit 403, ]7-69621 Villeurbanne Cedex, France 2 Mecanique des contacts U.R.A. 856 C.N.R.S., I.N.S.A. de LYON Bfit 113, F-69521 Villeurbanne Cedex, France 3 L.A.N., Math. U.R.A. 740 C.N.R.S., Universit5 LYON I, F-69622 Villeurbanne Cedex, France Received March 14, 1991/Revised version received March 20, 1992

Summary. A "mixed" finite difference method is analyzed for solving certain elliptic problems. This method, called L.P.D.E.M. (Locally exact Partial Differential Equation Method) was initially proposed in the frame of hydrodynamic lubrication. Convergence is obtained. Relations between this scheme and homogenization theory are also discussed. For a one-dimensional elliptic equation with no zero-order term and in conservative form, this method is an exact one. Some numerical results will also be given.

Mathematics Subject Classification (1991): 65N05, 65N10, 65N30, 76D08 1. Introduction We consider the problem: (1.1)

- div(D(x).gradu) + 7(x)u = f ( x ) u=O

xeO

on F

where (1.2)

D(x)=[~l(X) 0

c ~02 ( x])

withcq(x)>0

i= 1,2;7(x)>0

and Q is a bounded domain in IR2. For sufficiently smooth data, (1.1) has a unique solution u(x). The L.P.D.E.M. procedure, as proposed by Pan, Perlaman and Li [-7, 8] for hydrodynamic lubrication problems, appeared at first glance to be closed to a mixed finite element method by introducing the flux field variable (1.3)

Correspondenceto: M. Jai

q~=

--D(x) grad u

196

M. Jai et al.

Then an integration of the flux equation (1.3) along the xl axis allows us to obtain a relation between nodal values of the unknown u and the flux field 4~: xl +ht

(1.4)

u(xl + h l , x l ) -

U(Xl,X2)=-

S 4~l(s'x2)~l(s'x2) ds xt xl

(1.5)

U(Xl,X2)-u(xl-hl, X2)=-

~ ~ ( s , x2)~;l(s, x2)ds xt - h i

The primary approximation to be imposed is that the flux components are adequately represented by linear truncated Taylor expansions: c3 (~1(S, X2) = ~ I ( X 1 , X2) + (S - - X 1 ) ~ X l g l ) I ( X 1 , X 2 )

Substituting this equation in (1.4) and (1.5), we obtain a system of two equations 0 0 with two unknowns cb~ and ~ (~)1 which gives ~7-- q~l as a function ofu alone. The tl.x 1

same study of the second flux component, in the x2 direction, gives ~x-~q~: as a function of u only. Using the divergence equation, div q, + 7(x). u = j ; we generate a linear system where the only unknown is u. Due to this systematic splitting between u and q~, the L.P.D.E.M. method appears to be a finite difference scheme, in the final analysis. Using classical tools, we can establish that L.P.D.E.M. is second order accurate for sufficiently smooth u(x). An interesting one-dimensional property is that the exact point wise solution for y ( x ) = 0, is found independently of the number of nodes, if the analytical integration procedures are exact. Furthermore, certain stability properties are pointed out related with homogenization theory. For instance, the approximate L.P.D.E.M. solution for highly periodic coefficients is nothing other than the classical finite difference solution of the homogenized problem. Applying the L.P.D.E.M. approach to a diffusion-convection problem, a nonuniform difference scheme is obtained. But we can introduce something resembling a difference scheme which is uniformly first order accurate and is well suited to singular perturbation problems (see [5]).

2. Notation and generation of the difference scheme

2.1 Notation and basic properties For the sake of simplicity, we assume f2 to be rectangular and the discretization to 1 1 be regular with h a and h2 = - (NI and N 2 a r e non-negative N1 + 1 N2 + 1 integers, see Fig. 1). Set ~~h -~" {Qr = (ih~,jh2), 1 <_ i <- N~, 1 <=j <=N2}

Fh = {Q = (0,jh2), 1 =
A mixed finite difference method for elliptic problems

j Q ~ a

197

h

"-------'----~Q ~ rh

Qf

hl

Q2 IQ

" Ol

Q4 Fig. 1

X = (Xl, X2):

~ +h,

aT(x) = b]-(x)=

S

x,

al(t, x l ) '

~,+h~ ( t - x l ) I --dr; x,

x,

dt

OCl(t' X1)

a;(x) =

dt

S

x*-h, Oq(t, X2)

b;(x)=

~ xl -tll

( t - Xl) dt

2,(tzU)

Integrating with respect to the second variable x2, we find in a similar way a~, b~ Let a~ (x) Ck(X)-k = 1,2 rk(X) = (a~ .b ~ - a~ .b ~ )(x); r~(x)

dk(X)= _(a~ +a~)(x)., rk(x)

ek(x)--a-k(X) k= 1,2 r~(x)

Proposition 2.1. a~ , b ~s and dk are non-negative functions and b ~ , rk, ck and ek are negative ones and dk = -- (Ck + dk) for k = 1, 2. Proof It is obvious from (1.2) and the definition of a~, a ~ . .

[]

Let Q ~ f~h and Qi (i = t, 2, 3, 4) as in Fig. 1. The L.P.D.E.M. operator is defined by: EhU(Q) = q ( Q ) . u(Q3) + el(Q), u(Q1) + c2(Q), u(Q4) + e2(Q), u(Q2) + (dl + dz)(Q), u(Q)

Proposition 2.2. Assuming the solution u of(1.1) in C4(0), then we have for all Q in Qh, 1. Ehu(Q) + y(Q)u(Q) = f ( Q ) + gh(@)(Q) 2. leh(Cb)(O) ] < K . h z where eh( ~ ) is defined in the proof of the proposition and K is a constant independent ,from h. Proof 1. From (1.3), we have

(2.1)

Ou t~x-o = - e - t 4~i i = 1, 2 with 4~ = (4'1,4~2)

198

M. Jai et al.

Integrating equation (2.1), for i = 1, successively in the intervals [xl, xl + h~] and [Xa - hi, xl], we obtain; (2.2)

u(xl q- ]/1, X1) -- U(X1, X2) = --

(2.3)

u(xl, x2) - u(x~ - hi, x2) = -

xl +hi j" ~[)1(S, X 2 ) ~ t I (S, X2) ds xt xt

9 ~(s, x~)~-I (s, x2) ds x1 hi

Let us consider the Taylor expansion of the function s ~ ~ (s, x2): (2.4) qh(s, x2) = 4~1(Xl, x2) + (s - x1)6~Xl(Pl(X1,X2) (')2 X63X,~I(Xli + O ~ ( s - - x 1 ) , X 2 ) ,

q- ~(s - x1) 2

Xl < S < X l

+hi,

0<0~

< 1

On substituting for 4~ in (2.2) from (2.4), we get: (2.5) a +(x). ebl (x ) + b ~ (x). ~ - ebi (x ) = u(xl, x 2 ) - u(xl + hi, x 2 ) - I ~ (x) oxl where 1 xl

(2.6)

I]-(x)=~

+hi

~

02

( s - x l ) Z c ~ ; ' t s , xz). ( ? x ~ , t x ,

+ OT(s-x1),x2)ds

1

Similarly evaluating the right hand side of (2.3) gives: (2.7)

a l (x). ~l (X ) + b l (X). ~xl @i (x ) = u(x I -- h i , x 2 ) - - u(xl, x 2) -- [ l (X)

where I-~(x) is defined as 1~i(x) by integration in [x~ - hi, xl]. Equations (2.5) and (2.7) lead to a system of two equations with two unknowns O 4~1(x) and ~ ~l(x) whose determinant is rl(x) which is strictly negative from Proposition 2.1. The solution is: (2.8)

4~l(x) =

-bA(x)u(xl-hi,x2)+ b+ -~'b2 (x).u(xi,x2) rl

rt

b-; (x)u(xl + h i , x 2 ) + b ~i .I~ - b ~ .l~i (x ) rx rl

(2.9)

--~I(X)

0Xl

= C l ( X ) . U ( X 1 -- h i , x2) + d t ( x ) . u ( x 1 ,

+el(x).u(xl+hl,x2)-a

x2)

.,+ I ~ - a y . I . , + ( x ) F1

A mixed finite differencemethod for elliptic problems

199

The same study in the x2 direction for Eq. (2.1), with i = 2, gives b+2

(2.10)

b~ + b2

(D2(X) = - - - - ( X ) / A ( X 1 , X 2 - - h 2 ) "~ r2

r2

(x). u(x1, x2)

bf(x)u(x~, x2 + h2) + b~,12 - b 2 . I ~ (x) r2

r2

(2.11) (3x--~IJ)2(x) ~- CI(X)' tl(Xl'

X2

--

h2) + dl(x), u(xl, x2)

+ el(x).u(xl,x2 + h2) - a~.12 - a 2 . l ~ (x) r2

(1.3), (1.1), (2.9), (2.11) and the definition of the operator Eh induce (2.12)

E,,u(Q) + y(Q). u(Q) = f ( Q ) + eh(g')(Q) for each Q in Oh

with, 2

eh(qJ)(Q)= ~ a ~ - . I ~ - a ~ . I ~ ( Q ) k= 1

rk

So the first part of the proposition is proved. 2. Using the generalized mean theorem for the coefficients a~,b~, and l~(k = 1, 2)(see [5]), yields:

[a~ .I~ -a~..I~[ < Kh~,s and

Klh~ <=lrdx)l<-_Kzh~ k= 1,2 which gives

]~..(4,)(x)l __
[]

2.2 L.P.D.E.M. discretized problem The finite difference scheme for the elliptic problem is obtained by neglecting the terms of orders greater than 3 in the Taylor expansion of u. The L.PD.E.M. method, based on the same principle is obtained by neglecting the terms of orders greater than 2 in the Taylor expansion of the functions s--+cbl(s,x2) and S ~ (]~2(X1 , S).

Let un be the approximation of u. From (2.9) and (2.11), we have for all

Q, Q,(i = 1, 2, 3, 4) in Oh (Fig. 1). 0 4~1(Q) +

4~2(Q)= EhU(Q) + eh(q~)(Q)

Upon omitting O(h z) and higher order terms, eh(cl))(Q) vanishes and we obtain a new linear algebraic system: (2.13)

f E,,uh(Q) + ?(Q)uh(Q) = f ( Q ) (uh(Q)=O Q6rh

Q~Oh

200

M. Jai et al.

3. The one-dimensional case

Let us consider the problem: S-(c~(x)'u')'+7(x)'u=f(x)

(3.1)

O
1

N

( u ( 0 ) = u(1) = o ,

where a is a positive function in C2[0, l], ~ is a non-negative function in C~ and f is given in C~ l]. Define the grid x o = 0 < x l = h < . . .


1]

( h - N + l 1) ~

[0, 1]. By introducing the new variable r = -c~. u' and using the method described in the proof of Proposition 2.2(1.), it is proven that the solution u of (3.1) satisfies: a+ -r

a++

a-

(x).u(x) + E-(x).u(x + h) + ~;(x).u(x)

(x). u(x -- h)

r

r

= f ( x ) + gh(C~)(X) with a+(x)=

x+h

5 a-l(s)ds;

a-(x)=

x

b+(x)=

i

a-l(s)ds

x-h

x+h ~ (s-

x)c~-'(s)ds;

b-(x)=

x

i

(s- x)c~-l(s)ds

x-h

r(x) = a + . b - ( x ) -

a - . b +(x) ,

while the associated linear problem is: I

~

+ (Xi)" Uh(Xi- 1 )

=f(x,)

a+ - -- -{- a r

a-

(xi). uh(xi) + --(xi), uh(xi+ 1) + y(x3. uh(x~) r

l
I, uh(0) = uh(1) = 0 ,

Theorem 3.1. Let c~be a positive function in C2 [-0, 1] and f be in C~

1]. Then the

approximate L.P.D.E.M. solution of:

(3.2)

--(a(x).u')'=f',

u(0)=u(1)=0

is exact. Proof Applying the L.P.D.E.M. method by introducing the new variable

4) = - ~ . u' - f

(r is a constant according to (3.2))

we have u ' = - o ~ - I 4) - o~-1. f

A mixed finite difference method for elliptic problems

201

Integrating this equation successively in the interval [x, x + h] and [x - h, x], we obtain (q~ is a constant): x+h

u(x+h)-u(x)=-4).

x+h

~ :t-l(s)ds -

~ f . ~ - ' ( s ) ds

x

x

i

u(x)-u(x-h)=-q~, x

x

a-l(s)ds-

~ f.a

h

x

l(s)ds.

h

Eliminating q~ we obtain the following equation:

a + ( x ) . u ( x - h ) - ( a + + a ) . u ( x ) + a - ( x ) . u ( x +h) =

(/:)(xi ) ~z-1 .

o~-1

-

_

f.o~-i

x

.

~-1

x

h

)

.

As no approximation was used in deriving the difference scheme, the solution is nothing other than the exact one. []

Remark 3.1. Considering the variational formulation of (3.1), it can be shown that the results of Theorem 3.1 are still valid if c~is assumed to be only measurable and f i n Hi(0, 1).

4. Comparison with other procedures In order to compare the L.P.D.E.M. with the generalized finite element method [1] and the finite volume element method [6], we consider the one-dimensional model equation: (4.1)

-- (~(x).u')' = i f ( x ) ,

0 < x < 1;

u(0) = u(l) = 0 .

4.1 Generalized finite element method Since the standard finite element method, with piecewise polynomials of degree 1, gives a zero order accuracy for equations with coefficients assumed to be only measurable, Babuska and Osborn [1] introduced the generalized finite element method which gives first order accuracy. Applying this method for equation (4.1), we obtain:

a+(x).uh(x- h ) - (a + + a-).uh(x) + a-(x).ua(x + h) x+h

I f(s)dsx

i f(s)ds x-h

So both the L.P.D.E.M. and the generalized finite element methods have the same /eft hand sides, but the right hand sides are different. We refer to Sect. 8 for numerical tests.

202

M. Jai et al.

4.2 Finite volume element method

Let xi be the nodes of the basic grid, and consider another grid {ci}l__
~(ci- 1). u,,(xi- ~) + (c~(ci- 1) + ~(ci)). uh(xl) - c~(e3, u,,(xi + ~ ) = h. (f(ci) - f ( c i - 1 ))

Now, both the left and the right hand sides, and the coefficients of the linear system, differ from those of the L.P.D.E.M. scheme.

5. Discrete maximum principle

Here we will summarise some results concerning the properties of the operator Eh especially the discrete maximum principle [3]. These results are very close to the known results for the finite difference approximation for elliptic equation [4]. Further details may be found in Jai [5]. Lemma 5,1. There exists one function v > O, sufficiently smooth on O, depending on ctt and ~2, such that, for h <=ho: EhV ~ |

on

~ .

Lemma 5.2. (Discrete maximum principle). Let f and ,q be two re qular fimctions defined respectively on ~ and F and consider the discrete problem EhUh(Q) + y(Q).uh(Q) = f ( Q ) uh(Q)= g ( Q ) Qerh

Qef2h

l f f ( Q ) < O,for Q~f2h, and g(Q) < O,Jor Q in Fh, then Uh(Q) <=O,for Q in f2h W Fh.

Theorem 5.1. Let f and 9 be two regular functions defined respectively on ~2h and Fh, and let Uh be the solution of EhudQ) + 7 ( Q ) . u h ( Q ) = f ( Q ) Uh(Q) = g(Q) QeFh

Q~f2h

I f there exist two constants A and B, > O, such that If(Q)l < A, for Qe~2h, and I 9 ( Q ) I < B, for Q in Fh, then

[Uh(Q) [ < = B + A . I I v H ~ ,

for Q in Qh ~ Fh .

(where v is the function defined in Lemma 5.1).

From Varga [10] (Corollary 1, p. 85), we have the following Lemma 5.3. IrA = (al, j) is a real, irreducibly diagonally dominant n x n matrix with ai.j < Ofor all i ~e j, and ai, i > Ofor all 1 <- i <_ n, then ,4 -1 > O.

A mixed finite difference method for elliptic problems

203

6. Nodal error estimates Theorem 6.1. There is a unique solution Uh of(2.13). P r o ~ Consider the matrix Lh of the linear system (2.13). Combining Proposition 2.1 and L e m m a 5.3, we conclude that Lh is nonsingular and LiT 1 > 0, so the result follows. [~

Theorem 6.2. l f u is the exact solution of(l.1) and uh the solution of(2.13), then we have the error estimates: l u ( Q ) - uh(Q)l _-< K . h 2,

Q6~2h W Fh ,

where K is independent of h. Proof F r o m (2.13) and Proposition 2.2(1.), we obtain with l), = u - uh: Ehrh(Q) + ~'(Q).rh(Q) = e,h(~)(Q) rh(Q)=O Q e F a

Q~Qh

Applying T h e o r e m 5.1 to rh, with B = 0 and A = max [c,h(rb)(Q)], we obtain Qsf2h

Irh(O)[ < IIv[l~ . m a x [e,h(qb)(Q)[

VQ6f'2 h

Q~2

and using Proposition 2.2(2.), the result follows.

E3

7. Links with homogenization theory The features of the L.P.D.E.M. coefficient ai-+ and numerical experiments (see Sect. 8, Example 3) suggest some links with the homogenization theory.

7.1 Uniform periodicity Let us consider the problem -

-

+

on f2

(7.1) u=0

on F

where ~ ( i = 1, 2) are positive periodic functions in C2(O). Let e be the period, assumed to be small. Using directly the classical numerical methods (finite element, finite difference) to solve problem (7. l) gives very bad numerical results for small e, Searching for a way [2, 9] an asymptotic expansion of u ~ is introduced, and the homogenized problem satisfied by the first term of the expansion (u0) is then considered. Introduce y~ = xi/e~ and assume ~ to be written as a function ofyi only a~(x)=~i(yi)

i = 1,2

M. Jai et al.

204

We shall write: US(X) = Uo(X ) -t- g. R I ( X , y ) "+- / 3 2 U 2 ( X , y ) -1- ' ' "

It is well known that for small e, the solution u s tends toward the solution u0 of the corresponding homogenized problem: 1

(7.2)

I

~2

1 02

~z* c3x2 Uo Uo = 0

on

~, axzZUo + ?'(x). Uo = f ( x )

on ~2 ,

/" ,

dy

where c~* = o~/i ~

i = 1,2.

Theorem 7.1. Assume that the discretization meshes satisfy

hi = ml.e, (mi and m2eIN*); then the L.P.D.E.M. solution U~hO,/'(7.1) is nothing other than the classical finite difference approximation of the homogenized problem (7.2). Proof. We assume a rectangular domain f2. Using the L.P.D.E.M. scheme for solving problem (7.1) (see Sect. 2.2), we obtain the linearized problem (2.13) which is nOW:

fE~u~(Q)+?(Q).u~(Q)=f(Q)

(7.3)

(u~((?)

= o

0erh

QeY2h

,

where E~, is defined from the coefficients +

+

a T , b T , ri 9 . 9

which depend upon e. Using periodicity of a~ and performing a few classical change of variables, the dependence of these coefficients upon e are explicitly given by: a~(x)=aT(x)=e.c~*,

ri=--(mi.e)3.(o~*) 2 i = 1,2

so that (7.3) may now be written: -

-

1 u~(xl,l-i,Xz,j)-

-

-

- -

2 . U ~ ( X l , i , XZ,j) + U~h(XI,i+I,X2,j)

.

~*

(ml .a)2

1 u~,(xl,l, x 2 , j - 1 ) - 2.u~,(xi,i, x20) + u~(xl,,,x2,j+i) a* (m2 .e)2

-

-

~

+7(xi,l, x2,j).u~(xi,i, x2 d ) = f ( x l , i , x 2 , j )

l_
1 =
which is the classical approximate 5-point finite difference scheme, with hi = ml. e, i = 1, 2 of the homogenized problem (7.2). []

A mixed finite difference method for elliptic problems

205

7.2 Non-uniform periodicity The results of Theorem 7.1 can be extended to the more general non-uniform homogenized problem [2, 9]. In this case the right hand side of (7.1) becomes: ~

x

and the coefficients ~ are replaced by fl~(x, xi/e) with the non-uniform periodicity properties: For each x, fli(x, xl/e) are e-periodic with respect to the variable x;. The following result is established in Jai [5]. Theorem 7.2. Assume hi=mi.e ( i = 1,2 and ml,m2~N*).

The approximate L.P.D.E.M. u~, is a first order approximation of the 5-point finite difference solution of the corresponding homogenized problem. 8. Numerical results

In this section we illustrate the improvement offered by the L.P.D.E.M over the classical methods (finite difference, finite element, and so on) by considering five examples.

Example 1. Lubrication problem (Reynolds equation) We consider problem (8.1) to illustrate Theorem 3.l. Both L.P.D.E.M., finite difference (F.D.) and finite element (with piecewise polynomials of degree 1} are used to obtain approximate solution to the problem d (H3(x)p') ' = H ' 0 < x < 2 7 r where ' acts for dx (8.1)

{

p(0) = p(2~) = 0 ,

-

where H(x) -- 1 + ~.cos x, 0 < ~ < 1 represents the film thickness and p the pressure in a lubricated device. The results are shown in Fig. 2 (a = 0.8) and in Fig. 3 (~ = 0.95). For all values of ~ and mesh length h, the LP.D.E.M. approximation gives the exact solution (as predicted by the theory in Theorem 3.1). The result is particularly interesting for high values of a. In this case (Fig. 3), to obtain a good approximation of the solution by L.PD.E.M. we can take h = ~z/5, whereas a mesh size less than ~t/50 is needed for the finite difference and finite element methods.

Example 2. Equation with zero-order term Whereas the "exact approximation" properties have been proven for the L.P.D.E.M. approximation only for equations without zero order term, the following Table 1 exhibits the good behaviour of LP.D.E.M. versus finite difference (F.D.), finite element (F.E.) and generalized finite element (G.F.E.) (with piecewise polynomials of degree 1). The last of these was introduced by Babuska and Osborn [1] to generalize the standard finite element method. The model problem is: ~" - ( ( 1 + ~ . c o s ( 1 0 ~ x ) ) . u ' ) ' + ( 1 + x ) Z u = l + x O < x < 2 n ( u(0)=u(2n)=0 ~=0.95 As no exact solution is known in closed form, the value of the exact solution u at any grid point x~ is taken to be the value at x~ of the L.P.D.E.M. approximate solution obtained with h = 1/320.

(8.2)

206

M. Jai et al.

-

2.5

-

Exact Solution

---*----

2

F.E

/.

1.5

\

1t

O.fi 0

o c5

~

c::::;

--

--

~

c-4

~

,-,5

Fig. 2. H ( x ) = 1 + 0.8 c o s ( x ) h = 0.314

16"r

Exact Solution

14+

t3------ LP.D.E.M

12+ --

,t

F.D

- -

*

F.E

10+ 8+ 6+ 4+ 2+

Oi

I

i

"1"

Fig. 3.

H(x) =

~

I

I

I

I

I

1 + 0.95 c o s ( x ) h = 0.314

Figs. 2, 3. Example 1: C o m p a r i s o n of the Reynolds equation solution by the various methods

We give in Table l, the L ~ error = - m a x i m u m over i = 1,2 . . . . . N of lu(xi) - uh(xi)l (uh is the a p p r o x i m a t e solution obtained by the various methods).

Example 3. Comparison with the homogenized solution In order to illustrate the results of T h e o r e m 7.1, p r o b l e m (8.3) was run for various methods. Table 1. The discrete L ~~ error of the approximate solution by the various methods Method

h 1/10

L.P.D.E.M. F.D. F.E. G.F.E.

8• -4 1.8x 10 - t 1.8 x 10 -1 1.2• -2

1/20

1/40

4.8 x 10 -4 9.7• -2 1.2• i 4 x 1 0 -3

2.4• 2.7• 4.8• 1.2•

-5 -2 3

1/80

1/160

8• -s 2• -3 1.1• 3::<10-4

1.8 x 10 - s 2.5• -5 3 x 10-2 7.6•

A mixed finite difference method for elliptic problems

207

Table 2. Comparison of the various methods. Relative error with respect to the exact homogenized solution h

Method

h = c = 1/20 h=2.t:=l/10 h=5.c/2

(8.3)

L.P.D.E.M.F.D.

F.E.

H.F.D.

7.4x 10 -3 6.5x10 2 1.7x10 z

4.2x 10-1 4.2x10 1 4.2x10 1

3x10 4 1.2x10 3 2.9x10 4

1.32 6.3x10-1 3.6x10-1

+x

J" - ( ( 1 + a . c o s ( 2 7 z x / e ) ) . u ' ) ' + u = l [u(0)=u(1)=0 c~=0.8 0 < ~ : < 1

O
The homogenized solution Uo is associated with the following equation: 1 a* U or + Uo = l + x,

1 Uo(0) = Uo(1) = O;

a*

N/1 __ ~2

whose exact solution is

Uo(X) = C.(e -j~x - e

,j,*x) + I - e -,/"*x

where C-

e - \.a* _ 1 e V'a* -- e-

\'a*

Denoting by u~, uj,VD and u rE, the solutions obtained by the L.P.D.E.M., finite difference (F.D.) and finite element (F.E.) schemes respectively, and by u,. o the finite difference approximation (H.F.D.) of Uo, we will consider the relative error in L z discrete norm in comparison with the L.P.D.E.M. method, to the homogenized solution. T h o u g h Theorem 7.1 predicts correct behaviour of the L.P.D.E.M. solution only when h is a multiple of e, it appears in Table 2 that this result is still valid for other values of h as h = 5e/2.

Example 4. Equation with discontinuous coe~cients It is well k n o w n that for an equation with coefficients assumed to be only measurable, and for the standard finite element method with piecewise polynomials of degree 1, the following discretization errors hold: IlU -- UhilH, <--_K . h ~

Ilu -

u~llL2 <= K . h ~

In order to increase the accuracy order, Babuska and O s b o r n introduced in [1 ] the generalized finite element method, which gives a first order accuracy. Although the preceding theory has been established for smooth coefficients, the basic L.P.D.E.M. a p p r o a c h is also well suited to discontinuous coefficients problems [5]. The problem on which the numerical experiments are conducted has been introduced by Babuska and O s b o r n [1]:

-(b(x).u')'=

1,

0
1,

u(0)=u(l)=0

with b(x) =

1

O
100

1/2
=

1/2 1

208

M. Jai et al.

T a b l e 3. T h e discrete L ~ e r r o r o f the a p p r o x i m a t e s o l u t i o n by the v a r i o u s m e t h o d s Method

h

L.P.D.E.M. F.D. F.E. G.F.E.

1/5

1/21

1/41

1 . 4 x 10 -13 1.93• -2 1 . 8 4 x 10 - 2 3 . 9 x 10 - 3

lxl0 -t3 5.4• -3 5.3 x 10 - 3 2.6 x 10 - 4

3.7• -12 2.8 x 10 3 2.8 x 10 - 3 7x10 5

Table 3 exhibits the same standard features of exact integration of the L.P.D.E.M. scheme (see Remark 3.1).

Example 5. Results in two dimensional case

We now consider problem (8.4), to evaluate the L.P.D.E.M. for two dimensional calculations. In order to obtain information on the accuracy of this method, we compare the results with those obtained by a finite difference scheme. We consider the problem: =1 [u = 0

in O = ] 0 , 1 [ x ] 0 , 1 [

on F

where ~k(k = 1, 2) are highly varying coefficients. Direct comparison of L.P.D.E.M. and finite difference method is shown in Figs. 4 to 7 for various mesh lengths h. The results are most stable in Figs. 4 and 5, where each C~k(X)is a function of Xk only. This is a consequence of the homogenization procedure implicitly included in the L.P.D.E.M.

9. Conclusion and comments

It seems fair to state that the agreement between theory and experiment is rather good, especially the exact results for a one dimensional equation without zeroorder term. For the elliptical equation with zero-order term, the numerical results seem to suggest that the L.P.D.E.M. scheme is better than the others, by a factor of five. For rough coefficients (homogenization case), the numerical results obtained in one and two dimensions are in agreement with the theoretical analysis (Sect. 7) of the link between the L.P.D.E.M. method and the homogenization theory. The computation time needed to implement the L.P.D.E.M. is exactly the same as in the finite element method with piecewise polynomials of degree 1. We derived, in [5], a variant of the L.P.D.E.M. method for the diffusionconvection equation in conservative form. This variant gives an exact discretized

A mixed finite difference m e t h o d for elliptic problems

209

L.P.D.E.M

9

D

F.D

0.14 0.12 ~.

0.1

o,o8 0.06 ~" 0 . 0 4 0,02 0 5x5

t

t

I

I

10xl 0

20x20

30x30

40x40

l/hlxl/h2

Fig. 4. ~ = 0.5

--

0.7

L.P.D.E.M

[]

F.D

0.6 0.5

~ o.4 "~ 0.3 0.2 0.1 0 5•

I

I

I

I

I

10x10

20x20

30x30

40x40

50x50

l/hlxl/h2

Fig.

5. ~ =

0.9

Figs. 4, 5. The solution of the Reynolds equation. C o m p a r i s o n of the L.P.D.E.M. and the finite difference (FD.). ~k(X)= 1 + ct. COS(20~Xk) (k = 1, 2)

problem for the equation - (a(x). u')' + (b(x). u)' = f ( x )

and gives a uniform difference scheme useful for the singular perturbations problems.

210

M. Jai et al.

=

0.1

L.P.D.E,M

[]

F.D

0.08 ~.~ 0 . 0 6

/

e.. c,q

-

0.04 0.02 0 I 5x5

i

i

i

i

i

10xl 0

20x20

40x40

80x80

1 6 0 x l 60

1/hlx1/h2

Fig. 6 9

0.4

L.P.D.E.M

F.D

0.35 0.3 0.25 r

0.2 " E 0.15

0.1 0.05

0 5x5

I

I

I

I

I

10xl 0

20x20

40x40

80x80

1 6 0 x l 60

1/hlxl/h2

Fig. 7 Figs. 6, 7. The solution of the Reynolds equation. Comparison of the LP.D.E.M. and the finite difference (F.D.). ~k(x) = (1 + ~.cos(20nxx)). (1 + c~,cos(20nxl)) k = 1, 2

References I. Babuska, I., Osborn, J. (1983): Generalized finite element methods, their perforlnance and their relation to mixed methods. SIAM J. Numer. Anal. 20, 510-536 2. Bayada, G., Faure, J.B. (1989): A double scale analysis approach of the Reynolds roughness. Comments and application to the journal bearing. J. Tribol. 111, 323 329 3. Ciarlet, P.G. (1970): Discrete maximum principle for finite difference operators. Aequationes Mathematicae 2, 338-352 4. Godounov, S. (1973): Equation de la physique math6matique. Edition MIR, Moscou 5. Jai, M. (1990): Analyse num6rique d'un sch6ma aux differences finies "mixtes" L.P.D.E.M. et d'un algorithme de cavitation. Thdse Math., Lyon I 6. Mc Cormick, S., Cai, Z. (1990): On the accuracy of the finite volume element method for diffusion equations on composite grids. SIAM J. Numer. Anal. 27, 636-655

A mixed finite difference method for elliptic problems

211

7. Pan, C.H.T., Perlman, A., Li, W. (1987): A new numerical technique for the analysis of lubricating films. 1, Incompressible, isoviscous lubricant. Proc. 13th Leeds, Lyon, Symposium on Tribology, Leeds, pp. 8 12 8. Pan, CH.T., Li, W. (1987): A new numerical technique for the analysis of lubricating films. 2. Precision and convergence. Inter. Mech. Engr., pp. 47-52. 9. Sanchez Palencia, E. (1980): Non-homogeneous media and vibration theory. Lecture Notes in Physics, Vol. 127. Springer, Berlin Heidelberg New York 10. Varga, R.S. (1962): Matrix iterative analysis. Prentice-Hall, Englwood Cliffs, NJ