The gradient discretisation method for nonlinear porous media equations

The gradient discretisation method for nonlinear porous media equations

Jérôme Droniou School of Mathematical Sciences, Monash University, Clayton, Victoria 3800, Australia.  and  Kim-Ngan Le School of Mathematical Sciences, Monash University, Clayton, Victoria 3800, Australia.
July 2, 2019

The gradient discretisation method (GDM) is a generic framework for designing and analysing numerical schemes for diffusion models. In this paper, we study the GDM for the porous medium equation, including fast diffusion and slow diffusion models. Using discrete functional analysis techniques, we establish a strong -convergence of the approximate gradients and a uniform-in-time convergence for the approximate solution, without assuming non-physical regularity assumptions on the data or continuous solution. Being established in the generic GDM framework, these results apply to a variety of numerical methods, such as finite volume, (mass-lumped) finite elements, etc. The theoretical results are illustrated, in both fast and slow diffusion regimes, by numerical tests based on mass-lumped conforming finite elements.

2010 Mathematics Subject Classification:
65M08, 65M12, 65M60, 76S05

labelkeyrgb0.6,0,1 \definecolorlabelkeyrgb0.6,0,1 \definecolorvioletrgb0.580,0.,0.827 \definecolorshadecolorgray0.92 \definecolorTFFrameColorgray0.92 \definecolorTFTitleColorrgb0,0,0 \DeclareDocumentCommand\RPiD OD O,0 Π_#1(X_#1#2)

1. Introduction

In this paper, we study nonlinear porous media equations of the form


where with , is an open bounded domain of () with boundary , and .

In the case , the equation (1) is the standard model of diffusion of a gas in porous media, which is also called the slow diffusion model. The case describes the flow of an ideal gas in porous media while that of a diffusion of a compressible fluid through porous media. Other choices of exponent appear in different physical situations, such as thermal propagation in plasma () or plasma radiation (). The slow diffusion equation is not uniformly elliptic as it degenerates at unknown points where . An interesting feature is that, if is compactly supported, then at any the support of has a free boundary with a finite speed of propagation, see e.g. [5, 18].

In the fast diffusion case , the equation (1) is relevant in the description of plasma physics, the kinetic theory of gas or fluid transportation in porous media [3, 6, 19]. Since the modulus of ellipticity blows up whenever vanishes, the fast diffusion equation is a singular equation. An essential difference between the fast and slow diffusion cases is that for fast diffusion the solution decays to zero in some finite time depending on the initial data while for slow diffusion the solution decays to zero in infinite time like an inverse power of , see e.g. [4, 7].

If we formally represent (1) as


where , then the above equation is a classical nonlinear diffusion equation. The term in (2) induces the nonlinearity that raises many challenges in the analysis of the porous media equations. These problems have been studied extensively both in theory, see e.g. [4, 17, 19], numerical analysis, see e.g. [11, 12, 13, 1], as well as in numerical approximations (without convergence proof) [14].

In particular, authors in [11] study equation (2) with variable exponent of nonlinearity, i.e. is replaced by where (i.e. in our case). In order to deal with the degeneracy in the problem, an approximate regularized problem is investigated. A space-time discretization scheme using the finite element method in space and the discontinuous Galerkin method in time is proposed for the regularized model (not the original problem). Furthermore, error estimates are obtained with strong regularity assumptions on the solution of the regularized model, which are not expected to be satisfied by classical solutions to (1) such as Barrenblatt solutions [15]. A space and time dependent exponent is studied in [1] using the same method as in [11]. The slow diffusion case () is studied in [12] where a fully discrete scheme is proposed and error estimates are proved with strong assumptions on the solution and the pressure . In [13], the author studies the fast and slow diffusion cases in which a fully discrete Galerkin approximation is considered for the transformed Richards model

where and . Error estimates in non-standard quasi norms and rates of convergence are proved with strong regularity assumptions on the solution of (1). In a recent preprint [16], a general nonlinear diffusion equation on the whole space is considered. The authors propose monotone schemes of finite difference type on Cartesian meshes and obtain -convergence for these schemes.

In this paper, we use the gradient discretisation method to approximate (1) and discrete functional analysis techniques to obtain an -convergence result without assuming non-physical regularity assumptions on the data. The gradient discretisation method (GDM) is a generic framework for the design and analysis of numerical schemes for diffusion models. Using only a few discrete elements (a space, and a function and gradient reconstruction), it describes a variety of numerical schemes – such as finite volumes, finite elements, discontinuous Galerkin, etc. – and identifies three key properties of the discrete elements that ensure the convergence for linear and nonlinear models. Analyses carried out within the GDM therefore yield convergence theorems that directly apply to all the methods that fit into this framework. We refer to the monograph [10] for a detailed presentation of the GDM, and of the methods it contains.

The paper is organised as follows. Section 2.1 recalls notations of the GDM and introduce an implicit-in-time gradient scheme for (1). The definition of weak solutions and the main convergence result are stated in Section 2.2. We provide a priori estimates for the approximate solution in Section 3. Section 4 contains the initial weak convergence of the gradient scheme. The main convergence result, including the uniform-in-time convergence, is proved in Section 5. Numerical results are provided in Section 6 to illustrate the generic convergence results; to this purpose, we choose a particular gradient scheme, the mass-lumped conforming method, and we evaluate its accuracy when approximating the Barrenblatt solution. These tests are presented for a variety of exponents , including both fast diffusion and slow diffusion ranges. A brief conclusion is presented in Section 7, and technical results are gathered in an appendix (Section 8).

2. The gradient discretisation method and the main convergence result

2.1. Gradient scheme

We recall here the notions of the gradient discretisation method. The idea of this general analysis framework is to replace, in the weak formulation of the problem, the continuous space and operators by discrete ones; the set of discrete space and operators is called a gradient discretisation (GD), and the scheme obtained after substituting these elements into the weak formulation is called a gradient scheme (GS). The convergence of the obtained GS can be established based on only a few general concepts on the underlying GD. Moreover, different GDs correspond to different classical schemes (finite elements, finite volumes, etc.). Hence, the analysis carried out in the GDM directly applies to all these schemes, and does not rely on the specificity of each particular method.

Definition 2.1.

is a space-time gradient discretisation for homogeneous Dirichlet boundary conditions, with piecewise constant reconstruction, if

  1. the set of discrete unknowns is a finite dimensional real vector space,

  2. the linear map is a piecewise constant reconstruction operator in the following sense: there exists a basis of and a family of disjoint subsets of such that, for all , it holds , where is the characteristic function of ,

  3. the linear mapping gives a reconstructed discrete gradient. It must be chosen such that is a norm on ,

  4. is an interpolation operator,

  5. .

We then let and .

For any , we define the piecewise-constant-in-time functions , and by: For , for any , for a.e.

Note that is defined everywhere, including at . This will be required to state a uniform-in-time convergence result on this reconstructed function.

If and satisfies , we define . The piecewise constant feature of then shows that


Once a GD has been chosen, an implicit-in-time gradient scheme for (1) is defined the following way.

Algorithm 2.2 (GS for (1)).

Set and let satisfy:


for all ‘test function’ .

In order to establish the stability and convergence of the GS (4), sequences of space-time gradient discretisations are required to satisfy consistency, limit-conformity and compactness properties [10]. The consistency is slightly adapted here to account for the nonlinearity we consider. In the following, we let .

Definition 2.3 (Consistency).

A sequence of space-time gradient discretisations in the sense of Definition 2.1 is said to be consistent if

  • for all , letting

    we have as ,

  • for all , in as

  • as .

Definition 2.4 (Limit-conformity).

A sequence of space-time gradient discretisations in the sense of Definition 2.1 is said to be limit-conformity if, for all , letting

we have as .

Definition 2.5 (Compactness).

A sequence of space-time gradient discretisations in the sense of Definition 2.1 is said to be compact if


and has been extended by outside .

A sequence of GDs that is compact or limit-conforming also satisfies another important property: the coercivity [10, Lemmas 2.6 and 2.10].

Lemma 2.6 (Coercivity of sequences of GDs).

If a sequence of space-time gradient discretisations in the sense of Definition 2.1 is compact or limit-conforming, then it is coercive: there exists a constant such that

2.2. The main convergence result

To state the main result of this paper, we introduce a space of continuous functions for the weak topology of , and we define a weak solution to (1).

For a given , we set . We also set

Definition 2.7 (Weak solution to (1)).

Assume that , and . A weak solution to (1) is a function such that

  1. and in ,

  2. , ,

  3. , and for any


The existence of a weak solution to (1) will be obtained as a by-product of our convergence analysis.

Theorem 2.8 (Convergence of the gradient scheme).

Let be a sequence of gradient discretisations that is consistent, limit-conforming and compact. Then, for each there exists solution to the gradient scheme (4) with . Moreover, there exists a weak solution to (1) in the sense of Definition 2.7 such that, up to a subsequence as ,

  • strongly in ,

  • strongly in .

3. A priori estimates

We first provide a priori estimates for the solution to (4), and then deduce its existence in Corollary 3.2. For legibility, we drop the index in sequences of gradient discretisations, and we simply write instead of .

Lemma 3.1.

Let be defined by (5), and be a solution of (4). Then, for any integer number , we have


Consequently, there exists a constant depending only on , and such that


We choose the test function in (4) to write


For any and , we estimate the first term in the left hand side of (9), starting from


Since is increasing, is convex and thus above its tangent line, which implies for all . Applying this inequality with and , it follows from (10) that


Plugging that in (9) and using a telescopic sum yields (7). The a priori estimates (8) follow from this relation, by writing

Noting that for all and recalling that


we obtain the estimates in (8). ∎

The following corollary guarantees the existence of a solution to (4).

Corollary 3.2.

If is a gradient discretisation in the sense of Definition 2.1 then there exists at least one solution to the gradient scheme (4).


The result is obtained using (8), the properties and , and the same arguments as in [9, Corrollary 4.2]

In the following lemma, we estimate the dual norm of the discrete time derivative . The dual norm on is defined by

Lemma 3.3.

There exists a constant depending on , and such that


We first fix and then, for any , take in (4) such that and for all . Using the Cauchy–Schwarz inequality and the definition of we deduce

where . Taking the supremum over all satisfying , we deduce

Divide by , square, and sum over (and use Jensen or Cauchy–Schwarz to estimate the term involving ) to get

This together with (8) complete the proof of the lemma. ∎

To deal with the lack of global Lipschitz estimates on , we introduce cutoff functions. The definition of these functions is different in the case (fast diffusion) and (slow diffusion), because each of these cases correspond to a certain breakdown in the Lipschitz properties of : for fast diffusion, is not Lipschitz-continuous at 0, whereas for slow diffusion, the Lipschitz constant of explodes at . We therefore define and for any as



The cutoff functions and are globally Lipschitz continuous with respective Lipschitz constants and . Our goal here is to estimate the time-translates of . To achieve this, we will be using the cutoff functions, as well as the following relation.

Lemma 3.4.

Recalling that , for any we have


By noting that for and

we obtain

The above inequalities imply that for any ,

or, equivalently,


Together with the global Lipschitz property of and the monotonicity of , this yields

which completes the proof. ∎

We can now estimate the time translates of .

Lemma 3.5.

Let be a gradient discretisation in the sense of Definition 2.1 and be a solution to scheme (4). Then, there exists a constant depending only on , , and such that:

  • for ,

  • for ,


In this proof, denotes a generic constant that has the same dependencies as in the lemma. For and any we have




To estimate , define

The commutativity property (3) shows that and thus, when ,


By the Hölder and Chebyschev inequalities, (3) and (8), we estimate :


Similarly, we also have estimates for :


We estimate using the same arguments in [9, Lemma 4.4]. Thanks to Lemma 3.4 we have



For any there exists such that . We rewrite by expressing as the sum of its jumps in time, and use the definition of dual semi-norm to infer

Together with the Cauchy–Schwarz inequality (on the sums and the integrals) and (8), this implies


For any