ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ, 2022, том 58, № 7, с.947-961
ЧИСЛЕННЫЕ МЕТОДЫ
УДК 519.63+517.957
КОМПАКТНЫЕ И МОНОТОННЫЕ РАЗНОСТНЫЕ СХЕМЫ
ДЛЯ ОБОБЩЁННОГО УРАВНЕНИЯ ФИШЕРА
© 2022 г. П. П. Матус, Б. Д. Утебаев
Для обобщённого уравнения Фишера с нелинейной конвекцией строятся и изучаются на
стандартных шаблонах монотонные и компактные разностные схемы 4 + 1 и 4 + 2 поряд-
ков аппроксимации при нестандартных соотношениях на временной и пространственный
шаг. Получены априорные оценки разностного решения в нелинейном случае на основе
установленных двусторонних оценок сеточного решения. Эти результаты обобщаются на
двумерное обобщённое уравнение Фишера с нелинейной конвекцией. Проведённый вычис-
лительный эксперимент иллюстрирует эффективность рассматриваемых методов.
DOI: 10.31857/S037406412207007X, EDN: CELASR
Введение. Разработка численных методов исследования математических моделей реакци-
онно-диффузионных систем открывает новые возможности для изучения нелинейных явлений
в сложных физических, химических, биологических и других системах. К числу таких моде-
лей принадлежит популяционная модель Фишера или Колмогорова-Петровского-Пискунова
(КПП) (см. работы [1, 2])
∂u
∂2u
=
+ u(1 - u),
(1)
∂t
∂x2
в которой учитываются основные механизмы, характерные для популяционной динамики.
Уравнение Фишера (или КПП) и его модификации встречаются в различных задачах, на-
пример, в теории горения, в теории фазовых переходов, в физике плазмы и др. [3]. В книге [4,
с. 454] было предложено обобщение уравнения Фишера (обобщённое уравнение Бюргерса-
Фишера с нелинейной конвекцией)
∂u
∂h(u)
∂2u
+
=
+ u(1 - u)
∂t
∂x
∂x2
для учёта различных видов конвективного движения.
В последнее время много усилий было направлено на разработку компактных схем высоко-
го порядка точности, которые используют только узлы сетки, непосредственно примыкающие
к центральному узлу. Под компактными разностными схемами понимаются разностные схе-
мы повышенных порядков аппроксимации и/или точности, записывающиеся на стандартных
для данного уравнения шаблонах [5, с. 11; 6]. В работе [7] построены компактные разностные
схемы для уравнений конвекции-диффузии с дивергентными и недивергентными конвектив-
ными слагаемыми, которые базируются на использовании экспоненциальных схем. Следует
отметить также статью [8], в которой рассматриваются компактные и монотонные разност-
ные схемы (4 + 1), т.е. четвёртого порядка по пространственной переменной и первого по
временной для уравнения Фишера (1).
В работе А.А. Самарского [9] предложены безусловно монотонные схемы (нет соотноше-
ний между сеточными шагами и коэффициентами уравнения), равномерно сходящиеся со ско-
ростью O(h2) для эллиптического уравнения второго порядка типа конвекции-диффузии.
Идея построения монотонных разностных схем высокого порядка аппроксимации на стан-
дартных шаблонах была развита в статьях В.К. Полевикова [10, 11], в которых построены
компактные разностные схемы четвёрт√о порядка т√ности при соотношениях на шаги сет-
ки, удовлетворяющих неравенствам 1/
5 ≤ h1/h2 ≤
5.
Монотонные разностные схемы играют важную роль при математическом моделировании
прикладных задач, так как они позволяют получить численные решения без нефизических
947
6∗
948
МАТУС, УТЕБАЕВ
осцилляций (см. [12-14]). Существует достаточно много определений монотонности [15, с. 169;
16; 17, с. 229]. Наиболее часто используют определение, связанное с выполнением принципа
максимума на сеточном уровне [15, с. 169; 18, с. 27; 19]. В линейном случае оно достаточно
близко к определению монотонности схемы по Фридрихсу [20], т.е. в канонической форме
записи разностного уравнения все коэффициенты неотрицательны, а их сумма равна единице
(требование аппроксимации). В нелинейном случае такие схемы позволяют получить не только
априорные оценки, но и двусторонние оценки разностного решения для начально-краевых
задач для параболических уравнений [21], которые согласованы с аналогичными оценками в
дифференциальном случае, полученными О.А. Ладыженской [22, с. 22].
Настоящая работа посвящена построению и исследованию компактных и монотонных раз-
ностных схем на стандартных шаблонах для обобщённого уравнения Фишера с нелинейной
конвекцией. В п. 1 приводятся вспомогательные результаты, необходимые для дальнейшего ис-
следования разностных схем, в частности, лемма о двусторонних оценок разностного решения
и компактные безусловно монотонные схемы четвёртого порядка точности для одномерного
стационарного уравнения конвекции-диффузии.
В п. 2 построены компактные и монотонные разностные схемы с (4 + 1) порядком ап-
проксимации для обобщённого уравнения Фишера с квазилинейной конвекцией. Доказана мо-
нотонность таких схем при нестандартных соотношениях на временной и пространственный
шаг τ ≥ c1h2. Получены априорные оценки разностного решения в нелинейном случае на
основе установленных двусторонних оценок сеточного решения. Также строятся компактные
схемы (4 + 2). Для реализации нелинейной разностной схемы построен итерационный метод,
сохраняющий второй порядок аппроксимации по временной переменной при специальном вы-
боре начального приближения. Приведены результаты тестовых расчётов, подтверждающие
повышенный порядок точности разностных схем на стандартных шаблонах.
В п. 3 полученные результаты обобщаются на двумерное обобщённое уравнение Фишера с
нелинейной конвекцией.
1. Вспомогательные результаты.
1.1. Двусторонние оценки. Пусть задано начальное число точек - сетка ωh. Окрест-
ностью точки x называется множество M′(x) = M(x)\x, M(x) - шаблон. Пусть заданы
функции A(x), B(x, ξ), F (x), определённые при x ∈ ωh и принимающие вещественные зна-
чения. Далее каждой точке x ∈ ωh сопоставим одно и только одно уравнение вида [15, с. 226]
∑
A(x)y(x) =
B(x, ξ)y(ξ) + F (x), x ∈ ωh,
(2)
ξ∈M′(x)
называемое канонической формой записи разностной схемы. В соответствии с [15, с. 226] точка
x называется граничным узлом сетки, если в ней задано условие Дирихле:
y(x) = μ(x), x ∈ γh,
(3)
где γh - множество граничных узлов. Отметим, что при аппроксимации граничных условий
второго или третьего рода сетка может не содержать граничных узлов, т.е. все точки сетки
будут являться только внутренними узлами. Будем предполагать, что выполняются обычные
условия положительности коэффициентов:
A(x) > 0, B(x, ξ) > 0 для всех ξ ∈ M′(x),
(4)
∑
D(x) = A(x) -
B(x, ξ) > 0 для всех ξ ∈ M′(x).
(5)
ξ∈M′(x)
В соответствии с монографией А.А. Самарского [15, с. 169] разностные схемы (2), (3),
удовлетворяющие условиям (4), (5), будем называть монотонными.
В дальнейшем будем пользоваться следующими сеточными нормами:
∥ · ∥¯C = max| · |,
∥ · ∥C = max| · |.
x∈ωh
x∈ωh
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
КОМПАКТНЫЕ И МОНОТОННЫЕ РАЗНОСТНЫЕ СХЕМЫ
949
Сформулируем основные результаты, позволяющие установить двусторонние оценки се-
точного решения через входные данные задачи при незнакоопределённых входных данных
задачи F (x).
Лемма (см. [16]). Пусть выполнены условия положительности коэффициентов (4), (5).
Тогда максимальное и минимальное значения решения разностной схемы (2) принадлежит
интервалу изменения входных данных
F (x)
F (x)
min
≤ y(x) ≤ max
(6)
x∈ωh D(x)
x∈ωh D(x)
Следствие 1 (см. [15, с. 230]). Пусть выполнены условия леммы 1. Тогда для решения
разностной схемы (2) справедлива оценка в сеточном аналоге нормы C
∥y∥C = max|y(x)| ≤ ∥F/D∥C .
x∈ωh
1.2. Компактные и безусловно монотонные схемы четвёртого порядка аппрок-
симации для одномерного уравнения конвекции-диффузии. Отметим, что рассмат-
риваемые ниже конечно-разностные методы являются следствием разностных схем, постро-
енных В.К. Полевиковым для двумерного стационарного уравнения конвекции-диффузии в
статье [11].
Рассмотрим краевую задачу
(
)
d
du
du
k(x)
- r(x)
- q(x)u = -f(x),
0<x<l,
(7)
dx
dx
dx
u(0) = μ1, u(l) = μ2,
(8)
где
0 < k1 ≤ k(x) ≤ k2,
|r(x)| ≤ k3, q(x) ≥ 0.
Для упрощения дальнейших исследований ниже мы будем использовать следующие обо-
значения:
(
)
du
d
du
L(k,r)u = Lu - r(x)
,
Lu =
k(x)
dx
dx
dx
Здесь и далее предполагаем, что решение рассматриваемых дифференциальных задач суще-
ствует, единственно и обладает всеми непрерывными производными, необходимыми по ходу
изложения.
При построении монотонных схем четвёртого порядка точности на трёхточечном шаблоне
(xi-1, xi, xi+1) с использованием равномерной сетки
ωh = ωh
⋃γh = {xi = ih, i = 0,N, h =
= 1/N}, ωh = {xi = ih, i = 1, N - 1}, γh = {x0 = 0, xN = l} будем ориентироваться на
соотношение [15, с. 170]
L(k,r)u ∼ Λ(a,b)u = κ(aux)x - b+aux - b-a(+1)ux
с коэффициентами
(
)
1
4
1
r±
a=6
+
+
> 0, b± =
+ O(h4),
k(-1)
k(-0.5)
k
k
1
h|r|
r+ = 0.5(r + |r|) ≥ 0, r- = 0.5(r - |r|) ≤ 0, κ =
> 0, R =
> 0.
(9)
1+R+R2 +R3
2k
Для аппроксимации коэффициента k(x) используется шаблонный функционал, предложен-
ный А.А. Самарским в работе [23]. Полученная ниже схема не обобщается на квазилинейные
уравнения параболического типа ввиду неопределённости коэффициента в точке k(0.5).
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
950
МАТУС, УТЕБАЕВ
Здесь используются стандартные обозначения теории разностных схем (см. [15]):
(-1)
u-u
u(+1) - u
ux =
,
ux =
,
u = u(x), u(+1) = u(x + h), u(-1) = u(x - h),
h
h
k = k(x), k(0.5) = k(x + 0.5h), k(-0.5) = k(x - 0.5h).
На равномерной сетке ωh дифференциальную задачу (7), (8) заменим неконсервативной
разностной схемой
ãΛ(a,b)y - qy +f˜ = 0,
(10)
y0 = μ1, yN = μ2,
(11)
где
1
Λ(a,b)y = κ(ayx)x - b+ayx - b-a(+1)yx,
κ=
,
R=h|r|,
1+R+R2 +R3
2k
2
h
(r)2
r dk
d
(r)
ã=1+
p1,
p1 = p1 + σ1,
b=r
+ O(h4), p1 =
+
-2
,
12
k
k
k2 dx
dx k
)
(
)
2
1
(r
h
d
(q)
r
( r)
r=
+
p2
,
p2 = p2 + 2
+σ1
,
p2 = L(k,r)
,
k k
12
dx k
k
k2
[
) ]
[
) ]
2
h
(q)
(q
h2
(f)
(q
q=q+
L(k,r)
+
+σ1
q ,
f =f+
L(k,r)
+
+σ1
f
(12)
12
k
k
12
k
k
В соответствии с работами [23; 24, с. 70] покажем, что разностная схема (10), (11) имеет
порядок аппроксимации O(h4), т.е. для её невязки
ψ = ãΛ(a,b)u - L(k,r)u - qu +f˜+ qu - f
имеет место априорная оценка
∥ψ∥C ≤ Mh4, M = const > 0.
Разностный оператор Λ(a,b) аппроксимирует соответствующий дифференциальный опера-
тор L(k,r)u со вторым порядком. С учётом предположений (9) в отношении коэффициента a
справедливы следующие асимптотические разложения:
√
(
)
(
)
du
h
h2
k d
1
h3
1
aux = k
-
Lu +
√
Lu
-
L
Lu
+ O(h4),
dx
2
6
dx
k
24
k
√
(
)
(
)
du
h
h2
k d
1
h3
1
a(+1)ux = k
+
Lu +
√
Lu
+
L
Lu
+ O(h4).
(13)
dx
2
6
dx
k
24
k
С помощью разложений (13) приходим к соотношению
2
h
Λ(a,b)u = L(k,r)u +
L1u + O(h4),
(14)
12
где
(
(
)
(
))
√
1
r
d
1
L1u = L
Lu
-2
k
√ Lu
k
k
dx
k
Далее, чтобы избавиться от производных высоких порядков, не поддающихся разностной
аппроксимации на минимальном шаблоне, упростим выражение L1u. Так как
r
du
Lu = L(k,r)u +
k
,
k
dx
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
КОМПАКТНЫЕ И МОНОТОННЫЕ РАЗНОСТНЫЕ СХЕМЫ
951
то получим
(
)
)
1
( r
du
L1u = L(k,r)
L(k,r)u
+L(k,r)
k
+
k
k2
dx
)
(
)
r
d
(1
r
du
r√
d
1
+
k
L(k,r)u +
k
-2
k
√ Lu
=
k
dx k
k2
dx
k
dx
k
(
)
1
du
=L(k,r)
L(k,r)u
+p2k
- p1Lu,
(15)
k
dx
где p1, p2 определяются из формул (12).
Подставив в (15) L(k,r)u = qu - f из уравнения (7) и добавив параметр регуляризации
-σ1, равный тождественно нулю на точном решении u(x), приходим к соотношению
(
)
1
du
L1u = L(k,r)
(qu - f)
+p2k
- p1Lu - σ1(L(k,r)u - qu + f) =
k
dx
(
)
)
q
du
(f
=L(k,r)
u
- p1Lu + p2k
- σ1L(k,r)u - (f - qu)σ1 - L(k,r)
(16)
k
dx
k
Так как
(
)
(
)
(
)
)
q
q
r
d
q
(q
q
du d
(q)
L(k,r)
u
=L
u
-
k
u
=L(k,r)
u+
L(k,r)u + 2k
=
k
k
k
dx k
k
k
dx dx k
)
(q
q
du d
(q)
=L(k,r)
u+
(q - f) + 2k
,
k
k
dx dx k
то
(q)
q
du d
(q)
L1u = L(k,r)
u+
(q - f) + 2k
-
k
k
dx dx k
)
du
(f
- p1Lu + p2k
- σ1L(k,r)u - (f - qu)σ1 - L(k,r)
=
dx
k
(
[
) ]
d
(q))
du
(q)
(q
= -p1Lu + p2 + 2
k
-σ1L(k,r)u + L(k,r)
+
+σ1
q u-
dx k
dx
k
k
[
) ]
(f)
(q
- L(k,r)
+
+σ1
f
(17)
k
k
С учётом выражений (9), (12), (14) и (17) получим
2
h
h2
du
ψ = ã(Λ(k,r)u - L(k,r)u) +
p1Lu -
p2k
-(q-q)u+
f - f) + O(h4) =
12
12
dx
2
h
h2
= ã(Λ(k,r)u - L(k,r)u) +
L1u + O(h4) =
p1(Λ(k,r)u - L(k,r)u) + O(h4) = O(h4).
12
12
Параметр регуляризации σ1 = O(1) в выражении (16) будем выбирать из условий ã ≥ 1
и q≥ 0, т.е.
{
1
(q)}
σ1 = σ1(x) = - min p1,
L(k,r)
q
k
Для исследования монотонности разностной схемы (10), (11) запишем её в каноническом
виде:
Ciyi = Aiyi-1 + Biyi+1 + Fi, i = 1,N - 1,
y0 = μ1, yN = μ1,
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
952
МАТУС, УТЕБАЕВ
в котором
(
)
(
)
ã
ã
Ai = ai
(κ + hb+)
> 0, Bi = ai+1
(κ - hb-)
> 0,
h2
h2
ãκ
ã
Ci = (ai+1 + ai)
+
(b+ai -b-ai+1) + q > 0, Fi
f, Di = Ci - Ai - Bi > 0.
(18)
h2
h
Так как выполнены все условия положительности коэффициентов (18), то разностная схема
(10), (11) является безусловно монотонной.
Итак, согласно сформулированной выше лемме получаем двустороннюю оценку решения
разностной схемы (10) при произвольных незнакопостоянных входных данных задачи
{
}
{
}
Fi
Fi
min μ1, μ2, min
≤ yi ≤ max μ1,μ2, max
1≤i≤N-1 Di
1≤i≤N-1 Di
Также из следствия 1 получаем, что разностная схема (10) устойчива относительно правой
части и граничных условий, и для её решения справедлива априорная оценка
{
}
Fi
∥y∥¯C ≤ max
|μ1|, |μ2|,
Di
C
2. Компактные разностные схемы (4+1) для обобщённого уравнения Фишера с
нелинейной конвекцией. Целью данной работы является исследование компактных и моно-
тонных разностных схем (в смысле отсутствия ограничений типа Куранта на временной шаг),
для реализации которых не требуется применение итерационных методов [8]. При построе-
ний таких схем на стандартном трёхточечном шаблоне на верхнем слое будем пользоваться
разностной схемой (10).
В области
QT =Ω×[0,T],
Ω= {0 ≤ x ≤ l},
Ω= Ω⋃Γ, Γ - граница области, рассмотрим
начально краевую задачу для обобщённого уравнения Фишера с нелинейной конвекцией
∂u
∂2u
∂u
=
- h′(u)
+ u(1 - u),
0<x<l,
0<t≤T,
(19)
∂t
∂x2
∂x
с начальным и граничными условиями
u(x, 0) = u0(x), u0(x) ≥ 0, x ∈Ω,
(20)
u(0, t) = μ1(t) ≥ 0, u(l, t) = μ2(t) ≥ 0, t ∈ (0, T ],
(21)
т.е. априори предполагаем неотрицательность входных данных.
Введём стандартные обозначения теории разностных схем [15, с. 260]:
v-v
v = vni = v(xi,tn),
v = vn+1i = v(xi,tn+1), vt =
τ
На равномерной сетке ω = ωh × ωτ = {(xi, tn) ∈QT }, ωh = {xi = ih, i = 0,N, hN = l},
ωτ = {tn = nτ, n = 0, N0, τN0 = T },
ωτ = ωτ
⋃ {tN0 = T } дифференциальную задачу
(19)-(21) заменим разностной схемой
(1 + λ)yt = σ(κŷxx - b+ ŷx - b- ŷx) + (1 - σ)(κyxx - b+yx - b-yx) +
2
h
+
(κ(y - yŷ)¯xx - b+(y - yŷ)¯x - b-(y - yŷ)x) +
12
+ ã(κŷxx -b+ ŷx -b- ŷx) + (y - yŷ)(1 + λ),
(22)
y(x, 0) = u0(x), x ∈ ωh
ŷ0 = μ1,
ŷN = μ2,
(23)
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
КОМПАКТНЫЕ И МОНОТОННЫЕ РАЗНОСТНЫЕ СХЕМЫ
953
где
2
h
h2
λ=
σ1 ≥ 0, σ = 1 -
,
12
12τ
1
h|b|
1
R=h|b|
κ=
,
R=
,
κ=
,
,
1+R+R2 +R3
2
1+R+R2 +R3
2
b = b(y) = h′(yni ) + O(h4), b+ = 0.5(b + |b|) ≥ 0, b- = 0.5(b - |b|) ≤ 0,
b=r2 +σ1b,
b+ = 0.5(b + |b|) ≥ 0,
b- = 0.5(b - |b|) ≤ 0,
r1 + σ1
2
h
r2 = (b(y))¯xx - b(y)(b(y))x,
ã=
(r1 + σ1), r1 = (b(y))2 - 2(b(y))x.
(24)
12
Параметр регуляризации σ1 в выражении (24) будем выбирать из условия ã > 0, т.е.
⎧
⎨- min
{r1(y(xi, tn-1))}, если r1 < 0,
0≤i≤N
σ1 = σ1(yn
i
1≤n≤N0
)=⎩
0,
если r1 > 0.
Покажем, что разностная схема (22), (23) имеет четвёртый порядок аппроксимации по
пространству и первый по времени, т.е. для её невязки
ψ = -(1 + λ)ut + σ(κûxx - b+ûx - b-ûx) + (1 - σ)(κuxx - b+ux - b-ux) +
2
h
+
(κ(u - uû)¯xx - b+(u - uû)¯x - b-(u - uû)x) +
12
+ ã(κûxx -b+ ûx -b- ûx) + (u - uû)(1 + λ)
имеет место оценка
∥ψ∥C ≤ M1(h4 + τ), M1 = const > 0.
(25)
С учётом того, что
2
∂2u
∂u
h
κ(u)uxx - b+(u)ux - b-(u)ux ∼
- h′(u)
+
Eu + O(h4),
∂x2
∂x
12
где
(
)
)
∂3u
∂2u
(∂2u(1 - û)
∂
Eu =
- h′(u)
-
- h′(u)
(u(1 - û))
-
∂x2∂t
∂x∂t
∂x2
∂x
((
(
)∂2û
∂2
∂
)∂û)
- (h′(u))2 - 2∂
h′(u)
-
h′(u) + h′(u)
h′(u)
-
∂x
∂x2
∂x2
∂x
∂x
)
(∂2û
∂û
∂u
-σ1
- h′(u)
+ u(1 - û) -
∂x2
∂x
∂t
и
∂3u
∂2u
- h′(u)
∼ κuxxt - b+uxt - b-uxt,
∂x2∂t
∂x∂t
∂2u(1 - û)
∂
- h′(u)
(u(1 - û)) ∼ κ(u(1 - û))¯xx - b+(u(1 - û))¯x - b-(u(1 - û))x,
∂x2
∂x
(
(
)∂2û
∂2
∂
)∂û
(h′(u))2 - 2∂
h′(u)
-
h′(u) + h′(u)
h′(u)
+
∂x
∂x2
∂x2
∂x
∂x
)
(∂2û
∂û
+σ1
- h′(u)
∼ ã(u)(κ(u)ûxx -b+(u)ûx -b-(u)ûx),
∂x2
∂x
приходим к оценке (25).
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
954
МАТУС, УТЕБАЕВ
2.1. Двусторонняя оценка. Для получения априорных оценок сеточного решения запи-
шем разностную схему (22) в каноническом виде:
Ciyn+1i = Aiyn+1i-1 + Biyn+1i+1 + Fni, i = 1,N - 1,
в котором
[(
)
]
2
τ
h
h2
Ai =
1-
-
yn
i-1
(κ + hb+) + ãi(κ + hb+) ,
h2
12τ
12
[(
)
]
τ
h2
h2
Bi =
1-
-
yn
(κ - hb-) + ãi(κ - hb-) ,
i+1
h2
12τ
12
[(
)
]
τ
h2
h2
Ci =
1-
-
yn
(κ + hb+) + ãi(κ + hb+) +
i
h2
12τ
12
[(
)
]
τ
h2
h2
+
1-
-
yn
(κ - hb-) + ãi(κ - hb-) + (1 + λ)(1 + τyni ),
i
h2
12τ
12
1
1
Fi =
(κ + hb+)(1 + τ)yni-1 +
(κ - hb-)(1 + τ)yni+1 +
12
12
[
]
1
1
+ (1 + τ) (1 + λ) -
(κ + hb+) -
(κ - hb-) yni.
12
12
Пусть
m1 = min
{min{μ1(t), μ2(t)}, u0(x)}, m2 = max
{max{μ1(t), μ2(t)}, u0(x)}.
(x,t)∈QT
(x,t)∈QT
Имеет место следующее утверждение.
Теорема 1. Пусть выполнены следующие условия:
2
12
h
h<h0, h20 =
eT , τ > τ0, τ0 =
(26)
m2
12 - h2eT m2
Тогда разностное решение y(x, t), (x, t) ∈ ω, неотрицательно и для всех i = 0, N , k = 0, N0
имеет место двусторонняя оценка
0≤yki ≤etkm2.
(27)
Доказательство. Очевидно, что 0 ≤ m1 ≤ y0i = u0(xi) ≤ m2. По индукции предпола-
гаем, что оценка (27) имеет место для всех k = 1, n. Докажем теперь справедливость этих
неравенств и для k = n + 1. Действительно, в силу условий теоремы (26) и предположения
индукции
τ
τm2
τ
1
Ani ≥ r -
max
yni ≥ r -
eλtn > 0, r =
-
,
12
1≤i≤N-1
12
h2
12
Bni > 0, Cni > 0, Dni > 0, i = 1,N - 1.
В соответствии с двусторонней оценкой (6) имеет место неравенство
mn+11 ≤ yn+1i ≤ mn+12,
где
Fni
Fni
mn+11 = min
≥ 0, mn+12 = max
≥ 0.
0<i<N Dni
0<i<N Dn
i
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
КОМПАКТНЫЕ И МОНОТОННЫЕ РАЗНОСТНЫЕ СХЕМЫ
955
Заметим, что
Fni
(1 + τ) max yni
max
≤
≤ eτmaxyni.
i Dni
1 + τ minyni
i
Следовательно,
yn+1i ≤ max{mn+12,eτ maxyni} ≤ m2etn+1 .
Теорема доказана.
Следствие 2. Для решения разностной схемы (22), (23) имеет место априорная оценка
∥yn∥C ≤ m2etn , n = 0, N0.
2.2. Вычислительный эксперимент. Чтобы продемонстрировать точность настоящего
метода, рассмотрим уравнение
∂u
∂u
∂2u
+u
=
+ u(1 - u),
(28)
∂t
∂x
∂x2
сочетающее механизмы реакции, конвекции и диффузии, с начальным и граничными усло-
виями
1
1
(x)
u(x, 0) =
-
th
,
2
2
4
)
)
1
1
( 5t
1
1
( 5t
l
u(0, t) =
+
th
,
u(l, t) =
+
th
-
(29)
2
2
8
2
2
8
4
Это уравнение известно как уравнение Бюргерса-Фишера из-за свойств конвективного явле-
ния от уравнения Бюргерса и характеристики диффузионного переноса и реакции типа урав-
нения Фишера. Соответствующее точное решение этого уравнения имеет следующий вид [25]:
)
1
1
( 5t
x
u(x, t) =
+
th
-
2
2
8
4
Для определения порядка скорости сходимости по временной и пространственной перемен-
ным (4+1) в нормах L∞ = C и L2 воспользуемся правилом Рунге
∥z(2h, τ)∥L∞(L2)
∥z(h, 2τ)∥L∞(L2)
ph
= log2
,
pτ
= log2
(30)
L∞(L2)
L∞(L2)
∥z(h, τ)∥L∞(L2)
∥z(h, τ)∥L∞(L2)
Для проверки скорости сходимости вдоль пространственного и временного направлений
выберем шаги h и τ так, чтобы выполнялись условия h4 ≥ τ для определения порядка
скорости сходимости разностного решения по пространственной переменной ph∞ и τ ≥ h4
для определения порядка скорости сходимости сеточного решения по временной переменной
pτ∞ (см. [8]).
В табл. 1, 2 приведены порядки скорости сходимости по пространственному и временному
направлениям в нормах L∞ и L2. Величины, представленные в таблицах, соответствуют
моменту времени T = 1.
Таблица 1. Скорость сходимости по пространственному направлению
в нормах L∞ и L2
h = 0.1
τ
∥z∥L∞
ph∞
∥z∥L
phL2
2
h
τ
0.00202609
-
0.00148301
-
h/2
τ /24
0.00012658
3.99714
9.21E-05
4.00917
h/4
τ /44
7.70E-06
4.03757
5.60E-06
4.03892
h/8
τ /84
2.75E-07
4.80735
1.97E-07
4.82915
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
956
МАТУС, УТЕБАЕВ
Таблица 2. Скорость сходимости по временному направ-
лению в нормах L∞ и L2 (h = 0.001)
τ = 0.1
∥z∥L∞
pτ∞
∥z∥L
pτL2
2
τ
0.00203821
-
0.00148323
-
τ /2
0.00101756
0.99716
0.00074039
1.00233
τ /4
0.00050808
0.99145
0.00036966
1.00214
τ /8
0.00025383
1.00569
0.00018467
1.00125
Из результатов, представленных в табл. 1, 2, видно, что указанная разностная схема имеет
четвёртый порядок аппроксимации по пространству и первый по времени.
2.3. Компактные разностные схемы (4 + 2) для обобщённого уравнения Фише-
ра с нелинейной конвекцией. В данном пункте строятся компактные схемы (4 + 2). Для
реализации нелинейной разностной схемы построен итерационный метод, сохраняющий вто-
рой порядок аппроксимации по временной переменной при специальном выборе начального
приближения.
В области
QT рассмотрим начально краевую задачу для обобщённого уравнения Фишера
с нелинейной конвекцией (19)-(21).
На равномерной сетке ω дифференциальную задачу (19)-(21) заменим разностной схемой
(1 + λ)yt = σ(κŷxx - b+ ŷx - b- ŷx) + (1 - σ)(κyxx - b+yx - b-yx) +
2
h
+
(κ(y(0.5) - yŷ)¯xx - b+(y(0.5) - yŷ)¯x - b-(y(0.5) - yŷ)x) +
12
+ ã(κ(y(0.5))¯xx -b+(y(0.5))¯x - b-(y(0.5))x) + (y(0.5) - yŷ)(1 + λ),
(31)
где
2
h
1
h2
ŷ+y
λ=
σ2 ≥ 0, σ =
-
,
y(0.5) =
,
12
2
12τ
2
1
h|b|
1
κ=
,
R=
,
κ=
,
R=h|b|,
1+R+R2 +R3
2
1+R+R2 +R3
2
h′(ŷ) + h′(y)
b = b(y) =
,
b+ = 0.5(b + |b|) ≥ 0, b- = 0.5(b - |b|) ≤ 0,
2
b=r2 +σ2b,
b+ = 0.5(b + |b|) ≥ 0,
b- = 0.5(b - |b|) ≤ 0,
r1 + σ2
2
h
r2 = (b(y))¯xx - b(y)(b(y))x,
ã=
(r1 + σ2), r1 = b(y)2 - 2(b(y))x.
12
Здесь параметр регуляризации σ2 = σ2(y) определяется из условия ã > 0.
Разностная схема (31) является нелинейной. Для её реализации можно использовать сле-
дующий итерационный метод:
k
k+1
k
k
k+1ŷ - y
k
+
(1 +
λ)
= σ(κ
k+1ŷ¯xx - b
k+1ŷx) + (1 - σ)(κy¯xx - b+y¯x - b-yx) +
ŷ x -b-
τ
(
)
)
)
)
2
k+1
k
k+1
k
k+1
h
k
( kŷ + y
( kŷ + y
( kŷ + y
+
κ
-y
ŷ
-b+
-y
ŷ
-b-
-y
ŷ
+
12
2
2
2
xx
x
x
(
)
)
k
k
k+1
k
k
(k+1ŷ + y)
(k+1ŷ + y)
(k+1ŷ + y)
( kŷ + y
k
b+
b-
+
ã
κ
-
-
+
-y
ŷ
(1 +
λ),
2
2
2
2
xx
x
x
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
КОМПАКТНЫЕ И МОНОТОННЫЕ РАЗНОСТНЫЕ СХЕМЫ
957
с выбором начального приближения
⎧
⎨
μ1(tn+1),
i = 0,
n+1
0
yi
2yni - yn-1i, i = 1, N - 1,
=⎪⎩
μ2(tn+1),
i=N.
Так как начальное приближение использует два слоя, то нам нужно предварительно найти
решение также и на первом слое. Для этого используем дифференциальное уравнение (19)
при t = 0, т.е.
∂u
(x, 0) = u1(x), u1(x) = u′′0 - h′(u0)u′0 + u0(1 - u0).
∂t
Аппроксимировав его с порядком O(τ2) относительно точки n + 1/2, получим
y1 = y0 + τu1(x).
Аналогично рассмотренному выше случаю легко показать, что для погрешности аппрок-
симации разностной схемы (31) имеет место оценка
∥ψ∥C ≤ M2(h4 + τ2), M2 = const > 0.
2.4. Вычислительный эксперимент. Рассмотрим уравнение Бюргерса-Фишера (28) с
начальным и граничными условиями (29).
При расчёте по итерационной схеме потребуем выполнения условия
(k+1)
(k)
max
|yn+1i - yn+1i| ≤ ε, ε = 10-9,
0≤i≤N
где k - номер итерации, k = 0, 1, 2, . . .
Для определения порядка скорости сходимости по временной и пространственной перемен-
ным (4 + 2) в нормах L∞ = C и L2 воспользуемся правилом Рунге (30).
В табл. 3, 4 приведены порядки скорости сходимости по пространственному и временному
направлениям в нормах L∞ и L2.
Таблица 3. Скорость сходимости по пространственному направлению
в нормах L∞ и L2
h = 0.1
τ = 0.01
∥z∥L∞
ph∞
∥z∥L2
phL2
k
h
τ
6.96E-05
-
7.38E-05
-
6
h/2
τ /22
4.98E-06
3.80658
5.23E-06
3.81834
4
h/4
τ /42
3.20E-07
3.95658
3.36E-07
3.95920
3
h/8
τ /82
2.02E-08
3.98687
2.12E-08
3.98745
2
h/16
τ /162
1.26E-09
3.99809
1.32E-09
3.99840
2
Таблица 4. Скорость сходимости по временному направ-
лению в нормах L∞ и L2 (h = 0.001)
τ = 0.1
∥z∥L∞
pτ∞
∥z∥L
pτL2
k
2
τ
6.96E-05
-
7.38E-05
-
6
τ /2
1.90E-05
1.86749
2.01E-05
1.87559
5
τ /4
4.98E-06
1.93923
5.23E-06
1.94311
4
τ /8
1.27E-06
1.97088
1.33E-06
1.97267
3
τ /16
3.20E-07
1.98591
3.36E-07
1.98680
3
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
958
МАТУС, УТЕБАЕВ
3. Компактные разностные схемы (4 + 1) для двумерного обобщённого урав-
нения Фишера с нелинейной конвекцией. Далее построим монотонные и компактные
разностные схемы с 4 + 1 порядком аппроксимации на 9-точечном шаблоне для двумерного
обобщённого уравнения Фишера с нелинейной конвекцией.
QT
Пусть в области
= {(x, t) : x ∈Ω, 0 ≤ t ≤ T}, x = (x1,x2),
Ω = {0 ≤ xα ≤ lα,
Ω=Ω⋃∂Ω,где∂Ω-граница,требуетсянайтифункциюu(x,t),удовлетворяю-
α = 1,2},
щую начально-краевой задаче для двумерного обобщённого уравнения Фишера с нелинейной
конвекцией
)
∑
∂u
(∂2u
∂u
=
- h′α(u)
+ u(1 - u), x ∈ Ω, t ∈ (0, T ],
(32)
∂t
∂x2α
∂xα
α=1
с граничными условиями Дирихле
u(x, t) = μ(x, t), (x, t) ∈ ∂Ω × (0, T ]
(33)
и начальным условием
u(x, 0) = u0(x), x ∈Ω.
(34)
В области
QT рассмотрим равномерную сетку ω = ωh × ωτ , где
ωτ = {tn = nτ, n = 0,N0, τ = T/N0,
ωh = ωh
⋃γh},
множество внутренних узлов пространственной сетки определяется соотношением
ωh = {x = (xi11,x22 ), xαα = iαhα, iα = 1,Nα - 1, hαNα = lα, α = 1,2},
через γh обозначено множество её граничных узлов. На сетке ω дифференциальную задачу
(32)-(34) аппроксимируем разностной схемой
∑
∑
h2α
(1 + λ)yt =
(Λ(b)ασα ŷ + Λ(b)α(1 - σα)y) +
Λ(b)αΛ(b)3-αŷ +
12
α=1
α=1
)
∑
2
(σ∗0h2
α
+
Λ(b)α ŷ + ãαΛ(b)αŷ +h
Λ(b)αy(1 - ŷ)
+ (1 + λ)y(1 - ŷ),
(35)
12
12
α=1
y(x, t) = μ(x, t), x ∈ γh, t ∈ ωτ , (x, t) ∈ ωh × ωτ ,
y(x, 0) = u0(x), x ∈ ωh,
(36)
где
2
h
λ=
(σ∗0 + σ∗1),
|h|2 = h21 + h22, σα = 1 -hα
,
12
12τ
,
Λ(b)αy = καyxαxα - bα(y)yxα - bα(y)yxα , Λαb)yŷ = κα(yŷ)xαxα-bα(y)(yŷ)xα-bα(y)(yŷ)
xα
1
hα|bα|
Λ(b)αŷ = καŷxαxα - bα(y)ŷxα - bα(y)ŷxα, κα =
,
Rα =
,
1+Rα +R2α +R3α
2
1
κα =
,
Rα =hα|bα|,
Λ(b)αŷ = καŷxαxα -bα(y)ŷxα -bα(y)ŷxα,
1+Rα +R2α +R3α
2
Λ(b)αΛ(b)3-α = κα(κ3-α ŷx3-αx3-α - b3-α(y)ŷx3-α - b3-α(y)ŷx3-α ) -
- b+α(y)(κ3-αŷx3-αx3-α - b3-α(y)ŷx3-α - b3-α(y)ŷx3-α) -
- b-α(y)(κ3-αŷx3-αx3-α - b3-α(y)ŷx3-α - b3-α(y)ŷx3-α),
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
КОМПАКТНЫЕ И МОНОТОННЫЕ РАЗНОСТНЫЕ СХЕМЫ
959
)
2
h
(h2α
ãα =
r1α + σ∗
1
,
r1α = b2α - 2(bα)
,
xα
12
h2
2
bα =h
(r2α + bσ∗1), r2α = Λ(b)α(bα), b = h′α(yni
).
1i2
12ã
Параметр регуляризации σ∗1 будем выбирать из условия ãα > 0, т.е.
{
-min{r1α}, если r1α < 0,
α
σ∗1 = σ∗1(yni
)=
1i2
0,
если r1α > 0.
В соответствии с работами [23; 24, с. 70] нетрудно показать, что разностная схема (35), (36)
аппроксимирует исходную задачу (32)-(34) четвёртым порядком по пространству и первым по
времени, т.е. для её невязки
∑
∑
h2α
ψ = -(1 + λ)ut + (Λ(b)ασαû + Λ(b)α(1 - σα)u) +
Λ(b)αΛ(b)3-αû +
12
α=1
α=1
)
∑
2
(σ∗0h2
α
+
Λ(b)αû + ãαΛ(b)αû +h
Λ(b)αu(1 - û)
+ (1 + λ)u(1 - û)
12
12
α=1
имеет место априорная оценка
∥ψ∥C ≤ M3(|h|4 + τ), M3 = const > 0.
Для применения принципа максимума разностную схему (36) приведём к каноническому
виду (2) и проверим достаточные условия на коэффициенты (4), (5). После элементарных
преобразований находим
[(
)
]
∑
1
h2α
h2
A = (1 + λ)(1 + τynij) + τ
σα -
ynij +
σ∗
θα + ãα θα -θαθ3-αhα
,
0
h2α
12
12
12h2
α=1
3-α
((
)
)
(
)
τ
h21
h2
1
1
B1 =
σ1 -
yn
i-1j
ξ1 + ã1
ξ1
+τξ1
σ∗0 -
θ(-11)2 -
θ2
,
h21
12
12h21
12h22
12h2
1
((
)
)
(
)
τ
h21
h2
1
1
B2 =
σ1 -
yn
η1 + ã1η1
+τη1
σ∗0 -
θ(+11)2 -
θ2
,
i+1j
h21
12
12h21
12h22
12h2
1
((
)
)
(
)
τ
h22
h2
1
1
B3 =
σ2 -
yn
ξ2 + ã2
ξ2
+τξ2
σ∗0 -
θ1 -
θ(-12)
,
ij-1
1
h22
12
12h22
12h22
12h2
1
((
)
)
(
)
τ
h22
h2
1
1
B4 =
σ2 -
yn
η2 + ã2η2
+τη2
σ∗0 -
θ1 -
θ(+12)
,
ij+1
1
h22
12
12h22
12h22
12h2
1
τ
τ
τ
τ
B5 =
ξ1ξ-112 +
ξ-121ξ2 > 0, B6 =
η1ξ+112 +
η-121ξ2 > 0,
12h22
12h21
12h22
12h2
1
τ
τ
τ
τ
B7 =
ξ1η-112 +
ξ+121η2 > 0, B8 =
η1η+112 +
η+121η2 > 0,
12h22
12h21
12h22
12h2
1
[(
)
]
1
F = (1 + τ)
1+λ-
(θ1 + θ2) ynij + ξ1yni-1j + η1yni+1j + ξ2ynij-1 + η2yn
ij+1
,
12
τ
D = (1 + λ)(1 + τynij) +
(ξ1yni-1j + η1yni+1j + ξ2ynij-1 + η2ynij+1 - (θ1 + θ2)ynij ),
12
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
960
МАТУС, УТЕБАЕВ
где
θα = 2κα + hαb+α - hαb-α,
θα = 2κα + hαbα - hαbα,
ξα = κα + hαb+α,
ξα = κα + hαb+α, ηα = κα - hαb-α,
ηα = κα - hαbα.
Параметр регуляризации σ∗0 выбираем из условия неотрицательности коэффициентов Bα,
α = 1,4 (см. работы [10, 11]), т.е.
1
h1
√
σ∗0 ≥ max{2/h2α},
√
≤
≤
5.
(37)
α
5
h2
Итак, имеет место следующее утверждение.
Теорема 2. Пусть выполнены условия (37) и
{
}
12
h2α
max{hα} < h0, h20 =
eT , τ > τ0, τ0 = max
α=1,2
m3
α=1,2
12 - h2αeT m3
Тогда разностное решение y(x, t), (x, t) ∈ ω, неотрицательно и для всех iα = 0, Nα, α = 1, 2,
k = 0,N0 имеет место двусторонняя оценка
0≤yki
≤etkm3,
1i2
где m3 = max {max{μ1(t), μ2(t), μ3(t), μ4(t)}, u0(x)}.
(x,t)∈QT
Замечание. При выполнении условий теоремы 2 разностная схема (35) является монотон-
ной.
Следствие 3. Для решения разностной схемы (35), (36) имеет место априорная оценка
∥yn∥C ≤ m3etn , n = 0, N0.
СПИСОК ЛИТЕРАТУРЫ
1. Fisher R.A. The wave of advance of advantageous genes // Ann. Hum. Genetic. 1937. V. 7. № 4. P. 353-
369.
2. Kolmogorov A.N., Petrovskii I.G., Piskunov N.S. A study of the equation of diffusion with increase in
the quantity of matter, and its application to a biological problem // Moscow University Mathematics
Bull. 1937. V. 1. P. 1-26.
3. Шаповалов А.В., Трифонов А.Ю. Метод разложения Адомиана для одномерного нелокального
уравнения Фишера-Колмогорова-Петровского-Пискунова // Изв. вузов. Физика. 2019. Т. 62. № 4.
C. 135-143.
4. Murray J.D. Mathematical Biology: I. An Introduction. New York, 2002.
5. Толстых А.И. Компактные разностные схемы и их применение в задачах аэрогидродинамики. М.,
1990.
6. Матус П.П., Хоанг Тхи Киеу Ань. Компактные разностные схемы на трёхточечном шаблоне для
гиперболических уравнений второго порядка // Дифференц. уравнения. 2021. Т. 57. № 7. C. 963-975.
7. Утебаев Б.Д. Компактные разностные схемы для уравнений конвекции-диффузии // Весцi НAH
Беларусi. Сер. фiз.-мат. навук. 2021. Т. 57. № 3. С. 311-318.
8. Матус П.П., Утебаев Б.Д. Компактные и монотонные разностные схемы для параболических урав-
нений // Мат. моделирование. 2021. Т. 33. № 4. С. 60-78.
9. Самарский А.А. О монотонных разностных схемах для эллиптических и параболических уравнений
в случае несамосопряжённого эллиптического оператора // Журн. вычислит. математики и мат.
физики. 1965. Т. 5. № 3. C. 548-551.
10. Полевиков В.К. Схема повышенного порядка точности для задач высокоинтенсивного тепломассо-
обмена // Современные проблемы тепловой гравитационной конвекции. Минск, 1974. C. 84-88.
11. Полевиков В.К. Монотонная разностная схема повышенного порядка точности для двумерных урав-
нений конвекции-диффузии // Журн. Белорус. гос. ун-та. Математика. Информатика. 2019. № 3.
C. 71-83.
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
КОМПАКТНЫЕ И МОНОТОННЫЕ РАЗНОСТНЫЕ СХЕМЫ
961
12. Gaspar F.G., Lisbona F.J., Matus P., Tuyen V.T.K. Monotone finite difference schemes for quasilinear
parabolic problems with mixed boundary conditions // Comput. Methods in Appl. Math. 2016. V. 16.
№ 2. P. 231-243.
13. Matus P., Hieu L.M., Vulkov L.G. Analysis of second order difference schemes on nonuniform grids for
quasilinear parabolic equations // J. of Comput. and Appl. Math. 2017. V. 310. P. 186-199.
14. Matus P., Lemeshevsky S. Stability and monotonicity of difference schemes for nonlinear scalar
conservation laws and multidimensional quasi-linear parabolic equations // Comput. Methods in Appl.
Math. 2009. V. 9. № 3. P. 253-280.
15. Самарский А.А. Теория разностных схем. М., 1983.
16. Матус П.П., Хиеу Л.М., Волков Л.Г. Принцип максимума для разностных схем с незнакопостоян-
ными входными данными // Докл. НАН Беларуси. 2015. Т. 59. № 5. С. 13-17.
17. Галанин М.П., Савенков Е.Б. Методы численного анализа математических моделей. М., 2010.
18. Самарский А.А., Вабищевич П.Н. Численные методы решения задач конвекции-диффузии. М.,
1999.
19. Вабищевич П.Н., Самарский А.А. Монотонные разностные схемы для задач конвекции-диффузии
на треугольных сетках // Журн. вычислит. математики и мат. физики. 2002. Т. 42. № 9. C. 1368-
1382.
20. Friedrichs K.O., Hyers D.H. Symmetric hyperbolic linear differential equations // Comm. on Pure and
Appl. Math. 1954. V. 7. № 2. P. 345-392.
21. Матус П.П., Поляков Д.Б. О согласованных двусторонних оценках решений квазилинейных пара-
болических уравнений и их аппроксимаций // Дифференц. уравнения. 2017. Т. 53. № 7. С. 991-1000.
22. Ладыженская О.А., Солонников В.А., Уральцева Н.Н. Линейные и квазилинейные уравнения па-
раболического типа. М., 1967.
23. Самарский А.А. Схемы повышенного порядка точности для многомерного уравнения теплопровод-
ности // Журн. вычислит. математики и мат. физики. 1963. Т. 3. № 5. C. 812-840.
24. Берковский Б.М., Полевиков В.К. Вычислительный эксперимент в конвекции. Минск, 1988.
25. Wang Xinyi, Lu Yuekai. Exact solutions of the extended Burgers-Fisher equation // Chinese Phys. Lett.
1990. V. 7. № 4. P. 145-147.
Институт математики НАН Беларуси,
Поступила в редакцию 09.03.2022 г.
г. Минск,
После доработки 09.03.2022 г.
Католический университет имени Иоанна-Павла II,
Принята к публикации 25.05.2022 г.
г. Люблин, Польша
7
ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ том 58
№7
2022
