Journal of Applied Mathematics and Physics, 2014, 2, 88-93
Published Online April 2014 in SciRes. http://www.scirp.org/journal/jamp
http://dx.doi.org/10.4236/jamp.2014.25012
How to cite this paper: Cao, H.T. (2014) A Numerical Approach to a Nonlinear and Degenerate Parabolic Problem by Regu-
larization Scheme. Journal of Applied Mathematics and Physics, 2, 88-93. http://dx.doi.org/10. 4236/jam p. 2014.25012
A Numerical Approach to a Nonlinear and
Degenerate Parabolic Problem by
Regularization Scheme
Haitao Cao
Department of Mathematics and Physics, Changzhou Campus, Hohai University, Changzhou, China
Received October 2013
Abstract
In this work we propose a numerical scheme for a nonlinear and degenerate parabolic problem
having application in petroleum reservoir and groundwater aquifer simulation. The degeneracy of
the equation includes both locally fast and slow diffusion (i.e. the diffusion coefficients may ex-
plode or vanish in some point). The main difficulty is that the true solution is typically lacking in
regularity. Our numerical approach includes a regularization step and a standard discretization
procedure by means of C0-piecewise linear finite elements in space and backward-differences in
time. Within this frame work, we analyze the accuracy of the scheme by using an integral test
function and obtain several error estimates in suitable norms.
Keywords
Nonlinear Degenerate Parabolic Equation, Finite Element Method, Regularization, Implicit Scheme
1. Introduction
Denoting
n
RΩ⊂
(n ≥ 1) as the domain occupied by the porous medium with Lipschitz boundary and (0, T] (0
< T < +∞) as the time interval. In this work we propose a numerical method to the following nonlinear and de-
generate parabolic equation.
( ]
t
b(u)div(A(x)u + g(x, b(u))) = 0, (x, t) Ω 0,T∂ −∇∈×
(1.1)
with the boundary and initial conditions
u = 0 on Ω (0,T],u(x,0)=u00 in Ω,∂× ≥
(1.2)
where the function b(s) is uniformly bounded, continuous and strictly increasing in s. But the degeneracy condi-
tions, b′(s) = 0 and b′(s) = +∞ for some point, are included. It means that the equation includes both locally fast
and slow diffusion.
Throughout this paper we assume, without loss of generality, that b′(s) = 0, s → +∞ and b′(0+) = +∞. Thus the
problem (1.1)-(1.2) often used to describe the flows in porous media, infiltration and multi-phase change see [1]-

H. T. Cao
89
[11] etc. The existence of variational solution of (1.1)-(1.2) has been studied in [5] [6] etc. Here we point out
that the Equation (1.1) is usually obtained after using the Kirchhoff transformation. Such problems can be inves-
tigated through a parabolic regularization of function b(u) or by perturbing the boundary and initial data such
that the corresponding solutions do not take the degenerate points. In the past several decades, there are numer-
ous papers about discussing this problem. For example, in [4] [7], the authors consider a rather general class of
singular parabolic problems; two-phase Stefan problem and porous medium equation are included. By using a
regularization scheme and a parabolic duality technique, a fully discrete numerical approach is proposed and
analyzed. In [8] [9], Richards’ equation is analyzed by a numerical approach also consisting in a regularization
procedure and discretization by means of C0-piecewise linear finite elements (or mixed finite elements) in space
and backward-differences in time. In [10], basing on the maximum principle, the authors considered the porous
medium equation by perturbing the boundary and initial data to overcome the degeneracy. On the other hand,
there are many other papers dealing with such problem by some (linear) relaxation schemes or imposing some
suitable conditions to deal with the degeneracy and nonlinearity [11] 1213].
However, to our knowledge, most of papers have considered one of the degenerate cases, i. e. b′(s) = 0 for
some point but b is Lipchtize continuous (such as the Richards’ equation) or b′(s) > 0 but discontinuous at one
point 0 (such as Stefan problem or porous-medium equation). In our paper, we will deal with the two degenerate
cases simultaneously and present the error estimates of several unknowns by using an integral test function. Due
to the singularity and degeneracy of b, solutions of (1.1) may not be classical; therefore they must be understood
in the sense of distribution [5]. The proposed method in the paper, we first replace b by a regular function bϵ
(whose derivative is bounded by two values depending on regularized parameter ϵ > 0). The second step is a
standard discretization procedure by means of C0-piecewise linear finite elements in space and backward-dif-
ferences in time. Within this frame work, we analyze the accuracy of the scheme by using an integral test func-
tion and obtain several error estimates in suitable norms.
The layout of the paper is as follows. In Section 2, we give out the numerical formulation and main result.
In Section 3, we will prove a prior estimate and the main result of our paper.
Notations:
2 1-1
T0
Ω = Ω (0, T), L(0, T;H () ,H()× ΩΩ
is de dual space of
1
0
H
. Other Soblev space can be
referred to [13]. From now on, C will denote a generic positive constant which is independent of ϵ.
2. Problem Setting and Main Result
For the problem (1.1)-(1.2), we make the following hypotheses upon the data.
Assumption:
(H1) The function b mapping
[0, ] [0,1]+∞ 
is continuous and strictly increasing. And b′(s) = 0, s → +∞,
b′(0) = +∞.
(H2)
(())(1,) :
nn
ij
AAxi jnRM=≤≤ 
is continuous and satisfies
22
1 212
0, ||||,.
n
AR λλ λξξξλξξ ∃<≤≤⋅ ≤∈
(H3)
nn
i
g = (g(x, s))(1 in) : R R R≤≤× →
is continuous in s and fulfill the following identities:
2
|g(y, b(u))g(y, b(v))| C|b(u) b(v)||u v|. − ≤−−
(H4)
2
u0 L()().H
∞
∈ Ω∩ Ω
Remark 2.1. Due to the maximum principle, the solution of the problem (1.1)-(1.2) is great or equal to zero
[9] [10]. In physical model, the function b(s) usually presents the enthalpy of Stefan model or denotes the re-
duced saturation of fluid in porous media.
According to Alt and Luckhaus [5], we have at least that
2-12 1
t0
b(u)L(0, T;H()), uL(0, T;H()).∂ ∈Ω∈Ω
Because b is uniformly bounded, therefore we conclude that
0 -1
b(u) C(0, T;H()).∈Ω
This gives us b(u(·,t))
point wise for every t
∈ (0, T].
Let
H
T
be a decomposition of Ω into a regular conforming finite element mesh with maximal element di-
ameter H. Denote
H
KT
K
∈
Ω=∪
. Then we define finite element space
H
X
as:
0
{( );|0}.
H
XCislinearfor allKand χχ χ
∂Ω
=∈Ω =
Let N be an integer, τ = T/N,
τ iti=
. Our numerical method reads as: find
, ()
HH
i Hii
v Xbv
ε
η ∈=
(i = 1,···,
N), for all
H
X χ ∈

H. T. Cao
90
100
(,)(()( ,),)0,,
HH
ii ii H
AxvgxvI u
ηη χ ηχ
τ
−
−+ ∇+∇==
(2.3)
where
H
I
is the C0-piecewise linear interpolant operator to
H
X
and satisfies: for
2
()uH∈Ω
2() ()
,1 2,
s
s
HLH
I uuCHus
ΩΩ
− ≤≤≤
(2.4)
and
()bs
ε
is defined as
1
1
1
,[0,( )]
()
() .
( ),((),]
ss b
b
bs
bsss b
ε
εε
ε
εε
−
−
−
∈

=
+ ∈+∞

(2.5)
In the above terms, ϵ is small parameter. It is easy to get
1
| ()()|max{,||},'(),
()
bsbssbsb
εε
ε
εε εε
−
−≤≤ ≤
(2. 6)
The Equation (2.3) is a general nonlinear elliptic problem. It can be numerically solved by relaxation iteration
scheme or to apply a linearization scheme first. We point out that the regularization (2.5) of a degenerate prob-
lem is not necessary in implementing numerical analysis with suitable assumptions [11].
Our main result of the numerical approach is present as the following:
Theorem 2.2. Let
, () u bu θ =
and
H
k
v
be the solution of (1.1) and (2.3) respectively. Suppose (H1) -(H4)
hold, then there exist
11
2
2
2222
011
()
()
((),)()() ()
n kk
kk
nn
t tt
HH HH
kk kk
tt
kk
L
L
bv uvuvuvCH θ τε
−−
Ω== Ω
Ω
−−+∇−+ ∇−≤++
∑∑
∫ ∫∫
where the constant C is independent of τ, ϵ,H.
Remark 2.3. Due to the lacking regularity of solution, the result of Theorem 2.2 is not optimal with respect to
time discretization. On the other hand, if b′ is positive and bounded, we can get the classical results.
3. Proof of Main Result
Because of the lacking regularity, the solution of (1.1)-(1.2) must be understood in terms of distributions, as
proposed in [5]. Firstly, we formulate the weak solution of (1.1)-(1.2) as follows:
Definition3.1. We say that u is a weak solution of problem (1.1)-(1.2), if it satisfies the following two identi-
ties:
(1)
2T
b(u) L(Ω )∈
and
2 -1
t
b(u) L(0,T;H() ) ∂∈Ω
with
0
tt
b(u) dxdt(b(u(t))-b(u)) dxdt0,
TT
ϕϕ
ΩΩ
∂ +∂=
∫∫
(3. 7)
for every
2 112
(0,;())(0,;())L THHTL ϕ ∈ Ω∩ Ω
with φ(T) = 0.
(2) For
21
L(0,T;H() ) ϕ ∈Ω
t
b(u) dxdt(()(,()))dxdt0,
TT
Axu gxbu ϕϕ
Ω
Ω
∂+∇ +⋅∇=
∫∫
(3.8)
According to [5], it follows that
Lemma3.2. Assuming (H1)-(H4) hold, if u is a solution of (3.7)-(3.8), we have
12 21
()( )
(0, ;())
[0, ]
max (())()
T
t
LL
L TH
tT
Butbuu C
−
ΩΩ
Ω
∈+∂+∇ ≤
where
()
0
( ())()( ())()
ut
But utbutbsds= −
∫
.
For the discrete scheme (2.3), we also have
Lemma3.3. Assuming (H1)-(H4) hold, if for all
1iN=
,
H
i
v
solves problem (2.3), for any
0kN<<
,
we have
2
2
1()
( ()).
k
HH
ki
iL
B bvdxvC
ε
τ
Ω=Ω
+∇ ≤
∑
∫

H. T. Cao
91
Proof: Let
H
i
v χ =
in (2.3), we have
1
(,)(()() ,)0
H HH
iiiii i
vAx vgxv ηη τη
−Ω Ω
−+∇+∇ =
Summing
1ik=
, we get
001 1
11
(,) (,)(,)(()(,),)0
kk
H HHHHH
kki iiiii
ii
vvvvAxvgxv ηη ητη
ΩΩ −−ΩΩ
= =
−−−+∇+∇ =
∑∑
Denote the above terms by
123
0.III++=
1
22
1 000
1
22
23
11
() ()
(,)(,)()(( ))(()),
, ||.
2
H
i
H
i
kv
HHH H
kk k
v
i
kk
HH
ii
ii
LL
IvvbsdsdxBbv dxBbv dx
IvIC v
εεε
ηη
τ
τ
−
ΩΩ
Ω ΩΩ
=
= =
ΩΩ
≥−−= +
≥ ∇≤+∇
∑∫∫ ∫∫
∑∑
Combining all the above terms, we get the conclusion.
Proof of Theorem 2.2: Let
(( ))bu t θ =
, we integrate (1.1) over time from
1i
t
−
to
i
t
, for every
H
X χ ∈
, to
get
11
1
(,)((),)((),) 0,
ii
ii
tt
ii tt
A xudtgdt θθ χχθχ
−−
−ΩΩ Ω
−+∇∇+∇ =
∫∫
(3.9)
Subtract (2.3) from (3.9), we get
11
11
((),)(()(),) (()(),)0
ii
ii
tt
H
ii iiii
tt
A xuvdtggdt θη θηχχθηχ
−−
−− ΩΩΩ
−−−+∇−∇+−∇=
∫∫
Summing the equality from
1 i=
to
k
, we obtain
11
00 11
((),)(()(),)(() (),)0.
ii
ii
kk
tt
H
kki i
tt
ii
A xuvdtggdt θη θηχχθηχ
−−
Ω ΩΩ
= =
− −−+∇−∇+−∇=
∑∑
∫∫
Summing k from 1 to n again, we get
11
00
1 1111
((),)(()(),)(()(), )0.
ii
ii
n nknk
tt
H
kki i
tt
k kiki
A xuvdtggdt θη θηχχθηχ
−−
Ω ΩΩ
== == =
− −−+∇−∇+−∇=
∑ ∑∑∑∑
∫∫
Setting
1
k
k
tH
Hk H
t
I uv dtX χ
−
= −∈
∫
and
0 000
,( ),()
uH H
eu IububIu
ε
θη =−= =
, to gives
1 11
11
11
00
1 11
11
00
1
(( ),())(()(),())
(()() ,())
((),)(()(),
k ik
k ik
ik
ik
ki
ki
n nk
t tt
H HH
kkki k
t tt
k ki
nk tt
H
ik
tt
ki
ntt
H
kk ui
tt
k
uv dtAxuvdtuvdt
ggdtuvdt
e dtAxuvdt
θη θη
θη
θη θη
− −−
−−
−−
ΩΩ
== =
Ω
= =
Ω
=
− −−−+∇−∇−
+− ∇−
=− −−+∇−
∑ ∑∑
∫ ∫∫
∑∑ ∫∫
∑∫∫
1
11
11
11
)
(()() ,)
k
k
ik
ik
nk t
u
t
ki
nk tt
iu
tt
ki
e dt
ggdte dt θη
−
−−
Ω
= =
Ω
= =
+−
∑∑ ∫
∑∑ ∫∫
Denoting the above formulation by
123 456
,IIIIII++=++
we will deal with each terms in the following.
11
11 1
10011 12
11
11 11 1
111 112 113
(,) (,),
((),) ((),)(,)
,
kk
k
kk k
nn
kk
HH
kk kk
tt
kk
nn n
k kkt
HHH HH
kkkk ksk
t ttt
kk k
Iu vdtu vdtII
Ibvuv dtbvuv dtdsuv dt
III
θη θη
θ ηθ
−−
−− −
ΩΩ
= =
Ω ΩΩ
= ==
=−− −−−=+
=− −+−−+∂−
=++
∑∑
∫∫
∑∑ ∑
∫ ∫∫∫
It is easy to see that
111
0.I≥
Use (2.6) and Lemma 3.3, it follows that

H. T. Cao
92
211
22
22
2
22
112 ()
11 1
() ()
||()( )()( )
kk
kk
nn n
tt
HH H
kk k
tt
L
kk k
LL
I CvuvCuv δετδδεδ
−−
Ω
= ==
ΩΩ
≤+∇−≤ +∇−
∑∑ ∑
∫∫
Considering that
2-121
t0
b(u)L(0, T;H()), uL(0, T;H())∂ ∈Ω∈Ω
, we get
1
1
10
113 ()
()
|| .
k
k
tH
tk
tH
H
Iu vdtC τθ τ
−
−
Ω
Ω
≤∂ −≤
∫
So we get
01
2
2
2
11 1()
((),)( )()().
nk
k
n
tt
HH H
kk k
tt
kL
Ibv uvdtCuv θδετδ
−
Ω=Ω
≥ −−−+−∇−
∑
∫∫
Because
00 0000
|||() ()||() ()|
H
bu bubu bIu
ε εε
θη −≤ −+−
, we have
22
12
12
2
22
24
120 0
() ()
11
()
2
24
1()
| |()()()
( )()().
k
k
k
k
nn
tH
k
LH
t
kk
L
ntH
k
t
kL
I CuHuuv
C Huv
δ τεδ
δε δ
−
−
ΩΩ
== Ω
=Ω
≤++ ∇−
≤+ +∇−
∑∑
∫
∑∫
Using the following equality
22
11 11
2( )(),
nk nn
ki k k
ki kk
aa aa
= ===
= +
∑∑ ∑∑
we can get
2
11
2
2
11
2
22
2()
11
()
22
()
11
()
11
() ()() ()
22
(()() ).
kk
kk
kk
kk
nn
tt
HH
kk
tt
L
kk
L
nn
tt
HH
kk
tt
L
kk
L
IA xuvdtA xuvdt
Cu vdtu vdt
−−
−−
Ω
= =
Ω
Ω
= =
Ω
=∇−+∇−
≥∇−+ ∇−
∑∑
∫∫
∑∑
∫∫
For the term
3
I
, noticing
()
H
ii
bv
ε
η =
, we use (H3), Cauchy inequality and the fact (2.6) to get
111 1
11
22
311 11
22
311 1
() ()
(()(()),())(()(()),()),
||()()( ())()
ikik
ikik
ik
ik
nk nk
ttt t
HH HH
ikii k
ttt t
ki ki
nk n
tt
H HH
i kk
tt
kik k
LL
Iggbvdtuvdtg vg bdtuvdt
ICggb vdtuvdtv
θη
δ τθδετ
−−−−
−−
= == =
ΩΩ
= ==
ΩΩ
=−∇− +−∇−
≤−+∇−+
∑∑ ∑∑
∫∫∫ ∫
∑∑ ∑
∫∫
2
01
2
2
1()
22
11
()
( )((),)().
kk
k
n
L
nn
tt
HH H
ii k
tt
kk
L
Cb vuvdtuvdt δτθ δε
−
=Ω
= =
ΩΩ
≤−− +∇−+
∑
∑∑
∫∫
Similarly, the right hand terms
456
,,III
can be estimated as the following. Use the fact (2.4), we have
2 22
0
122()
22
0
11 22
212
2
4()(,; ())
00
()
22
22
5(,; ())
11
() ()
6
||( ),
||( )()( )(),
||(( ),)
k
n
k
L
kk
n
kk
nn
t
kkukku
L LttL
t
kk
L
nn
tt
HH
ku k
LttL
tt
kk LL
HH
ii
t
Ie dteCH
Iu vdtCeu vdtCH
ICbvuv
θη τθη
δδδδ
τθ
−
Ω
−−
ΩΩ
= =
Ω
Ω
= =
ΩΩ
≤−≤ −≤
≤ ∇−+∇≤ ∇−+
≤ −−
∑∑
∫
∑∑
∫∫
22
0
00
22 22
(,; ())
11
((),)().
kk
n
nn
tt
HH
u ii
LttL t
kk
dteCb vuvdtCH ε τθε
Ω
= =
ΩΩ
+∇+ ≤−−++
∑∑
∫∫
Combining all the terms and choosing
δ
properly, we have
0 11
2
2
0
22
11
()
()
22
1
((),) ()()
(( ),)().
n kk
kk
k
nn
t tt
HH HH
kk kk
t tt
kk
L
L
ntHH
ii
t
k
bvuvdtu vdtuvdt
Cb vuvdtCH
θ
τ θτε
−−
Ω== Ω
Ω
Ω
=
−−+∇−+ ∇−
≤−−+ ++
∑∑
∫ ∫∫
∑∫

H. T. Cao
93
Finally, we use Gronwall inequality to get the conclusion of the Theorem 2.2.
Acknowledgements
The author is supported by the Fundamental Research Funds for the Central Universities (NO.2013B10114).
References
[1] Alt, H.W. (1985) Nonsteady Fluid Flow through Porous Media. Applications and Theory, 3, 222-228.
[2] Damlamina , A. (1977) Some Result in the Multiphase Stefan Problem. Communications in Partial Differential Equa-
tions, 2, 1017-1044. http://dx.doi.org/10.1080/03605307708820053
[3] van Duijn C.J. and Peletier L.A. (1982) Nonstationary Filtration in Partially Saturated Porous Media. Archive for Ra-
tional Mechanics and Analysis, 78, 173-198. http://dx.doi.org/10.1007/BF00250838
[4] Nochetto, R.H. and Verdi, C. (1988) Approximation of Degenerate Parabolic Problems Using Numerical Intergration.
SIAM Journal on Numerical Analysis, 25, 784-814. http://dx.doi.org/10.1137/0725046
[5] Alt, H.W. and Luckhaus, S. (1983) Quasilinear Elliptic-Parabolic Differential Equations. Mathematische Zeitschrift,
183, 311-341. http://dx.doi.org/10.1007/BF01176474
[6] Ladyzenskaya O., Solonnikov V. and Ural’ceva N. (1968) Linear and Quasilinear Equations of Parabolic Type. Tra ns-
lations of Mathematical Monographs.
[7] Nochetto R .H. (1987) Error Estimates for Multidimensional Singular Parabolic Problems. Japan Journal of Industrial
and Applied Mathematics, 4, 111-138. http://dx.doi.org/10.1007/BF03167758
[8] Pop, I.S. (2002) Error Estimates for a Time Discretization Method for Richards’ Equation. Computers & Geosciences,
6, 141-160. http://dx.doi.org/10.1023/A:1019936917350
[9] Radu, F., Pop, I.S. and Knabner, P. (2004) Oder of Convergence Estimates for an Euler Implicit Mixed Finite Element
Discretization of Richards’ Equation. SIAM Journal on Numerical Analysis, 22, 1452-1478.
http://dx.doi.org/10.1137/S0036142902405229
[10] Pop, I.S. and Yong, W.A . (1997) A maximum Principle Based Numerical Approach to Porous Medium Equation.
Proceedings of the 14th Conference on Scientific Computing, 207-218.
[11] Arbogast, T., Wheeler, M.F. and Zhang, N.Y. (1996) A Nonlinear Mixed Finite Element Method for a Degenerate
Parabolic Equation Arising in Flow in Porous Media. SIAM Journal on Numerical Analysis, 33, 1669-1687.
http://dx.doi.org/10.1137/S0036142994266728
[12] Jager, W. and Kacur, J. (1995) Solution of Doublely Nonlinear and Degenerate Parabolic Problems by Relaxation
Schemes. Mathematical Modelling and Numerical Analysis, 29, 605-627.
[13] Adams, R.A. and Fournier, J.J.F. (2003) Sobolev Spaces. 2nd Edition, Academic Press, New York.