PUBLICATIONS DE L’INSTITUT MATHÉMATIQUE Nouvelle série, tome 91(105) (2012), 137–153 DOI: 10.2298/PIM1205137M

STABILITY AND CONVERGENCE OF THE DIFFERENCE

SCHEMES FOR EQUATIONS OF ISENTROPIC GAS DYNAMICS IN LAGRANGIAN COORDINATES

Piotr Matus and Dmitry Polyakov

Abstract. For the initial-boundary value problem (IBVP) if isentropic gas dynamics written in Lagrangian coordinates written in terms of Riemann invariants we show how to obtain necessary conditions for existence of global smooth solution using the Lax technique. Under these conditions we formulate the existence theorem in the class of piecewise-smooth functions. A priori estimates with respect to the input data for the difference scheme approximating this problem are obtained. The estimates of stability are proved using only restrictions on the initial and boundary conditions corresponding the differential problem. In the general case the estimates have been obtained only for the finite instant of time t < t_{0}. The monotonicity has been proved in the both cases. The uniqueness and convergence of the difference solution are also considered. The results of the numerical experiment illustrating theoretical statements are given.

1. Introduction

Gas dynamics equations play a key role in the mathematical description of the gas processes. The nonlinearity of these equations can generate various physical effects such as shock waves or boundary layers independent of the smoothness of the input data. The question then arises whether there are conditions on the input data which guarantee the absence of any irregularities of the solution.

For the Cauchy problem the necessary conditions for existence of global smooth solution have been obtained by Lax in [5]. Later it was proved in [18] that these conditions are sufficient. In [1], the wellposedness of solutions of 2 × 2 hyperbolic systems with boundary damping is studied under special restrictions on input data. In [10], the unique solvability of special class of hyperbolic IBVP has been proved. By now, the most complete results in investigating the stability of the difference

2010 Mathematics Subject Classification: 65N06; 65M12. Key words and phrases: Isentropic gas dynamics; Riemann invariants; Lagrangian coor-

dinates; initial-boundary value problem; well-posedness and blow-up; finite difference scheme; stability.

The work is supported by Belarusian Republican Foundation for Fundamental Research under grant No.F10R–121.

solution with respect to small perturbations of the input data have been obtained for linear problems of mathematical physics [15]. The main problem in studying the stability of nonlinear difference schemes is the necessity of estimating the difference derivatives. The linearized difference schemes constructed in this paper are monotone. That is why with the use of the maximum principle we are going to get all strong estimates imposing restrictions only on the initial and boundary conditions. Stability and monotonicity of difference schemes for nonlinear scalar conservation laws and multidimensional quasi-linear parabolic equations have been obtained in [6]. In [7, 8], the nonlinear stability of the difference scheme approximating the IBVP for isentropic gas in Eulerian coordinates written in terms of Riemann invariants is investigated.

In this paper we investigate stability of the difference scheme approximating the IBVP for isentropic gas in Lagrangian coordinates written in terms of Riemann invariants using the technique proposed in [8, 9]. Notice that the Lagrangian coordinates in an one-dimensional case do not allow to allocate areas of subsonic and supersonic flows. We prove the global stability with respect both to small perturbations of initial and boundary conditions and monotonicity. Conditions only on the input data, allowing to guarantee the absence of shock waves have been obtained. In this paper we show that the conditions given for difference schemes coincide with the necessary conditions of absence of gradient catastrophe in the differential case.

The paper is organized as follows: Section 2 is devoted to the statement of the IBVP for a gas dynamics system. For the approximation of the system of equations in Riemann invariants, a linearized difference scheme is used. The conditions for the initial and boundary data that guarantee stability of the difference scheme are introduced. In Section 3, we obtain the necessary conditions of absence of shock waves for the differential problem. In Section 4, the stability of the proposed difference scheme is studied. In Section 5, we investigate the monotonicity. Section 6 is devoted to the investigation of the convergence of the difference scheme. Section 7 presents the statements of the numerical experiment, which illustrate the theoretical results.

2. Problem statement

In the domain \overline{Q} = \overline{Ω} × [0, +∞), Ω = \{x : 0 < x < l\}, let us consider the

IBVP for the system of equations of the gas dynamics written in the Lagrangian coordinates [11, 14]:

\frac{∂v}{∂t} − \frac{∂u}{∂x} = 0, \frac{∂u}{∂t} + \frac{∂p(v)}{∂x} = 0, (x, t) ∈ Q,

p(v) = K^{2}v^{−γ} , γ = 1 + 2ε, ε = const > 0, K = const > 0, (2.1)

v(x, 0) = v_{0}(x), u(x, 0) = u_{0}(x), 0 ⩽ x ⩽ l, (2.2)

\begin{matrix}u(0, t) + c_{1} ((v(0, t))^{−ε} − 1) = μ_{1}(t), t ⩾ 0, \\ − u(l, t) + c_{1} ((v(l, t))^{−ε} − 1) = μ_{2}(t), t ⩾ 0. \end{matrix} (2.3)

Here v = v(x, t), u = u(x, t), p = p(x, t) denote the specific volume, the velocity, and the pressure respectively, c_{1} = \frac{K\sqrt{γ}}{ε} . Hereinafter c_{i} denote positive constants.

Later we will prove that under some assumptions on input data the following estimate is valid: 0 < v_{min} ⩽ v(x, t) ⩽ v_{max} , (x, t) ∈ \overline{Q}. Then the eigenvalues of the matrix of the system (2.1) are real and distinct:

λ_{1} = +\sqrt{−\frac{∂p(v)}{∂v}} > 0, λ_{2} = −\sqrt{−\frac{∂p(v)}{∂v}} < 0.

Therefore, system (2.1) is hyperbolic and boundary conditions (2.3) are well posed [4]. Suppose that the initial data satisfy the inequalities

\begin{matrix}−|u_{0}(x)| + c_{1} (v_{0}(x)^{−ε} − 1) & ⩾ 0, \\ |u_{0}(x)| + c_{1} (v_{0}(x)^{−ε} − 1) & ⩽ c_{2}, \end{matrix} 0 ⩽ x ⩽ l, (2.4)

\begin{matrix}−c_{3} ⩽ −u′_{0}(x) − c_{1}εv_{0}(x)^{−(ε+1)} v′_{0}(x) & ⩽ 0, \\ 0 ⩽ u′_{0}(x) − c_{1}εv_{0}(x)^{−(ε+1)} v′_{0}(x) & ⩽ c_{3}, \end{matrix} 0 ⩽ x ⩽ l (2.5)

The differential problem (2.1)–(2.3) in Riemann invariants [11]

r = u + c_{1}(v^{−ε} − 1), s = −u + c_{1}(v^{−ε} − 1),

has the following form:

\frac{∂r}{∂t} + a(s + r)\frac{∂r}{∂x} = 0, \frac{∂s}{∂t} − a(s + r)\frac{∂s}{∂x} = 0, (x, t) ∈ Q, (2.6)

\begin{matrix}r(x, 0) = r_{0}(x) = u_{0}(x) + c_{1} ((v_{0}(x))^{−ε} − 1) , 0 ⩽ x ⩽ l, \\ s(x, 0) = s_{0}(x) = −u_{0}(x) + c_{1} ((v_{0}(x))^{−ε} − 1) , 0 ⩽ x ⩽ l, \end{matrix} (2.7)

r(0, t) = μ_{1}(t), s(l, t) = μ_{2}(t), t ⩾ 0, (2.8)

where a(s + r) = c_{4}(s + r + 2c_{1})^{1+ \frac{1}{ε}} , c_{4} = ε/2(2c_{1})^{\frac{1}{ε}}. For the convenience let us define b = s + r. Suppose that the boundary conditions satisfy the following inequalities:

(2.9)  0 ⩽ μ_{1}(t) ⩽ c_{2},  0 ⩽ μ_{2}(t) ⩽ c_{2}, t ⩾ 0,
(2.10)  −c_{3}a_{max} ⩽ μ′_{1}(t) ⩽ 0, c_{3}a_{max} ⩽ μ′_{2}(t) ⩽ 0, t ⩾ 0,

where a_{max} = c_{4}(2c_{2} + 2c_{1})^{1+ \frac{1}{ε}} .

Let us assume the fulfilment of the conjugation conditions:

r_{0}(0) = μ_{1}(0), s_{0}(l) = μ_{2}(0), (2.11)

μ′_{1}(0) + a(r_{0}(0) + s_{0}(0))r′_{0}(0) = 0, μ′_{2}(0) − a(r_{0}(l) + s_{0}(l))s′_{0}(l) = 0. (2.12)

Then from conditions (2.4)–(2.5) for the Riemann invariants we get the inequalities:

(2.13)  0 ⩽ r_{0}(x) ⩽ c_{2}, 0 ⩽ s_{0}(x) ⩽ c_{2}, 0 ⩽ x ⩽ l,
(2.14)  −c_{3} ⩽ s′_{0}(x) ⩽ 0 ⩽ r′_{0}(x) ⩽ c_{3}, 0 ⩽ x ⩽ l.

3. Formation of shock waves for the differential problem

It is well known [14] that in general case the smooth solutions of initial value problems (or initial-boundary value problems) for quasilinear hyperbolic systems exist only locally in time [17] and singularities will appear in finite time, even if the initial data (and the boundary data) are sufficiently smooth and small. The theory of gradient catastrophe is constructed mainly for the Cauchy problem [5, 18]:

\frac{∂v}{∂t} − \frac{∂u}{∂x} = 0, \frac{∂u}{∂t} + \frac{∂p(v)}{∂x} = 0, x ∈ ℝ, t > 0, (3.1)

v(x, 0) = v_{0}(x), u(x, 0) = u_{0}(x), x ∈ ℝ. (3.2)

In [5] the necessary condition for existence of global smooth solution have been derived.

Proposition 3.1. [5] Let initial data (3.2) satisfy the inequalities

−|u_{0}(x)| + c_{1} (v_{0}(x)^{−ε} − 1) ⩾ 0,
|u_{0}(x)| + c_{1} (v_{0}(x)^{−ε} − 1) ⩽ c_{2}, x ∈ ℝ.

Then the condition

−c_{3} ⩽ −u′_{0}(x) − c_{1}εv_{0}(x)^{−(ε+1)} v′_{0}(x) ⩽ 0
⩽ u′_{0}(x) − c_{1}εv_{0}(x)^{−(ε+1)} v′_{0}(x) ⩽ c_{3}, ∀x ∈ ℝ,

is necessary for global smooth solution of the problem (3.1)–(3.2).

Remark 3.1. Actually this condition is also sufficient (see [18]).

Therefore for the Cauchy problem there is a criterion of formation of shock waves. In the case of IBVP the situation is more difficult. The existing results have been proved with essential restrictions on input data. In [1], for example, the initial data are supposed to have small magnitudes. Below we will obtain necessary conditions for the differential case, using the technique from [5] and the idea of differentiation with respect to time variable from [4]. We will need the following generalizations of the results from [5].

Proposition 3.2. Let z(t) be the solution of the initial-value problem

(3.3)  \frac{dz}{dt} = a(t)z^{2}, z(t_{0}) = m,

in the interval (t_{0}, T). Suppose that the function a(t) satisfies the inequality

0 < A < a(t), t_{0} ⩽ t ⩽ T,

and that m is positive. Then T < t_{0} + (mA)^{−1} .

Proposition 3.3. Suppose that a(t) satisfies the inequality |a(t)| < B. Then initial-value problem (3.3) has a solution for |t − t_{0}| < |mB|^{−1} .

In the domain \overline{Q} let us consider the following IBVP for two quasilinear equations in Riemann invariants:

ŕ = r_{t} + λ(r, s)r_{x} = 0, (x, t) ∈ Q, (3.4) s̀ = s_{t} − λ(r, s)s_{x} = 0, (x, t) ∈ Q, (3.5)

\begin{aligned}r(x, 0) & = r_{0}(x), s(x, 0) = s_{0}(x), 0 ⩽ x ⩽ l, \\ r(0, t) & = μ_{1}(t), s(l, t) = μ_{2}(t), t ⩾ 0, \end{aligned} (3.6)

supposing that conditions (2.9), (2.11), (2.13) are fulfilled and

(3.7)  0 < λ(r_{0}, s_{0}) ⩽ c_{λ}, 0 < λ_{r}(r_{0}, s_{0}) ⩽ c_{λ_{r}} , 0 < λ_{s}(r_{0}, s_{0}) ⩽ c_{λ_{s}} .

Let us note that from equations (3.4)–(3.5) we see that the values of invariants transfer along characteristics, so

0 < λ(r, s) ⩽ c_{λ}, 0 < λ_{r}(r, s) ⩽ c_{λ_{r}} , 0 < λ_{s}(r, s) ⩽ c_{λ_{s}} ,
∀ r(x, t), s(x, t), 0 ⩽ x ⩽ l, t ⩾ 0.

Remark 3.2. From conditions (2.13) we obtain that λ(r_{0}, s_{0}) > 0, so λ(r, s) > 0 for all r(x, t), s(x, t), 0 ⩽ x ⩽ l, t ⩾ 0. So the fulfillment of conditions (2.9), (2.11), (2.13) is sufficient for hyperbolicity of system (2.1).

Let us differentiate the equation (3.4) with respect to t:

(3.8)  r_{t t} + λr_{x t} + λ_{r}r_{t}r_{x} + λ_{s}s_{t}r_{x} = 0.

From (3.4) and (3.5) we get that 0 = s̀ = ś − 2λs_{x}, so s_{x} = \frac{ś}{2λ}. Therefore

(3.9)  s_{t} = \frac{1}{2}ś.

Abbreviating ω = r_{t} and substituting (3.9) in (3.8), we get the following equation:

(3.10)  ώ − \frac{λ_{r}}{λ} ω^{2} − \frac{λ_{s}}{2λ} śω = 0.

Denote by k = k(r, s) a function satisfying k_{s} = −\frac{λ_{s}}{2λ} . Using (3.4), we have

ḱ = k_{r} ŕ + k_{s} ś = −\frac{λ_{s}}{2λ} ś.

By substituting this equation into (3.10) one gets

ώ − \frac{λ_{r}}{λ} ω^{2} + ḱω = 0.

Multiplying by e^{k} and abbreviating z = e^{k}ω, we get

ź = e^{−k} \frac{λ_{r}}{λ} z^{2},

the first equation of (3.3) with a = e^{−k} λ_{r} /λ and initial condition

m = \{ \begin{matrix}− min_{0⩽x⩽l} e^{k(r_{0} ,s_{0} )} λ(r_{0}, s_{0})r′_{0}(x), & t_{0} = 0, 0 ⩽ x ⩽ l; \\ e^{k(r_{0} ,s_{0} )} μ′_{1}(t_{0}), & t_{0} ⩾ 0. \end{matrix}

Similarly for equation (3.5) we get equation (3.3) with a = e^{−k} λ_{s} /λ and the initial condition

m = \{ \begin{matrix}max_{0⩽x⩽l} e^{k(r_{0} ,s_{0} )} λ(r_{0}, s_{0})s′_{0}(x), & t_{0} = 0, 0 ⩽ x ⩽ l; \\ e^{k(r_{0} ,s_{0} )} μ′_{2}(t_{0}), & t_{0} ⩾ 0. \end{matrix}

Using propositions 3.2 and 3.3 and conditions (3.7) we obtain the corresponding conditions for the signs of derivations of the input data and requirement of their boundness. Finally, we derive

Theorem 3.1. Let inequalities (2.9), (2.11) and (2.13) be true. Then the conditions (2.5), (2.10) and (2.12) are necessary for existence of global smooth solutions of problem (2.1)–(2.3).

Moreover, using results from [10] we can formulate the following unique solvability theorem for the problem already written in Riemann invariants.

Theorem 3.2. If the above-formulated assumptions (2.9)–(2.14) are valid, then for an arbitrary given t_{0} > 0, there exists a unique continuous generalized solution of problem (2.6)–(2.8) in the class of functions with weak discontinuities and this solution transfers along the characteristics from the initial and boundary data.

Remark 3.3. Unfortunately the previous result is not valid for the boundary conditions not stated in Riemann invariants.

4. Stability of difference scheme

In the domain \overline{Q} we introduce a uniform grid \overline{ω}_{hτ} : \overline{ω}_{hτ} = \overline{ω}_{h} × \overline{ω}_{τ} : \overline{ω}_{h} = \{x_{i} = ih, i = \overline{0, N}, h = \frac{l}{N} \} , \overline{ω}_{τ} = \{t_{n} = nτ, n ∈ ℕ_{0}\} . On the grid \overline{ω}_{hτ} we approximate the differential problem in the Riemann invariants (2.6)–(2.8) by the linearized difference scheme

\begin{matrix}r_{ht,i} + a_{h,i} r̂_{hx̄,i} = 0, i = \overline{1, N}, \\ s_{ht,i} − a_{h,i} ŝ_{hx,i} = 0, i = \overline{0, N −1} , \end{matrix} (4.1)

\begin{matrix}r_{h,i}^{0} = r_{0,i} = u_{0,i} + c_{1} ((v_{0,i})^{−ε} − 1) , i = \overline{0, N}, \\ s_{h,i}^{0} = s_{0,i} = −u_{0,i} + c_{1} ((v_{0,i})^{−ε} − 1) , i = \overline{0, N}, \end{matrix} (4.2)

r̂_{h,0} = μ_{1}^{n+1} , ŝ_{h,N} = μ_{2}^{n+1} , n ∈ ℕ_{0}. (4.3)

Hereinafter we use standard notations of the difference schemes theory [15, 16]:

y = y_{i}^{n} = y(x_{i}, t_{n}), ŷ = y_{i}^{n+1} = y(x_{i}, t_{n+1} ),

y_{x̄,i} = \frac{y_{i} − y_{i−1}}{h} , y_{x,i} = \frac{y_{i+1} − y_{i}}{h} , y_{t,i} = \frac{ŷ_{i} − y_{i}}{τ} .

When investigating the stability of the difference problem (2.6)–(2.8), we will use the following canonical form [7] of the two-point difference scheme for the initial value problem

C_{i}y_{i} = A_{i}y_{i−1} + F_{i}, i = \overline{i_{0} + 1, i_{N}} , y_{i_{0}} = μ_{1}, (4.4)

C_{i}y_{i} = B_{i}y_{i+1} + F_{i}, i = \overline{i_{0}, i_{N} −1} , y_{i_{N}} = μ_{2}. (4.5)

We will need the following results [7].

Lemma 4.1. Let the conditions A_{i} ⩾ 0, D_{i} = C_{i} − A_{i} > 0, i = \overline{i_{0} + 1, i_{N}} , be met. Then for the solution of problem (4.4) the estimate

max_{i_{0} ⩽i⩽i_{N}} |y_{i}| ⩽ max\{ |μ_{1}|, max_{i_{0} <i⩽i_{N}} \frac{|F_{i}|}{D_{i}}\}

is valid.

Lemma 4.2. Let the conditions B_{i} ⩾ 0, D_{i} = C_{i} − B_{i} > 0, i = \overline{i_{0}, i_{N} −1} , be met. Then for the solution of problem (4.5) the estimate

max_{i_{0} ⩽i⩽i_{N}} |y_{i}| ⩽ max\{ |μ_{2}|, max_{i_{0} ⩽i<i_{N}} \frac{|F_{i}|}{D_{i}}\}

is valid.

Lemma 4.3. Let conditions A_{i} ⩾ 0, C_{i} > 0, i = \overline{i_{0} + 1, i_{N}} , be met. If F_{i} ⩾ 0, i = \overline{i_{0} + 1, i_{N}} , μ_{1} ⩾ 0 (F_{i} ⩽ 0, i = \overline{i_{0} + 1, i_{N}} , μ_{1} ⩽ 0), then for the solution of problem (4.4) the estimate y_{i} ⩾ 0 (y_{i} ⩽ 0), i = \overline{i_{0} + 1, i_{N}} , is valid.

Lemma 4.4. Let conditions B_{i} ⩾ 0, C_{i} > 0, i = \overline{i_{0}, i_{N} −1} , be met. If F_{i} ⩾ 0, i = \overline{i_{0}, i_{N} −1} , μ_{2} ⩾ 0 (F_{i} ⩽ 0, i = \overline{i_{0}, i_{N} −1} , y_{i_{N}} = μ_{2} ⩽ 0), then for the solution of problem (4.5) the estimate y_{i} ⩾ 0 (y_{i} ⩽ 0), i = \overline{i_{0}, i_{N} −1} , is valid.

Along with (4.1)–(4.3) we consider the perturbed problem:

\begin{matrix}r̃_{ht,i} + ã_{h,i} ˆ r̃_{hx̄,i} = 0, i = \overline{1, N}, \\ s̃_{ht,i} − ã_{h,i} ˆ s̃_{hx,i} = 0, i = \overline{0, N −1} , \end{matrix} (4.6)

\begin{matrix}r̃_{h,i}^{0} = r̃_{0,i} = ũ_{0,i} + c_{1} ((ṽ_{0,i})^{−ε} − 1) , i = \overline{0, N}, \\ s̃_{h,i}^{0} = s̃_{0,i} = −ũ_{0,i} + c_{1} ((ṽ_{0,i})^{−ε} − 1) , i = \overline{0, N}, \end{matrix} (4.7)

ˆr̃_{h,0} = μ̃_{1} ^{n+1} , ˆ s̃_{h,N} = μ̃_{2} ^{n+1} , n ∈ ℕ_{0}. (4.8)

Let the following inequalities analogous to (2.9), (2.11), (2.13) for the perturbed data (4.7) be satisfied:

0 ⩽ r̃_{0}(x) ⩽ c_{2}, 0 ⩽ s̃_{0}(x) ⩽ c_{2}, 0 ⩽ x ⩽ l, (4.9)

0 ⩽ μ̃_{1} (t) ⩽ c_{2}, 0 ⩽ μ̃_{2} (t) ⩽ c_{2}, t ⩾ 0, (4.10)

r̃_{0}(0) = μ̃_{1} (0), s̃_{0}(l) = μ̃_{2} (0). (4.11)

Subtracting from difference equations (4.6)–(4.8) the corresponding equations (4.1)–(4.3), we come to the problem for the perturbations δr_{i} = r̃_{h,i} − r_{h,i} , δs_{i} = s̃_{h,i} − s_{h,i} :

\begin{matrix}δr_{t,i} + ã_{h,i} δr̂_{x̄,i} + δa_{i} r̂_{x̄,i} = 0, i = \overline{1, N}, \\ δs_{t,i} − ã_{h,i} δŝ_{x,i} − δa_{i} ŝ_{x,i} = 0, i = \overline{0, N −1} , \end{matrix} (4.12)

\begin{matrix}δr_{i}^{0} = δr_{0,i} = δu_{0,i} + c_{1}δ (v_{0,i}^{−ε} ), i = \overline{0, N}, \\ δs_{i}^{0} = δs_{0,i} = −δu_{0,i} + c_{1}δ (v_{0,i}^{−ε} ), i = \overline{0, N}, \end{matrix} (4.13)

δr̂_{0} = δμ_{1}(t̂), δŝ_{N} = δμ_{2}(t̂), n ∈ ℕ_{0}, (4.14) where δa_{i} = ã_{h,i} − a_{h,i} . Hereinafter we will use the following grid norms:

∥y_{h}∥_{C_{h}^{+}} = max_{1⩽i⩽N} |y_{h,i} |, ∥y_{h}∥_{C_{h}^{−}} = max_{0⩽i⩽N−1} |y_{h,i} |,
∥y_{h}∥C̄_{h} = max_{0⩽i⩽N} |y_{h,i} |, ∥y_{h}^{n}∥_{C_{τ}} = max_{0⩽j⩽n} |y_{h}^{j}|.

From (4.12)–(4.13) it follows that we need to analyze stability before studying the behavior of all difference derivatives in the nonlinear terms of the difference equations [9, 12, 13]. Let us prove

Theorem 4.1. 1. Let conditions (2.9), (2.11), (2.13) be met. Then the following estimates:

∥r^{n}∥C̄_{h} ⩽ ∥r^{0}∥C̄_{h} ⩽ c_{2}, ∥s^{n}∥C̄_{h} ⩽ ∥s^{0}∥C̄_{h} ⩽ c_{2}, ∥a^{n}∥C̄_{h} ⩽ a_{max} = c_{4} (2c_{2} + 2c_{1})^{1+ \frac{1}{ε}} ,

b_{i}^{n} ⩾ 0, a_{i}^{n} ⩾ a_{min} = c_{4}(2c_{1})^{1+ \frac{1}{ε}} , i = \overline{0, N}, n ∈ ℕ_{0},

are valid.

2. Let conditions (2.9), (4.9), (4.11) be met. Then the following estimates:

∥r̃^{n} ∥C̄_{h} ⩽ ∥r̃^{0}∥C̄_{h} ⩽ c_{2}, ∥s̃^{n} ∥C̄_{h} ⩽ ∥s̃^{0}∥C̄_{h} ⩽ c_{2}, ∥ã^{n}∥C̄_{h} ⩽ a_{max} = c_{4}(2c_{2} + 2c_{1})^{1+ \frac{1}{ε}} ,

b̃_{i}^{n} ⩾ 0, ã_{i}^{n} ⩾ a_{min} = c_{4}(2c_{1})^{1+ \frac{1}{ε}} , i = \overline{0, N}, n ∈ ℕ_{0},

are valid.

Proof. We will prove only the first part of the theorem, because the second one can be proved similarly. Let us rewrite the difference scheme (4.1)–(4.3) in the canonical form:

C_{i}^{n}r_{h,i}^{n+1} = A_{i}^{n}r_{h,i−1}^{n+1} + F_{r,i}^{n} , i = \overline{1, N}, r_{h,0}^{n+1} = μ_{1}^{n+1} ,

C_{i}^{n}s_{h,i}^{n+1} = A_{i}^{n}s_{h,i+1}^{n+1} + F_{s,i}^{n} , i = \overline{0, N −1} , s_{h,N}^{n+1} = μ_{2}^{n+1} ,

A_{i}^{n} = \frac{a_{i}^{n}τ}{h}, C_{i}^{n} = 1 + A_{i}^{n}, F_{r,i}^{n} = r_{h,i}^{n} , F_{s,i}^{n} = s_{h,i}^{n} .

We will use the method of mathematical induction. From conditions (2.9), (2.13) we get ∥b^{0}∥C̄_{h} ⩽ 2c_{2}, ∥a^{0}∥C̄_{h} ⩽ a_{max} , b_{i}^{0} ⩾ 0, a_{i}^{0} ⩾ a_{min} , i = \overline{0, N}. On the first layer the coefficients satisfy the conditions: A_{i}^{0} > 0, i = \overline{0, N}. By lemmas 4.1 and 4.2 we obtain that

∥r^{1}∥C̄_{h} ⩽ max{|μ_{1}^{1}|, ∥r^{0}∥_{C_{h}^{+}} } ⩽ max{|μ_{1}^{1}|, ∥r^{0}∥C̄_{h} } ⩽ c_{2},
∥s^{1}∥C̄_{h} ⩽ max{|μ_{2}^{1}|, ∥s^{0}∥_{C_{h}^{−}} } ⩽ max{|μ_{2}^{1}|, ∥s^{0}∥C̄_{h} } ⩽ c_{2}.

Therefore, ∥b^{1}∥C̄_{h} ⩽ 2c_{2} è ∥a^{1}∥C̄_{h} ⩽ a_{max} . Using lemmas 4.3 and 4.4 we prove that r_{i}^{1} ⩾ 0, i = \overline{1, N}, r_{0}^{1} = μ_{1}^{1} ⩾ 0, s_{i}^{1} ⩾ 0, i = \overline{0, N −1} , s_{N}^{1} = μ_{2}^{1} ⩾ 0. So b_{i}^{1} ⩾ 0, a_{i}^{1} ⩾ a_{min} , i = \overline{0, N}. The proof can be completed by using the method of mathematical induction. □

Remark 4.1. In the previous theorem existence and boundness of the solution of problem (4.1)–(4.3) are proved.

Let us rewrite the equations (4.1) in the following form:

\begin{matrix}r_{h t,i} + a_{h,i} r_{h tx̄,i} = −a_{h,i} r_{h x̄,i} = \frac{a_{h,i}}{ǎ_{h,i}} r_{h t̄,i} , \\ s_{h t,i} − a_{h,i} s_{h tx,i} = a_{h,i} s_{h x,i} = \frac{a_{h,i}}{ǎ_{h,i}} s_{h t̄,i} . \end{matrix}

Abbreviating g = r_{h t}, f = s_{h t}, we obtain the difference problem:

\begin{matrix}C_{i}^{n+1} g_{i}^{n+1} = A_{i}^{n+1} g_{i−1}^{n+1} + F_{g,i}^{n+1} , i = \overline{1, N}, g_{0}^{n+1} = μ_{1 t}^{n+1} , n ∈ ℕ_{0}, \\ C_{i}^{n+1} f_{i}^{n+1} = A_{i}^{n+1} f_{i+1}^{n+1} + F_{f,i}^{n+1} , i = \overline{0, N −1} , f_{N}^{n+1} = μ_{2 t}^{n+1} , n ∈ ℕ_{0}, \\ A_{i}^{n+1} = \frac{a_{i}^{n+1} τ}{h}, F_{g,i}^{n+1} = \frac{a_{i}^{n+1}}{a_{i}^{n}} g_{i}^{n}, F_{f,i}^{n+1} = \frac{a_{i}^{n+1}}{a_{i}^{n}} f_{i}^{n}, C_{i}^{n+1} = 1 + A_{i}^{n+1} , \end{matrix} ; (4.15)

with the initial conditions:

\begin{matrix}C_{i}^{0}g_{i}^{0} = A_{i}^{0}g_{i−1}^{0} + F_{g,i}^{0} , i = \overline{1, N}, g_{0}^{0} = μ_{1 t}^{0}, \\ C_{i}^{0}f_{i}^{0} = A_{i}^{0}f_{i+1}^{0} + F_{f,i}^{0} , i = \overline{0, N −1} , f_{N}^{0} = μ_{2 t}^{0}, \\ A_{i}^{0} = \frac{a_{i}^{0}τ}{h}, F_{g,i}^{0} = −a_{i}^{0}r_{0x̄,i} , F_{f,i}^{0} = a_{i}^{0}s_{0x,i} , C_{i}^{0} = 1 + A_{i}^{0}. \end{matrix} ; (4.16)

Let us take into consideration the norm Q^{n} = c_{5} (∥g^{n}∥C̄_{h} + ∥f^{n}∥ C̄_{h} ) , where c_{5} = c_{4} (2c_{2} + 2c_{1}) ^{\frac{1}{ε}} (1 + \frac{1}{ε}) . Later we will use

Lemma 4.5. Let conditions (2.9), (2.11), (2.13) be met. Then

∥δa^{n}∥C̄_{h} ⩽ c_{5} (∥δs^{n}∥C̄_{h} + ∥δr^{n}∥C̄_{h} ) , n ∈ ℕ_{0}. Proof. Using the mean value theorem and Theorem 4.1, we obtain

ã_{i}^{n} − a_{i}^{n} = c_{4}( 1+ \frac{1}{ε}) (θ(s̃_{i} +r̃_{i})+(1 − θ)(s_{i}+r_{i})+2c_{1}) ^{\frac{1}{ε}} (δs_{i}+δr_{i}), θ ∈ (0, 1).  □

Corollary 4.1. Let conditions (2.9), (2.11), (2.13) be met. Then ∥a_{t}^{n}∥C̄_{h} ⩽ Q^{n}, n ∈ ℕ_{0}.

Theorem 4.2. Let conditions (2.9)–(2.14) be met. Then for the difference derivatives the estimate

c_{5}∥g^{n+1} ∥C̄_{h} ⩽ c_{6}a_{min} , c_{5}∥f^{n+1} ∥C̄_{h} ⩽ c_{6}a_{min} , c_{6} = \frac{c_{3}c_{5}a_{max}}{a_{min}} , n ∈ ℕ_{0},

is valid.

Proof. Using the positiveness of coefficients (4.15)–(4.16), conditions (2.9)– (2.14) and lemmas 4.3 and 4.4 we get that g_{i}^{n+1} ⩽ 0, f_{i}^{n+1} ⩽ 0, i = \overline{0, N}, n ∈ ℕ_{0}. Therefore ∥a^{n+1} /a^{n} ∥C̄_{h} ⩽ 1, n ∈ ℕ_{0}. Then, taking into account that coefficients (4.15)–(4.16) satisfy the conditions of lemmas 4.1 and 4.2, we obtain:

∥g^{n+1} ∥C̄_{h} ⩽ max\{|μ_{1 t}^{n+1} |, ∥ \frac{a^{n+1}}{a^{n}} ∥_{C_{h}^{+}} · ∥g^{n}∥_{C_{h}^{+}} \} ⩽ max\{c_{3}a_{max} , ∥g^{n}∥C̄_{h} \}
⩽ … ⩽ max{c_{3}a_{max} , ∥g^{1}∥C̄_{h} } ⩽ max{c_{3}a_{max} , ∥a^{0}∥C̄_{h} ∥r_{0x̄}∥_{C_{h}^{+}} } ⩽ c_{3}a_{max} .

Similarly: ∥f^{n+1} ∥C̄_{h} ⩽ ∥a^{0}∥C̄_{h} ∥s_{0x}∥_{C_{h}^{−}} ⩽ c_{3}a_{max} . □

Let us rewrite scheme (4.12)–(4.14) in the canonical form:

\begin{matrix}C_{i}^{n}δr_{h,i}^{n+1} = A_{i}^{n}δr_{h,i−1}^{n+1} + F_{r,i}^{n} , i = \overline{1, N}, δr_{h,0}^{n+1} = δμ_{1}(t̂), \\ C_{i}^{n}δs_{h,i}^{n+1} = A_{i}^{n}δs_{h,i+1}^{n+1} + F_{s,i}^{n} , i = \overline{0, N −1} , δs_{h,N}^{n+1} = δμ_{2}(t̂), \\ A_{i}^{n} = \frac{ã_{i}^{n}τ}{h}, C_{i}^{n} = 1 + A_{i}^{n}, F_{r,i}^{n} = δr_{h,i}^{n} + τδa_{i}^{n} \frac{g_{i}^{n}}{a_{i}^{n}} , F_{s,i}^{n} = δs_{h,i}^{n} + τδa_{i}^{n} \frac{f_{i}^{n}}{a_{i}^{n}} . \end{matrix} ; (4.17)

Let us prove

Theorem 4.3. Let conditions (2.9)–(2.14) and (4.9)–(4.11) be met. Then difference scheme (4.1)–(4.3) is stable with respect to the initial and boundary data for τ ⩽ c_{6}^{−1} and for its solution the estimate max{∥δs^{n+1} ∥C̄_{h} , ∥δr^{n+1} ∥C̄_{h} } ⩽ max{∥δμ_{1}^{n+1} ∥_{C_{τ}} , ∥δμ_{2}^{n+1} ∥_{C_{τ}} , ∥δs^{0}∥C̄_{h} , ∥δr^{0}∥C̄_{h} } is valid.

Proof. Taking into account the above estimates for the derivatives, the proof of Lemma 4.5 and the theorem conditions, we obtain for θ ∈ (0, 1):

1 + τ \frac{g_{i}^{n}}{a_{i}^{n}} c_{4}( 1 + \frac{1}{ε}) (θ(s̃_{i}^{n} + r̃_{i}^{n}) + (1 − θ)(s_{i}^{n} + r_{i}^{n}) + 2c_{1})^{\frac{1}{ε}} ⩾ 1 − c_{6}τ ⩾ 0.

Similarly

1 + τ \frac{f_{i}^{n}}{a_{i}^{n}} c_{4}( 1 + \frac{1}{ε}) (θ(s̃_{i}^{n} + r̃_{i}^{n}) + (1 − θ)(s_{i}^{n} + r_{i}^{n}) + 2c_{1}) ^{\frac{1}{ε}} ⩾ 1 − c_{6}τ ⩾ 0.

Let α = c_{6}τ, β = 1 − c_{6}τ. We will use the method of mathematical induction. On the zero layer we get the following estimates for the perturbations δr, δs:

max{∥δs^{0}∥C̄_{h} , ∥δr^{0}∥C̄_{h} } ⩽ max{∥δμ_{1}^{n+1} ∥_{C_{τ}} , ∥δμ_{2}^{n+1} ∥_{C_{τ}} , ∥δs^{0}∥C̄_{h} , ∥δr^{0}∥C̄_{h} } . Therefore on the first layer we obtain:

\begin{matrix}max & \{∥δs^{1}∥C̄_{h} , ∥δr^{1}∥C̄_{h} \} \\ ⩽ max\{∥δμ_{1}^{n+1} ∥_{C_{τ}} , ∥δμ_{2}^{n+1} ∥_{C_{τ}} , β∥δr^{0}∥C̄_{h} + α∥δs^{0}∥C̄_{h} , β∥δs^{0}∥C̄_{h} + α∥δr^{0}∥C̄_{h} \} . \end{matrix}

Using the inequalities

β∥δr^{0}∥C̄_{h} + α∥δs^{0}∥C̄_{h} ⩽ β max{∥δr^{0}∥C̄_{h} , ∥δs^{0}∥C̄_{h} } + α max{∥δr^{0}∥C̄_{h} , ∥δs^{0}∥C̄_{h} }
= (α + β) max{∥δr^{0}∥C̄_{h} , ∥δs^{0}∥C̄_{h} } = max{∥δr^{0}∥C̄_{h} , ∥δs^{0}∥C̄_{h} } ,
β∥δs^{0}∥C̄_{h} + α∥δr^{0}∥C̄_{h} ⩽ max{∥δr^{0}∥C̄_{h} , ∥δs^{0}∥C̄_{h} }

(4.18)

it is easy to prove that

max{∥δs^{1}∥C̄_{h} , ∥δr^{1}∥C̄_{h} } ⩽ max{∥δμ_{1}^{n+1} ∥_{C_{τ}} , ∥δμ_{2}^{n+1} ∥_{C_{τ}} , ∥δr^{0}∥C̄_{h} , ∥δs^{0}∥C̄_{h} } . Finally by induction we show that max{∥δs^{n+1} ∥C̄_{h} , ∥δr^{n+1} ∥C̄_{h} } ⩽ max{∥δμ_{1}^{n+1} ∥_{C_{τ}} , ∥δμ_{2}^{n+1} ∥_{C_{τ}} , ∥δs^{0}∥C̄_{h} , ∥δr^{0}∥C̄_{h} }.

□

Remark 4.2. Much in the same manner, the uniqueness of the solution of difference scheme (4.1)–(4.3) is proved. The proof is carried out by contradiction. We assume that there exist two solutions (r_{1}, s_{1}), (r_{2}, s_{2}) that satisfy the difference scheme with the same initial and boundary conditions. Then for the differences R = r_{2}−r_{1}, S = s_{2}−s_{1} we get the following difference equations with homogeneous initial and boundary conditions:

R_{t,i} + a_{2 i} R̂_{x̄,i} + A_{i} r̂_{x̄,i} = 0, i = \overline{1, N},
S_{t,i} − a_{2 i} Ŝ_{x,i} − A_{i} ŝ_{x,i} = 0, i = \overline{0, N −1} ,
R_{i}^{0} = 0, S_{i}^{0} = 0, i = \overline{0, N},
R̂_{0} = 0, Ŝ_{N} = 0, n = \overline{0, N_{0}−1} ,

where A_{i} = a_{2 i} − a_{1 i}. Using the introduced technique, we get the estimate max{∥R^{n+1} ∥C̄_{h} , ∥S^{n+1} ∥C̄_{h} } ⩽ 0. Therefore, the difference scheme has a unique solution.

In the theorem above the stability in the case when a shock wave is not generated is proved. To investigate the stability in the case when a shock wave arises we must prove the following statement.

Theorem 4.4. Let conditions (2.11), (2.13), (4.9) and

(4.19)  \begin{matrix}μ_{1}(t) = μ_{1}, μ_{2}(t) = μ_{2}, μ_{1}, μ_{2} = const, \\ |s′_{0}(x)| ⩽ c_{3}, |r′_{0}(x)| ⩽ c_{3}, 0 ⩽ x ⩽ l, \end{matrix}

be met. Then difference scheme (4.1)–(4.3) is stable with respect to the initial and boundary data for τ ⩽ c_{7}^{−1} and for the solution the following estimate

max{∥δs^{n+1} ∥C̄_{h} , ∥δr^{n+1} ∥C̄_{h} } ⩽ max{∥δμ_{1}^{n+1} ∥_{C_{τ}} , ∥δμ_{2}^{n+1} ∥_{C_{τ}} , ∥δs^{0}∥C̄_{h} , ∥δr^{0}∥C̄_{h} }

is valid for t_{n+1} ⩽ T < t_{0}, t_{0} = a_{min} /c_{5}a_{max} (∥r0 ; x̄∥_{C_{h}^{+}} + ∥s_{x}^{0}∥_{C_{h}^{−}} ) .

Proof. In this case the condition ∥a^{n+1} /a^{n} ∥C̄_{h} ⩽ 1, n ∈ ℕ_{0} is not met. Let us estimate this ratio for n = 0 assuming that 1 − t_{n+1} Q^{0}/a_{min} > 0. Using Corollary 4.1 we obtain 1 − τ∥a_{t}^{0}∥C̄_{h} /a_{min} ⩾ 1 − t_{n+1} Q^{0}/a_{min} > 0. Therefore

∥ \frac{a^{1}}{a^{0}} ∥C̄_{h} = \frac{1}{∥1 − τa_{t}^{0}/a^{0}∥ C̄_{h}} ⩽ \frac{1}{1 − τ∥a_{t}^{0}∥C̄_{h} /a_{min}} ⩽ \frac{1}{1 − τQ^{0}/a_{min}} .

Let us note that in the case of constant boundary values their difference derivatives with respect to time variable are equal to zero. So using lemmas 4.1 and 4.2 we get that Q^{1} ⩽ \frac{Q^{0}}{1−τQ^{0}/a_{min}} . Then

1 − \frac{τ∥a_{t}^{1}∥C̄_{h}}{a_{min}} ⩾ 1 − \frac{τQ^{1}}{a_{min}} ⩾ 1 − \frac{τQ^{0}/a_{min}}{1 − τQ^{0}/a_{min}} = \frac{1 − t_{2}Q^{0}/a_{min}}{1 − t_{1}Q^{0}/a_{min}} > 0.

By mathematical induction it is possible to show that Q^{n+1} ⩽ \frac{Q^{n}}{1−τQ^{n} /a_{min}} , n ∈ ℕ_{0} and using the discrete analogue of the Bihary inequality [8] we obtain the estimate

Q^{n+1} ⩽ \frac{Q^{0}}{1 − t_{n+1} Q^{0}/a_{min}} .

Since for (4.16) the conditions of lemmas 4.1 and 4.2 are fulfilled, then we get

Q^{0} ⩽ c_{5}a_{max} (∥r0 ; x̄∥_{C_{h}^{+}} + ∥s_{x}^{0}∥_{C_{h}^{−}} ). Therefore

Q^{n+1} ⩽ \frac{c_{5}a_{max} (∥r0 ; x̄∥_{C_{h}^{+}} + ∥s_{x}^{0}∥_{C_{h}^{−}} )}{1−t_{n+1} c_{5}a_{max} a_{min}^{−1} (∥r0 ; x̄∥_{C_{h}^{+}} +∥s_{x}^{0}∥_{C_{h}^{−}} )} ⩽ c_{7}a_{min} , c_{7} = \frac{2c_{5}a_{max} a_{min}^{−1} c_{3}}{1−2T c_{5}a_{max} a_{min}^{−1} c_{3}} .

Therefore, the derivatives are bounded in the norm Q for the time instant T < t_{0},

t_{0} = a_{min} /c_{5}a_{max} (∥r0 ; x̄∥_{C_{h}^{+}} + ∥s_{x}^{0}∥_{C_{h}^{−}} ) . Thus

c_{5}∥g^{n+1} ∥C̄_{h} ⩽ c_{7}a_{min} , c_{5}∥f^{n+1} ∥C̄_{h} ⩽ c_{7}a_{min} , c_{7} = \frac{c_{3}c_{5}a_{max}}{a_{min}} , n ∈ ℕ_{0}.

The end of proof is similar to Theorem 4.3, where estimates (4.18) with α = c_{7}τ, β = 1 − c_{7}τ. □

5. Investigation of monotonicity

The concept of monotonicity of difference schemes is important in the theory of numerical methods as it means absence of nonphysical oscillations in numerical calculations. Monotonicity for the linear difference scheme follows from the requirement of its coefficients positivity [3] or fulfillment of a grid maximum principle [15]. The most natural definition is given in [2]. Let us consider an abstract problem

(5.1)  L_{h}y = φ,

where L_{h}–the nonlinear difference operator defining the structure of the difference scheme, y–the difference solution, φ–problem input data. Perturbing the input data φ in (5.1), we get the equation for a perturbed solution

(5.2)  L_{h} ỹ = φ̃.

Subtracting (5.1) from (5.2), we obtain the following problem for δy = ỹ − y:

(5.3)  L̃_{h}(δy, y, ỹ) = δφ,

where δφ = φ̃ − φ.

Definition 5.1. Difference scheme (5.1) is called monotonic, if from the condition δφ ⩾ 0 (δφ ⩽ 0) follows the inequality δy ⩾ 0 (δy ⩽ 0).

Therefore the problem for investigation monotonicity coincides with the problem for investigation stability. Thus, investigation of monotonicity in a nonlinear case is reduced to the requirement of positiveness of coefficients or to the fulfillment of a grid maximum principle, but already for the problem for perturbations (5.3), which does not coincide with problem (5.1) or (5.2). In this case the requirement of positivity of coefficients leads to requirement of the monotonicity of difference derivatives.

Let us consider problem (4.12)–(4.14). The first theorem is proved in the case of global stability.

Theorem 5.1. Let conditions (2.9)–(2.14) and (4.9)–(4.11) be met and δr_{i}^{0} ⩾ 0, δs_{i}^{0} ⩾ 0, i = \overline{0, N}, δμ_{1}^{j} ⩾ 0, δμ_{2}^{j} ⩾ 0, j = \overline{0, n+1} (δs_{i}^{0} ⩽ 0, δr_{i}^{0} ⩽ 0, i = \overline{0, N}, δμ_{1}^{j} ⩾ 0, δμ_{2}^{j} ⩾ 0, j = \overline{0, n+1} ). Then for τ ⩽ c_{6}^{−1} the following estimates are true: δr_{i}^{n+1} ⩾ 0, δs_{i}^{n+1} ⩾ 0, i = \overline{0, N}, n ∈ ℕ_{0} (δs_{i}^{n+1} ⩽ 0, δr_{i}^{n+1} ⩽ 0, i = \overline{0, N}, n ∈ ℕ_{0}), i.e., difference scheme (4.1)–(4.3) is monotonic.

Proof. From (2.10) and (2.14) we obtain

r_{x̄,i}^{0} = \frac{1}{h} ∫_{x_{i−1}}^{x_{i}} r′_{0}(x)dx ⩾ 0, i = \overline{1, N}, s_{x,i}^{0} = \frac{1}{h} ∫_{x_{i}}^{x_{i+1}} s′_{0}(x)dx ⩽ 0, i = \overline{0, N −1} ,

g_{0}^{n} = \frac{1}{τ} ∫_{t_{n}}^{t_{n+1}} μ′_{1}(t)dt ⩽ 0, f_{N}^{n} = \frac{1}{τ} ∫_{t_{n}}^{t_{n+1}} μ′_{2}(t)dt ⩽ 0, n ∈ ℕ_{0}.

Let us consider the first case. Let δr_{i}^{0} ⩾ 0, δs_{i}^{0} ⩾ 0, i = \overline{0, N}. We have already proved in theorem 4.3 that all the coefficients of (4.17) are nonnegative. Since the boundary values perturbation is also nonnegative in this case we obtain the necessary result. The second case is considered similarly. □

The next theorem is proved for the case when a shock wave arises.

Theorem 5.2. Let conditions (2.11), (2.13), (4.9) and (4.19) be met. Assume δr_{i}^{0} ⩾ 0, δs_{i}^{0} ⩾ 0, i = \overline{0, N}, δμ_{1}^{j} ⩾ 0, δμ_{2}^{j} ⩾ 0, j = \overline{0, n+1} (δs_{i}^{0} ⩽ 0, δr_{i}^{0} ⩽ 0, i = \overline{0, N}, δμ_{1}^{j} ⩾ 0, δμ_{2}^{j} ⩾ 0, j = \overline{0, n+1} ). Then for τ ⩽ c_{7}^{−1} and for t_{n+1} ⩽ T < t_{0}, t_{0} = a_{min} /c_{5}a_{max} (∥r0 ; x̄∥_{C_{h}^{+}} + ∥s_{x}^{0}∥_{C_{h}^{−}} ) the following estimates are true: δr_{i}^{n+1} ⩾ 0,

δs_{i}^{n+1} ⩾ 0, i = \overline{0, N} (δs_{i}^{n+1} ⩽ 0, δr_{i}^{n+1} ⩽ 0, i = \overline{0, N}), i.e., difference scheme (4.1)–(4.3) is monotonic unless a smooth solution exists.

Proof. The proof is similar to the previous one using the proof of Theorem 4.4. □

6. Investigation of convergence

Let us note that since a priori estimates of the stability have been obtained, the investigation of the convergence of the difference solution to the differential one becomes much simpler since the problem for the error of the method can be written in a linear form.

Theorem 6.1. Suppose that there exists a solution r(x, t), s(x, t) ∈ C^{2,2}(\overline{Q}) of problem (2.6)–(2.8) and r_{0}(x) ⩾ 0, s_{0}(x) ⩾ 0, 0 ⩽ x ⩽ l. Then the solution of difference scheme (4.1)–(4.3) converges to the solution of differential problem (4.1)–(4.3) and the following a priori estimate holds

∥r_{h}^{n+1} − r^{n+1} ∥C̄_{h} + ∥s_{h}^{n+1} − s^{n+1} ∥C̄_{h} ⩽ c_{9}t_{n+1} (h + τ), n ∈ ℕ_{0}.
Proof. Let us define Δr = r_{h} − r, Δs = s_{h} − s, Δa = a_{h} − a. We have

\begin{aligned}Δr_{t,i} + a_{h}Δr̂_{x̄,i} + Δar̂_{x̄,i} & = ψ_{1,i} , i = \overline{1, N}, \\ Δs_{t,i} − a_{h}Δŝ_{x,i} − Δaŝ_{x,i} & = ψ_{2,i} , i = \overline{0, N −1} , \end{aligned}

\begin{aligned}Δr(x, 0) & = 0, Δs(x, 0) = 0, x ∈ ω̄_{h} , \\ Δr(0, t) & = 0, Δs(l, t) = 0, x ∈ ω_{τ}, \end{aligned}

\begin{matrix}ψ_{1} = ( \frac{∂r}{∂t} − r_{t}) + a( \frac{∂r}{∂x} − r̂x̄) , \\ ψ_{2} = ( \frac{∂s}{∂t} − s_{t}) − a( \frac{∂s}{∂x} − ŝ_{x} ) . \end{matrix}

Therefore max_{t∈ω_{τ}}∥ψ_{1}(t)∥C̄_{h} +max_{t∈ω_{τ}}∥ψ_{2}(t)∥C̄_{h} ⩽ c_{10} ·(h+τ). Let us write this difference

scheme in a canonical form:

(6.1)

\begin{matrix}C_{i}^{n}Δr_{h,i}^{n+1} = A_{i}^{n}Δr_{h,i−1}^{n+1} + F_{r,i}^{n} , i = \overline{1, N}, Δr_{h,0}^{n+1} = 0, \\ C_{i}^{n}Δs_{h,i}^{n+1} = A_{i}^{n}Δs_{h,i+1}^{n+1} + F_{s,i}^{n} , i = \overline{0, N −1} , Δs_{h,N}^{n+1} = 0, \end{matrix}

A_{i}^{n} = \frac{a_{h,i}^{n} τ}{h}, C_{i}^{n} = 1 + A_{i}^{n},

F_{r,i}^{n} = Δr_{h,i}^{n} + τΔa_{i}^{n} r̂_{x̄,i} + τψ_{1,i} , F_{s,i}^{n} = Δs_{h,i}^{n} + τΔa_{i}^{n} ŝ_{x,i} + τψ_{2,i} .

For the coefficients of difference scheme (6.1) the conditions of lemmas 4.1 and 4.2 are fulfilled. Therefore

∥Δr^{n+1} ∥C̄_{h} ⩽ ∥Δr^{n}∥C̄_{h} + τc_{5}c_{11}(∥Δr^{n}∥C̄_{h} + ∥Δs^{n}∥C̄_{h} ) + τ∥ψ_{1}^{n}∥C̄_{h} ,

∥Δs^{n+1} ∥C̄_{h} ⩽ ∥Δs^{n}∥C̄_{h} + τc_{5}c_{12}(∥Δr^{n}∥C̄_{h} + ∥Δs^{n}∥C̄_{h} ) + τ∥ψ_{2}^{n}∥C̄_{h} ,

where c_{11} = max_{(x,t)∈\overline{Q}_{T}} | \frac{∂r}{∂x} |, c_{12} = max_{(x,t)∈\overline{Q}_{T}} | \frac{∂s}{∂x} |. Summing these estimates, we obtain:

∥Δr^{n+1} ∥C̄_{h} + ∥Δs^{n+1} ∥C̄_{h}

⩽ (1 + τc_{5}(c_{11} + c_{12})) (∥Δs^{n}∥C̄_{h} + ∥Δr^{n}∥C̄_{h} ) + τ (∥ψ_{1}^{n}∥C̄_{h} + ∥ψ_{2}^{n}∥C̄_{h} )

⩽ e^{c_{5} (c_{11} +c_{12} )τ} (∥Δs^{n}∥C̄_{h} + ∥Δr^{n}∥C̄_{h} ) + τ (∥ψ_{1}^{n}∥C̄_{h} + ∥ψ_{2}^{n}∥C̄_{h} )
⩽ ∑_{k=0}^{n} τe^{t_{n−k} c_{5} (c_{11} +c_{12} )} (∥ψ_{1}^{k}∥C̄_{h} + ∥ψ_{2}^{k}∥C̄_{h} ) ⩽ c_{9}t_{n+1} (h + τ).  □

7. Numerical experiment

Let us consider differential problem (2.1) with initial and boundary conditions

\begin{matrix}v_{0}(x) = ((1 − cos( \frac{2πx}{l} )) c_{1}^{−1} + 1) ^{− \frac{1}{ε}} , u_{0}(x) ≡ 0, 0 ⩽ x ⩽ l, \\ μ_{1}(t) = μ_{2}(t) ≡ 0, t ⩾ 0, \end{matrix} ; (7.1)

\begin{matrix}v_{0}(x) ≡ 1, u_{0}(x) ≡ 0, 0 ⩽ x ⩽ l, \\ μ_{1}(t) = μ_{2}(t) = 1 − cost, t ⩾ 0. \end{matrix} ; (7.2)

In all experiments we use the following values of the parameters: c_{1} = 1, ε = \frac{1}{3}, l = 1, N = 1000, τ = h. The specified value of parameter ε corresponds to the case of monoatomic gas (γ = \frac{5}{3}) [14].

Let us note that all input data satisfy conjugation conditions (2.11)–(2.12). Input data (7.1) satisfy condition (2.4) on the initial data and conditions (2.9) and (2.10) on the boundary data, but do not satisfy the condition (2.5) on derivatives of the initial data, as in this case −u′_{0}(x) − c_{1}εv_{0}(x)^{−(ε+1)} v′_{0}(x) = \frac{2π}{l} sin( \frac{2πx}{l} ) > 0, x ∈ (0, \frac{l}{2}).

In the case of input data (7.2) both conditions (2.4)–(2.5) on the initial data and a condition (2.9) on the boundary data are satisfied, but condition (2.10) on the derivatives of the boundary data is not satisfied, as μ′_{1}(t) = μ′_{2}(t) = sint > 0 at 0 < t < π.

Thus in all the cases shock waves arise which is shown in Fig. 1–2.

(a)  (b)

Figure 1. Profiles of velocity and specific volume at the initial moment t = 0 (a); the solution of difference scheme (4.1)–(4.3) with input data (7.1) at time t = 0.1 (b).

(a)  (b)

Figure 2. Profiles of velocity and specific volume at the initial moment t = 0 (a); the solution of difference scheme (4.1)–(4.3) with input data (7.2) at time t = 1.25 (b).

Thus, necessity of conditions (2.10) and (2.14) for a stability of the difference scheme (4.1)–(4.3) is experimentally confirmed.

8. Conclusions

We obtained necessary conditions for the absence of shock waves for the IBVP for the system of equations for isentropic gas in Riemann invariants in Lagrangian coordinates and proved a priori estimates of stability for monotone difference schemes approximating this problem. In investigating the stability, we used restrictions only for the input data (initial and boundary conditions). On the basis of the investigations performed we can draw the following conclusions.

1. To get proper a priori estimates expressing stability of the difference scheme or its continuous dependence on the input data, first we need to prove the existence of a solution of the difference problem in strong norms. In our case we have proved necessary conditions, formulated the solvability result and have shown that this conditions are sufficient for the global stability of the difference scheme.

2. The uniqueness of the solution of the difference scheme follows from the stability.

3. Once a priori estimates of the stability have been obtained, the investigation of the convergence of the difference solution to the differential one becomes much simpler since the problem for the error of the method can be written in the linear form. To use the Lax theorem (from the approximation and stability convergence follows), it is necessary to prove the stability with respect to the right-hand side. Here we restricted ourselves to considering only homogeneous equations of the gas dynamics.

4. The system of equations for isentropic gas in Lagrangian coordinates is hyperbolic without vacuum only. In our case we proved for both differential and difference problems the absence of vacuum at any moment of time at any point of space.

5. We obtained conditions that guarantee the absence of shock waves. As was shown in the numerical experiment, the emergence of a shock wave is connected both with behavior of the initial and boundary conditions.

6. It is essential that the constructed difference scheme be monotone not only with respect to the approximated solution, but also to its derivatives. In that case, the maximum principle can be used to prove nonlinear stability of the difference scheme.

References

1. P. Bergeret, Classification of smooth solutions to 2 × 2 hyperbolic systems with boundary damping, Math. Methods Appl. Sci. 20 (1997), 1563–1598.

2. E. Godlewski and P.-A. Raviart, Hyperbolic Systems of Conservation Laws, Paris: Ellipses, 1991.

3. S. Godunov, Difference method for numerical calculation of discontinuous solutions of the fluid dynamics equations, Mat. Sb., N. Ser. 47(89):3 (1959), 273–306 (in Russian)

4. S. K. Godunov, Equations of Mathematical Physics, Nauka, Moscow, 1979, (in Russian).

5. P. Lax, Development of singularities of solutions of nonlinear hyperbolic partial differential equations, J. Math. Phys. 5:5 (1964), 611–613.

6. P. Matus, S. Lemeshevsky, Stability and monotonicity of difference schemes for nonlinear scalar conservation laws and multidimensional quasi-linear parabolic equations, Comput. Methods Appl. Math. 9:3 (2009), 253–280.

7. G. Marcinkiewicz, P. Matus, M. Chuiko, Stability of difference schemes in terms of Riemann invariants for a polytropic gas, Comput. Math. Math. Phys. 50:6 (2010), 1024–1037.

8. P. Matus, M. Chuiko, Investigation of the stability and convergence of difference schemes for a polytropic gas with subsonic flows, Differ. Equ. 45:7 (2009), 1074–1085.

9. P. Matus, A. Kolodinska, Nonlinear stability of the difference schemes for equations of isentropic gas dynamics, Comput. Methods Appl. Math. 8:2 (2008), 155–170.

10. A. D. Myshkis, A. M. Filimonov, On the global continuous solvability of the mixed problem for one-dimensional hyperbolic systems of quasilinear equations, Differ. Equ. 44:3 (2008), 413–427.

11. T. Nishida, J. Smoller, Mixed problems for nonlinear conservation laws, J. Differ. Equ. 23:2 (1977), 244–269.

12. V. Ostapenko, On the strong monotonicity of nonlinear difference schemes, Comput. Math. Math. Phys. 38:7 (1998), 1119–1133.

13. V. Ostapenko, Strong monotonicity of finite-difference schemes for systems of conservation laws, Comput. Math. Math. Phys. 39:10 (1999), 1619–1635.

14. B. L. Rozhdestvenskii and N. N. Yanenko, Systems of Quasilinear Equations and Their Applications to Gas Dynamics, Transl. Math. Monogr. 55, Am. Math. Soc., 1983.

15. A. A. Samarskii, The Theory of Difference Schemes, Pure Appl. Math. 240, Marcel Dekker, 2001.

16. A. Samarskii, P. Matus, and P. Vabishchevich, Difference Schemes with Operator Factors, Kluwer, Dordrecht/Boston/London, 2002.

17. Ta-tsien Li, Wen-ci Yu, Boundary Value Problems for Quasilinear Hyperbolic Systems, Duke Univ. Math. Ser. V, 1982.

18. Z. Yanchun, Global smooth solutions of the one dimentional gas dynamics systems, IMA preprint 545, University of Minnesota, 1989, 1–10.

Institute of Mathematics (Received 25 02 2011) NAS of Belarus (Revised 27 04 2011) 11 Surganov Str. 220072 Minsk Belarus

The John Paul II Catholic University of Lublin Al. Raclawickie 14 20-950 Lublin Poland

matus@im.bas-net.by

Institute of Mathematics NAS of Belarus 11 Surganov Str. 220072 Minsk Belarus

mitia87@gmail.com