Компактная схема для решения супердиффузионного уравнения

Полный текст

Введение

Уравнения в частных производных с запаздывающим аргументом широко применяются в математическом моделировании (см., например, [1, 2]), также в последнее время все чаще в математических моделях используются дробные производные. Численные алгоритмы решения уравнений с частными производными с эффектом запаздывания изучались, например, в [3, 4]. Литература по численным методам решения уравнений с дробными производными огромна, отметим лишь работы [5–8], в которых рассматривались численные методы решения уравнения супердиффузионного типа, т. е. с дробной производной по пространственной переменной порядка от 1 до 2. Для супердиффузионных уравнений с функциональным эффектом запаздывания в работе [9] построен и исследован численный метод порядка ${\Delta }^2+h\mathrm{,}$ где $\Delta$ — шаг разбиения по времени, а h — шаг разбиения по пространству, а в работе [10] построен и исследован численный метод порядка ${\Delta }^2+h^2\mathrm{.}$ В работе [11] на основе процедуры экстраполяции Ричардсона был построен метод порядка ${\Delta }^3+h^3\mathrm{,}$ однако при жестком условии пропорциональности шагов $\Delta$ и h.

Целью данной работы было построение для супердиффузионного уравнения с эффектом запаздывания численного метода порядка ${\Delta }^2+h^4\mathrm{,}$ основанного на идеях компактной разностной схемы. Компактная схема изучалась для разных дробных уравнений без запаздывания во многих работах, в частности, для субдиффузионных уравнений в [12, 13], для дробного по времени и пространству уравнения Блоха–Торри в [14], для субдиффузионных уравнений с постоянным запаздыванием в [15].

В данной работе для супердиффузионного уравнения с производной Рисса и эффектом запаздывания строится аналог метода Кранка–Никольсон по временной переменной и аналог компактной схемы для дробной пространственной переменной. Основной результат работы — доказательство сходимости в энергетической норме с порядком ${\Delta }^2+h^4\mathrm{,}$ Техника доказательства потребовала ограничить эффект запаздывания до уравнений с несколькими переменными запаздываниями. В силу переменности запаздываний, существенным элементом алгоритма является интерполяция и экстраполяция с заданными свойствами для функции-предыстории численной модели, в работе применяется кусочно-линейная интерполяция с экстраполяцией продолжением.

1. Постановка задачи и основные предположения

Рассмотрим уравнение супердиффузионного типа с несколькими переменными запаздываниями

$\frac{\partial u(x,t)}{\partial t}=d\frac{{\partial }^{\alpha }u}{\partial |x{|}^{\alpha }}+f(x,t,u(t-{\tau }_1\left(t\right)),u(t-{\tau }_2\left(t\right)),\dots ,u(t-{\tau }_K\left(t\right)\left)\right),$ (1.1)

где $0\le t\le T,0\le x\le X$ — независимые переменные, $u(x,t)$ — искомая функция решения, $d>0$ — коэффициент супердиффузии, величины запаздывания удовлетворяют ограничениям $0\le {\tau }_{k}t\le \tau$ для всех $k=\mathrm{1,}\dots ,K$ и $t\in \left[\mathrm{0,}T\right].$ Дробная производная Рисса по пространству порядка a,  $1<\alpha <\mathrm{2,}$ определяется формулой

$\frac{{\partial }^{\alpha }u(x,t)}{\partial |x{|}^{\alpha }}=-\frac{1}{2\cos (\alpha \pi /2)\Gamma (2-\alpha )}\frac{{\partial }^2}{\partial x^2}\underset{0}{\overset{X}{\int }}\frac{u(\xi ,t)}{|x-\xi {|}^{\alpha -1}}d\xi .$

Заданы граничные условия

$u\left(\mathrm{0,}t\right)=\mathrm{0,}u(X,t)=\mathrm{0,}0\le t\le T,$ (1.2)

и начальные условия

$u(x,s)=\phi (x,s),0\le x\le X,-\tau \le s\le 0.$ (1.3)

Будем предполагать, что решение $u(x,t)$ задачи (1.1)–(1.3) существует и единственно. Кроме того, при доказательстве сходимости численных методов будем предполагать необходимую гладкость решения $u(x,t).$

Дополнительно будем предполагать, что функция $f(x,t,u^1,u^2,\dots ,u^K)$ липшицева по $u^k\mathrm{,}$ т. е. существует постоянная $L_f$ такая, что для всех $x\in \left[\mathrm{0,}X\right],$$t\in \left[\mathrm{0,}T\right],$ $u^1\mathrm{,}{u}^2\mathrm{,}\dots \mathrm{,}{u}^K\mathrm{,}$ $v^1\mathrm{,}{v}^2\mathrm{,}\dots \mathrm{,}{v}^K$ выполняется

$\left|f\right(x,t,u^1,u^2,\dots ,u^K)-f(x,t,v^1,v^2,\dots ,v^K\left)\right|\le L_f\sum _{k=1}^K|u^k-v^k|.$ (1.4)

2. Дискретизация. Метод Кранка–Никольсон с кусочно-линейной интерполяцией и дробным оператором центральных разностей

Введем шаг по времени $\Delta =\frac{\tau }{M_0},$ где $M_0$ — натуральное, и пусть $M=\left[\frac{T}{\Delta }\right].$ Введем точки (узлы по времени) $t_j=j\Delta ,$ $j=-M_0,\dots ,M.$ Полуцелые узлы будем обозначать $t_{j+\frac{1}{2}}=t_j+\frac{\Delta }{2}.$

Разобьем отрезок $\left[\mathrm{0,}X\right]$ на части, введя шаг по пространству $h=\frac{X}{N},$ где N — целое. Введем точки (узлы по пространству) $x_i=ih,$ $i=\mathrm{0,}\dots ,N.$ Аппроксимацию функции $u(x_i,t_j)$ в узлах будем обозначать $U_j^i\mathrm{.}$

При всяком фиксированном $i=\mathrm{0,}\dots ,N$ введем дискретную предысторию к моменту $t_j\mathrm{,}$ $j=\mathrm{0,}\dots ,M:$ ${\left\{U_k^i\right\}}_j=\{U_k^i,j-M_0\le k\le j\}.$

Определение 2.1. Назовем отображение $I:{\left\{U_k^i\right\}}_j\to U_j^i(\cdot )\in C[t_j-\tau ,t_{j+\frac{1}{2}}]$ оператором интерполяции-экстраполяции дискретной предыстории.

Определение 2.2. Будем говорить, что оператор интерполяции-экстраполяции имеет порядок погрешности p на точном решении, если существуют константы $C_1$ и $C_2$ такие, что для всех i, j и $t\in [t_j-\tau ,t_{j+\frac{1}{2}}]$ выполняется неравенство

$\left|U_j^i\right(t)-u(x_i,t\left)\right|\le C_1\underset{j-M_0\le k\le j}{max}|U_k^i-u(x_i,t_k\left)\right|+C_2{\Delta }^p.$

Рассмотрим кусочно-линейную интерполяцию $I\left({\left\{U_k^i\right\}}_j\right)=U_j^i(\cdot ),$ задаваемую соотношениями

$\begin{array}{c}{U}_j^i(t_j+s)=\frac{1}{\Delta }\left(\right(t_k-t_j-s)U_{k-1}^i+(t_j+s-t_{k-1}\left)U_k^i\right),\\ t_{k-1}\le t_j+s\le t_k,j-M_0\le k\le j,\end{array}$ (2.1)

с экстраполяцией продолжением

$U_j^i(t_j+s)=\frac{1}{\Delta }\left(\right(-s)U_{j-1}^i+(\Delta +s\left)U_j^i\right),t_j\le t_j+s\le t_{j+\frac{1}{2}}.$ (2.2)

Кусочно-линейная интерполяция с экстраполяцией продолжением имеет порядок 2, если точное решение $u(x,t)$ дважды непрерывно дифференцируемо по t на промежутке $[-\tau ,T]$ (см. [16, c. 97]). Кусочно-линейная интерполяция обладает следующим свойством. Пусть $u\left(t\right)$ и $v\left(t\right)$ — результаты кусочно-линейной интерполяции по узлам ${\left\{u_k\right\}}_j$ и ${\left\{v_k\right\}}_j$ соответственно и $t\in [t_{k-1},t_k],$ $j-M_0\le k\le j.$ Тогда

$\left|u\right(t)-v(t\left)\right|\le |u_{k-1}-v_{k-1}|+|u_k-v_k|.$ (2.3)

Экстраполяция продолжением обладает следующим свойством. Пусть $u\left(t\right)$ и $v\left(t\right)$ являются результатами экстраполяции продолжением кусочно-линейной интерполяции по узлам ${\left\{u_k\right\}}_j$ и ${\left\{v_k\right\}}_j$ соответственно и $t\in [t_j,t_{j+\frac{1}{2}}].$ Тогда

$\left|u\right(t)-v(t\left)\right|\le \frac{1}{2}|u_{j-1}-v_{j-1}|+\frac{3}{2}|u_j-v_j|.$ (2.4)

Оба эти свойства проверяются по определению кусочно-линейной интерполяции и экстраполяции продолжением. Объединяя (2.3) и (2.4), получаем, что если $u\left(t\right)$ и $v\left(t\right)$ результаты кусочно-линейной интерполяции с экстраполяцией продолжением по узлам ${\left\{u_k\right\}}_j$ и ${\left\{v_k\right\}}_j$ соответственно и $t\in [t_j-\tau ,t_{j+\frac{1}{2}}],$ то

$\left|u\right(t)-v(t\left)\right|\le \frac{3}{2}\left(\right|u_{k-1}-v_{k-1}|+|u_k-v_k\left|\right),$ (2.5)

где $t\in [t_{k-1},t_k]$ в случае интерполяции и $k=j$ в случае экстраполяции. Назовем это свойство квазилипшицевостью с константой 3/2.

Для аппроксимации дробной производной Рисса выберем метод дробных центральных разностей [7]. Пусть

$g_k^{\alpha }=\frac{{(-1)}^k\Gamma (\alpha +1)}{\Gamma (\alpha /2-k+1)\Gamma (\alpha /2+k+1)},k=\mathrm{0,}\mp \mathrm{1,}\mp \mathrm{2,}\dots .$

Тогда, если функция $u(x,t)$ пять раз непрерывно дифференцируема по x на отрезке $\left[\mathrm{0,}X\right],$ то (см. [7])

$\frac{{\partial }^{\alpha }u}{\partial |x{|}^{\alpha }}=-h^{-\alpha }\sum _{k=-\frac{X-x}{h}}^{\frac{x}{h}}{g}_k^{\alpha }u(x-kh,t)+O\left(h^2\right),\left|O\right(h^2\left)\right|\le C_3{h}^2.$

Методом Кранка–Николсон назовем неявную разностную схему

$\begin{array}{c}\frac{U_{j+1}^i-U_j^i}{\Delta }=-\frac{1}{2}d{h}^{-\alpha }\sum _{k=-N+i+1}^{i-1}{g}_k^{\alpha }{U}_{j+1}^{i-k}-\frac{1}{2}d{h}^{-\alpha }\sum _{k=-N+i+1}^{i-1}{g}_k^{\alpha }{U}_j^{i-k}\\ +f(x_i,t_{j+\frac{1}{2}},U_j^i(t_{j+\frac{1}{2}}-{\tau }_1\left(t_{j+\frac{1}{2}}\right)),U_j^i(t_{j+\frac{1}{2}}-{\tau }_2\left(t_{j+\frac{1}{2}}\right)),\dots ,U_j^i(t_{j+\frac{1}{2}}-{\tau }_K\left(t_{j+\frac{1}{2}}\right)),\\ i=\mathrm{1,}\dots ,N-\mathrm{1,}j=\mathrm{0,}\dots ,M-\mathrm{1,}\end{array}$ (2.6)

с начальными условиями

$U_j^i\left(t\right)=\phi (x_i,t),-\tau \le t\le \mathrm{0,}i=\mathrm{0,}\dots ,N,$

и граничными условиями

$U_j^0=\mathrm{0,}{U}_j^N=\mathrm{0,}j=\mathrm{0,}\dots ,M,$

где $U_j^i\left(t\right)$ — результат действия оператора кусочно-линейной интерполяции (2.1) с экстраполяцией (2.2). Обозначим

${\delta }_x^{\alpha }{U}_j^i=h^{-\alpha }\sum _{k=-N+i+1}^{i-1}{g}_k^{\alpha }{U}_j^{i-k},$

${\tilde{f}}_{j+\frac{1}{2}}^i=f(x_i,t_{j+\frac{1}{2}},U_j^i(t_{j+\frac{1}{2}}-{\tau }_1\left(t_{j+\frac{1}{2}}\right)),U_j^i(t_{j+\frac{1}{2}}-{\tau }_2\left(t_{j+\frac{1}{2}}\right)),\dots ,U_j^i(t_{j+\frac{1}{2}}-{\tau }_K\left(t_{j+\frac{1}{2}}\right)\left)\right).$

Тогда метод (2.6) может быть записан в краткой форме

$\frac{U_{j+1}^i-U_j^i}{\Delta }=-\frac{1}{2}d{\delta }_x^{\alpha }{U}_{j+1}^i-\frac{1}{2}d{\delta }_x^{\alpha }{U}_j^i+{\tilde{f}}_{j+\frac{1}{2}}^i.$

Этот метод был изучен в работе [10] в более общем случае функционального запаздывания и двух пространственных переменных. Как следует из результатов работы [10], метод устойчив и сходится с порядком ${\Delta }^2+h^2\mathrm{.}$ Цель этой статьи состоит в построении алгоритма более высокого порядка сходимости на основе метода (2.6).

3. Компактная схема

Введем следующие операторы

${\delta }_x^2{U}_j^i=\frac{U_j^{i-1}-2U_j^i+U_j^{i+1}}{h^2},i=\mathrm{1,}\dots ,N-\mathrm{1,}$

$H{U}_j^i=U_j^i+\frac{\alpha h^2}{24}{\delta }_x^2{U}_j^i,i=\mathrm{1,}\dots ,N-1.$

Сконструируем метод, который по аналогии с подобными методами для других уравнений будем называть компактной схемой:

$\frac{H{U}_{j+1}^i-H{U}_j^i}{\Delta }=-\frac{1}{2}d{\delta }_x^{\alpha }{U}_{j+1}^i-\frac{1}{2}d{\delta }_x^{\alpha }{U}_j^i+H{\tilde{f}}_{j+\frac{1}{2}}^i.$ (3.1)

Перепишем метод в виде

$(H+\frac{\Delta }{2}d{\delta }_x^{\alpha })U_{j+1}^i=(H-\frac{\Delta }{2}d{\delta }_x^{\alpha })U_j^i+\Delta H{\tilde{f}}_{j+\frac{1}{2}}^i.$ (3.2)

Рассмотрим главную матрицу $A=H+\frac{\Delta }{2}d{\delta }_x^{\alpha }$ этой системы. Распишем матрицы

$H=\left(\begin{array}{cccccccc}1-\frac{\alpha }{12}& \frac{\alpha }{24}& \mathrm{...}& 0& 0& 0& \mathrm{...}& 0\\ \frac{\alpha }{24}& 1-\frac{\alpha }{12}& \mathrm{...}& 0& 0& 0& \mathrm{...}& 0\\ .& .& \mathrm{...}& .& .& .& \mathrm{...}& .\\ 0& 0& \mathrm{...}& \frac{\alpha }{24}& 1-\frac{\alpha }{12}& \frac{\alpha }{24}& \mathrm{...}& 0\\ .& .& \mathrm{...}& .& .& .& \mathrm{...}& .\\ 0& 0& \mathrm{...}& 0& 0& 0& \mathrm{...}& 1-\frac{\alpha }{12}\\ & & & & & & & \end{array}\right),$

${\delta }_x^{\alpha }=\left(\begin{array}{cccccccc}{g}_0& g_{-1}& \mathrm{...}& g_{2-i}& g_{1-i}& g_{-i}& \mathrm{...}& g_{2-N}\\ g_1& g_0& \mathrm{...}& g_{3-i}& g_{2-i}& g_{1-i}& \mathrm{...}& g_{3-N}\\ .& .& \mathrm{...}& .& .& .& \mathrm{...}& .\\ g_{i-1}& g_{i-2}& \mathrm{...}& g_{-1}& g_0& g_1& \mathrm{...}& g_{i+1-N}\\ .& .& \mathrm{...}& .& .& .& \mathrm{...}& .\\ g_{N-2}& g_{N-3}& \mathrm{...}& g_{N-i}& g_{N-i-1}& g_{N-i-2}& \mathrm{...}& g_0\\ & & & & & & & \end{array}\right).$

Матрица H симметричная и положительно определенная. В силу свойств коэффициентов $g_0>\mathrm{0,}$ $g_{-k}=g_k<\mathrm{0,}$ $k\ne \mathrm{0,}$ $\sum _{k=-\infty }^{\infty }{g}_k=\mathrm{0,}$ матрица ${\delta }_x^{\alpha }$ также симметричная и положительно определенная. Откуда следует, что матрица A — симметричная и положительно определенная, т. е. система (3.2) однозначно разрешима.

4. Невязка компактной схемы

Невязкой без интерполяции метода (11) назовем

${\psi }_j^i=H\frac{u(x_i,t_{j+1})-u(x_i,t_j)}{\Delta }+\frac{1}{2}d{\delta }_x^{\alpha }u(x_i,t_{j+1})+\frac{1}{2}d{\delta }_x^{\alpha }u(x_i,t_j)-H{f}_{j+\frac{1}{2}}^i,$

$f_{j+\frac{1}{2}}^i=f(x_i,t_{j+\frac{1}{2}},u(x_i,t_{j+\frac{1}{2}}-{\tau }_1\left(t_{j+\frac{1}{2}}\right)),u(x_i,t_{j+\frac{1}{2}}-{\tau }_2\left(t_{j+\frac{1}{2}}\right)),\dots ,u(x_i,t_{j+\frac{1}{2}}-{\tau }_K\left(t_{j+\frac{1}{2}}\right)).$

Лемма 4.1. Если точное решение $u(x,t)$ шесть раз непрерывно дифференцируемо по x и по t, причем смешанные частные производные до шестого порядка непрерывны, функция f дважды непрерывно дифференцируема по x, то невязка без интерполяции метода (3.1) представима в виде

${\psi }_j^i=O({\Delta }^2+h^4),\left|O\right({\Delta }^2+h^4\left)\right|\le C_4({\Delta }^2+h^4).$

Доказательство. Разложим по формуле Тейлора величины, входящие в определение невязки без интерполяции метода (3.1) в окрестности точки с координатами $(x_i,t_{j+\frac{1}{2}}).$ Для краткости будем опускать аргумент $(x_i,t_{j+\frac{1}{2}})$ у функции $u=u(x_i,t_{j+\frac{1}{2}})$ и ее производных. Тогда

${\psi }_j^i=H(\frac{u}{\Delta }+\frac{1}{2}{u}_t^{\text{'}}+\frac{\Delta }{4}{u}_{tt}^{\text{'}\text{'}}+O({\Delta }^2)-\frac{u}{\Delta }+\frac{1}{2}{u}_t^{\text{'}}-\frac{\Delta }{4}{u}_{tt}^{\text{'}\text{'}}+O({\Delta }^2\left)\right)$

$+\frac{d}{2}{\delta }_x^{\alpha }(u+\frac{\Delta }{2}{u}_t^{\text{'}}+O({\Delta }^2)+u-\frac{\Delta }{2}{u}_t^{\text{'}}+O({\Delta }^2\left)\right)-H{f}_{j+\frac{1}{2}}^i$

$=H{u}_t^{\text{'}}+d{\delta }_x^{\alpha }u-H{f}_{j+\frac{1}{2}}^i+O\left({\Delta }^2\right).$ (4.1)

Умножим уравнение (1.1) при $u=u(x_i,t_{j+\frac{1}{2}})$ на H, получим

$H{u}_t^{\text{'}}=dH\frac{{\partial }^{\alpha }u}{\partial |x{|}^{\alpha }}+H{f}_{j+\frac{1}{2}}^i.$ (4.2)

Выполняется (см. [14])

$H\frac{{\partial }^{\alpha }u}{\partial |x{|}^{\alpha }}=-{\delta }_x^{\alpha }u+O\left(h^4\right).$ (4.3)

Тогда из (4.1), (4.2), (4.3) вытекает утверждение леммы.

Определение 4.2. Невязкой метода (3.1) с кусочно-линейной интерполяцией и экстраполяцией продолжением назовем величину

${\widehat{\psi }}_j^i=H\frac{u(x_i,t_{j+1})-u(x_i,t_j)}{\Delta }+\frac{1}{2}d{\delta }_x^{\alpha }u(x_i,t_{j+1})+\frac{1}{2}d{\delta }_x^{\alpha }u(x_i,t_j)-H{\widehat{f}}_{j+\frac{1}{2}}^i,$ (4.4)

${\widehat{f}}_{j+\frac{1}{2}}^i=f(x_i,t_{j+\frac{1}{2}},\widehat{u}(x_i,t_{j+\frac{1}{2}}-{\tau }_1\left(t_{j+\frac{1}{2}}\right)),\widehat{u}(x_i,t_{j+\frac{1}{2}}-{\tau }_2\left(t_{j+\frac{1}{2}}\right)),\dots ,\widehat{u}(x_i,t_{j+\frac{1}{2}}-{\tau }_K\left(t_{j+\frac{1}{2}}\right)),$

где $\widehat{u}(x_i,t_j+s)$ — результат кусочно-линейной интерполяции с экстраполяцией продолжением точного решения:

$\widehat{u}(x_i,t_j+s)=\frac{1}{\Delta }\left(\right(t_k-t_j-s\left)u\right(x_i,t_{k-1})+(t_j+s-t_{k-1}\left)u\right(x_i,t_k\left)\right),$

$t_{k-1}\le t_j+s\le t_k,j-M_0\le k\le j,$

$\widehat{u}(x_i,t_j+s)=\frac{1}{\Delta }\left(\right(-s\left)u\right(x_i,t_{j-1})+(\Delta +s\left)u\right(x_i,t_j\left)\right),t_j\le t_j+s\le t_{j+\frac{1}{2}}.$

Лемма 4.2. При условиях леммы 4.1 невязка метода (3.1) с кусочно-линейной интерполяцией и экстраполяцией продолжением представима в виде

${\widehat{\psi }}_j^i=O({\Delta }^2+h^4),\left|O\right({\Delta }^2+h^4\left)\right|\le C_5({\Delta }^2+h^4).$

Доказательство. Невязка без интерполяции и невязка с кусочно-линейной интерполяцией и экстраполяцией продолжением связаны соотношением

${\widehat{\psi }}_j^i={\psi }_j^i+H{f}_{j+\frac{1}{2}}^i-H{\widehat{f}}_{j+\frac{1}{2}}^i.$ (4.5)

В силу определения оператора H, липшицевости (1.4) функции f и свойств оператора кусочно-линейной интерполяцией с экстраполяцией продолжением получаем

$\left|H\right(f_{j+\frac{1}{2}}^i-{\widehat{f}}_{j+\frac{1}{2}}^i\left)\right|\le \underset{i=\mathrm{1,}\dots ,N-1}{max}|f_{j+\frac{1}{2}}^i-{\widehat{f}}_{j+\frac{1}{2}}^i|$

$\le \text{}\underset{i=\mathrm{1,}\dots ,N-1}{max}{L}_f\sum _{k=1}^K\left|u\right(x_i,t_{j+\frac{1}{2}}\text{}-\text{}{\tau }_k\left(t_{j+\frac{1}{2}}\right))-\widehat{u}(x_i,t_{j+\frac{1}{2}}\text{}-\text{}{\tau }_k\left(t_{j+\frac{1}{2}}\right)\left)\right|\le K{L}_f{C}_2{\Delta }^2.$

Учитывая эту оценку, а также утверждение леммы 4.1, из (4.5) вытекает утверждение леммы 4.2.

5. Сходимость компактной схемы

Рассмотрим векторное пространство размерности $N-\mathrm{1,}$ в котором введем скалярное произведение: если $u=(u^1,u^2,\dots ,u^{N-1})$ и $v=(v^1,v^2,\dots ,v^{N-1}),$ то

$(u,v)=\sum _{i=1}^{N-1}h{u}^i{v}^i.$

Введем энергетическую норму

$\parallel u{\parallel }^2=(u,u),$

компактную норму и норму предыстории погрешности

$\parallel u{\parallel }_H^2=(u,Hu),\parallel {\left\{u_k\right\}}_j{\parallel }_{\tau }=\underset{j-M_0⩽k⩽j}{max}\parallel u_k{\parallel }_H.$

Выполняются оценки (см. [12])

$\frac{2}{3}\parallel u{\parallel }^2\le \parallel u{\parallel }_H^2\le \parallel u{\parallel }^2.$ (5.1)

Справедливы следующие вспомогательные оценки.

Лемма 5.1.

$\parallel Hu\parallel \le \parallel u\parallel .$ (5.2)

Доказательство. Наряду с вектором $u=(u^1,u^2,\dots ,u^{N-1})$ введем векторы смещения $u^{-}=(\mathrm{0,}{u}^1,\dots ,u^{N-2})$ и $u^{+}=(u^2,u^3,\dots \mathrm{,0}).$ По определению энергетической нормы выполняются соотношения $\parallel u^{-}\parallel \le \parallel u\parallel ,$ $\parallel u^{+}\parallel \le \parallel u\parallel .$

Вектор $Hu$ представим в виде

$Hu=u+\frac{\alpha h^2}{24}{\delta }_x^{2}u=u+\frac{\alpha h^2}{24}\frac{u^{-}-2u+u^{+}}{h^2}=\frac{1}{24}(\alpha u^{-}+(24-2\alpha )u+\alpha u^{+}),$

тогда

$\parallel Hu\parallel \le \frac{1}{24}(\alpha \parallel u^{-}\parallel +(24-2\alpha )\parallel u\parallel +\alpha \parallel u^{+}\parallel )\le \frac{1}{24}(\alpha \parallel u\parallel +(24-2\alpha )\parallel u\parallel +\alpha \parallel u\parallel )=\parallel u\parallel .$

Лемма 5.2. Если $a_1\mathrm{,}{a}_2\mathrm{,}\dots \mathrm{,}{a}_K$ — произвольные числа, то

${\left(\sum _{k=1}^K{a}_k\right)}^2\le 2^{K-1}\sum _{k=1}^K{a}_k^2.$

Доказательство. Проводится индукцией по K.

Введем погрешность ${\epsilon }_j^i=u(x_i,t_j)-U_j^i$ и послойную векторную погрешность метода (3.1) как вектор ${\epsilon }_j$ с координатами ${\epsilon }_j^i\mathrm{,}$ $i=\mathrm{1,}\dots ,N-1.$

Теорема 5.1. Пусть выполняются условия леммы 4.1. Тогда

$\parallel {\epsilon }_j\parallel \le C({\Delta }^2+h^4).$

Доказательство. Используя (3.1) и (4.4), выпишем уравнения для координат вектора погрешности

$H\frac{{\epsilon }_{j+1}^i-{\epsilon }_j^i}{\Delta }=-\frac{1}{2}d{\delta }_x^{\alpha }{\epsilon }_{j+1}^i-\frac{1}{2}d{\delta }_x^{\alpha }{\epsilon }_j^i+H({\tilde{f}}_{j+\frac{1}{2}}^i-{\widehat{f}}_{j+\frac{1}{2}}^i)+{\widehat{\psi }}_j^i.$ (5.3)

Введем векторы ${\tilde{f}}_{j+\frac{1}{2}}=({\tilde{f}}_{j+\frac{1}{2}}^1,{\tilde{f}}_{j+\frac{1}{2}}^2,\dots ,{\tilde{f}}_{j+\frac{1}{2}}^{N-1}),$ ${\widehat{f}}_{j+\frac{1}{2}}=({\widehat{f}}_{j+\frac{1}{2}}^1,{\widehat{f}}_{j+\frac{1}{2}}^2,\dots ,{\widehat{f}}_{j+\frac{1}{2}}^{N-1}),$ ${\widehat{\psi }}_{j+\frac{1}{2}}=({\widehat{\psi }}_{j+\frac{1}{2}}^1,{\widehat{\psi }}_{j+\frac{1}{2}}^2,\dots ,{\widehat{\psi }}_{j+\frac{1}{2}}^{N-1})$ и запишем (5.3) в векторном виде

$H\frac{{\epsilon }_{j+1}-{\epsilon }_j}{\Delta }=-\frac{1}{2}d{\delta }_x^{\alpha }{\epsilon }_{j+1}-\frac{1}{2}d{\delta }_x^{\alpha }{\epsilon }_j+H({\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}})+{\widehat{\psi }}_j.$

Умножим это уравнение скалярно на вектор ${\epsilon }_{j+\frac{1}{2}}=\frac{1}{2}({\epsilon }_j+{\epsilon }_{j+1}),$ получим

$\left(H\frac{{\epsilon }_{j+1}-{\epsilon }_j}{\Delta },{\epsilon }_{j+\frac{1}{2}}\right)=-\left(d{\delta }_x^{\alpha }{\epsilon }_{j+\frac{1}{2}},{\epsilon }_{j+\frac{1}{2}}\right)+\left(H({\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}}),{\epsilon }_{j+\frac{1}{2}}\right)+\left({\widehat{\psi }}_j,{\epsilon }_{j+\frac{1}{2}}\right).$ (5.4)

Левая часть этого уравнения преобразуется к виду

$\left(H\frac{{\epsilon }_{j+1}-{\epsilon }_j}{\Delta },{\epsilon }_{j+\frac{1}{2}}\right)=\frac{1}{2\Delta }(\parallel {\epsilon }_{j+1}{\parallel }_H^2-\parallel {\epsilon }_j{\parallel }_H^2).$ (5.5)

Первое слагаемое в правой части (5.4) в силу положительности d и положительной определенности ${\delta }_x^{\alpha }$ оценивается следующим образом

$-\left(d{\delta }_x^{\alpha }{\epsilon }_{j+\frac{1}{2}},{\epsilon }_{j+\frac{1}{2}}\right)\le 0.$ (5.6)

Оценим второе слагаемое в правой части (5.4), учитывая (5.2) и (5.1),

$\text{}\left(H({\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}}),{\epsilon }_{j+\frac{1}{2}}\right)\le \parallel H({\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}})\parallel \parallel {\epsilon }_{j+\frac{1}{2}}\parallel$

$\le \parallel {\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}}\parallel \frac{3}{2}\parallel {\epsilon }_{j+\frac{1}{2}}{\parallel }_H\le \parallel {\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}}\parallel \frac{3}{4}(\parallel {\epsilon }_{j+1}{\parallel }_H+\parallel {\epsilon }_j{\parallel }_H).$ (5.7)

Оценим компоненты вектора ${\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}},$ учитывая липшицевость (1.4) функции f 

$\text{}|{\tilde{f}}_{j+\frac{1}{2}}^i\text{}-{\widehat{f}}_{j+\frac{1}{2}}^i|=\left|f\right(x_i,t_{j+\frac{1}{2}},U_j^i(t_{j+\frac{1}{2}}\text{}-{\tau }_1(t_{j+\frac{1}{2}}\left)\right),U_j^i(t_{j+\frac{1}{2}}\text{}-{\tau }_2(t_{j+\frac{1}{2}}\left)\right),\dots ,U_j^i(t_{j+\frac{1}{2}}\text{}-{\tau }_K(t_{j+\frac{1}{2}}\left)\right))$

$-f(x_i,t_{j+\frac{1}{2}},\widehat{u}(x_i,t_{j+\frac{1}{2}}-{\tau }_1\left(t_{j+\frac{1}{2}}\right)),\widehat{u}(x_i,t_{j+\frac{1}{2}}-{\tau }_2\left(t_{j+\frac{1}{2}}\right)),\dots ,\widehat{u}(x_i,t_{j+\frac{1}{2}}-{\tau }_K\left(t_{j+\frac{1}{2}}\right)\left)\right|$

$\le L_f\sum _{l=1}^K\left|U_j^i\right(t_{j+\frac{1}{2}}-{\tau }_l\left(t_{j+\frac{1}{2}}\right))-\widehat{u}(x_i,t_{j+\frac{1}{2}}-{\tau }_l\left(t_{j+\frac{1}{2}}\right)\left)\right|.$ (5.8)

Учитывая свойство квазилипшицевости (2.5) оператора кусочно-линейной интерполяции с экстраполяцией продолжением, получаем

$\left|U_j^i\right(t_{j+\frac{1}{2}}-{\tau }_l\left(t_{j+\frac{1}{2}}\right))-\widehat{u}(x_i,t_{j+\frac{1}{2}}-{\tau }_l\left(t_{j+\frac{1}{2}}\right)\left)\right|\le \frac{3}{2}\left(\right|{\epsilon }_{k\left(l\right)-1}|+|{\epsilon }_{k\left(l\right)}\left|\right),$ (5.9)

где $t_{j+\frac{1}{2}}-{\tau }_l\left(t_{j+\frac{1}{2}}\right)\in \left[t_{k\left(l\right)-1}{,}_{k\left(l\right)}\right]$ или $k\left(l\right)=j$ в случае $t_{j+\frac{1}{2}}-{\tau }_l\left(t_{j+\frac{1}{2}}\right)>t_j.$ Таким образом, из (5.8) и (5.9) следует

$|{\tilde{f}}_{j+\frac{1}{2}}^i-{\widehat{f}}_{j+\frac{1}{2}}^i|\le L_f\frac{3}{2}\sum _{l=1}^K\left(\right|{\epsilon }_{k\left(l\right)-1}|+|{\epsilon }_{k\left(l\right)}\left|\right).$

По определению энергетической нормы получаем

$\parallel {\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}}{\parallel }^2=\sum _{i=1}^{N-1}h|{\tilde{f}}_{j+\frac{1}{2}}^i-{\widehat{f}}_{j+\frac{1}{2}}^i{|}^2\le {\left(L_f\frac{3}{2}\right)}^2\sum _{i=1}^{N-1}h{\left(\sum _{l=1}^K\right(\left|{\epsilon }_{k\left(l\right)-1}\right|+\left|{\epsilon }_{k\left(l\right)}\right|\left)\right)}^2.$

Из полученного соотношения, используя оценку леммы 5.2, следует

$\parallel {\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}}{\parallel }^2\le {\left(L_f\frac{3}{2}\right)}^2\sum _{i=1}^{N-1}h{2}^{2K-1}\sum _{l=1}^K\left(\right|{\epsilon }_{k\left(l\right)-1}{|}^2+\left|{\epsilon }_{k\left(l\right)}{|}^2\right)$

$={\left(L_f\frac{3}{2}\right)}^2{2}^{2K-1}\sum _{l=1}^K(\parallel {\epsilon }_{k\left(l\right)-1}{\parallel }^2+\parallel {\epsilon }_{k\left(l\right)}{\parallel }^2).$

Извлекая корень, получаем

$\parallel {\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}}\parallel \le L_f\frac{3}{2}{2}^{K-\frac{1}{2}}{\left(\sum _{l=1}^K\right(\parallel {\epsilon }_{k\left(l\right)-1}{\parallel }^2+\parallel {\epsilon }_{k\left(l\right)}{\parallel }^2\left)\right)}^{\frac{1}{2}}\le 3L_f{2}^K\sum _{l=1}^K(\parallel {\epsilon }_{k\left(l\right)-1}\parallel +\parallel {\epsilon }_{k\left(l\right)}\parallel ).$

Оценивая каждое слагаемое суммы в правой части этого соотношения, получаем неравенство

$\parallel {\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}}\parallel \le 3L_f{2}^{K+1}K\underset{j-M_0\le k\le j}{max}\parallel {\epsilon }_k\parallel .$ (5.10)

Из (5.7) и (5.10) вытекает оценка второго слагаемого в правой части (5.4)

$\left(H({\tilde{f}}_{j+\frac{1}{2}}-{\widehat{f}}_{j+\frac{1}{2}}),{\epsilon }_{j+\frac{1}{2}}\right)\le 9L_f{2}^{K-1}K(\parallel {\epsilon }_{j+1}{\parallel }_H+\parallel {\epsilon }_j{\parallel }_H)\underset{j-M_0\le k\le j}{max}\parallel {\epsilon }_k\parallel$

$\le 27L_f{2}^{K-2}K(\parallel {\epsilon }_{j+1}{\parallel }_H+\parallel {\epsilon }_j{\parallel }_H)\underset{j-M_0\le k\le j}{max}\parallel {\epsilon }_k{\parallel }_H.$ (5.11)

Оценим третье слагаемое в правой части (5.4). Имеем

$\left({\widehat{\psi }}_j,{\epsilon }_{j+\frac{1}{2}}\right)\le \parallel {\widehat{\psi }}_j\parallel \parallel {\epsilon }_{j+\frac{1}{2}}\parallel .$ (5.12)

Согласно лемме 4.2 выполняется

$\parallel {\widehat{\psi }}_j{\parallel }^2=\sum _{i=1}^{N-1}h|{\widehat{\psi }}_j{|}^2\le C_5^2{({\Delta }^2+h^4)}^{2}X,$

то есть

$\parallel {\widehat{\psi }}_j\parallel \le C_5\sqrt{X}({\Delta }^2+h^4).$ (5.13)

Учитывая оценку

$\parallel {\epsilon }_{j+\frac{1}{2}}\parallel \le \frac{3}{2}\parallel {\epsilon }_{j+\frac{1}{2}}{\parallel }_H\le \frac{3}{4}(\parallel {\epsilon }_{j+1}{\parallel }_H+\parallel {\epsilon }_j{\parallel }_H),$

из (5.12) и (5.13) получаем

$\left({\widehat{\psi }}_j,{\epsilon }_{j+\frac{1}{2}}\right)\le \frac{3}{4}{C}_5\sqrt{X}({\Delta }^2+h^4)(\parallel {\epsilon }_{j+1}{\parallel }_H+\parallel {\epsilon }_j{\parallel }_H).$ (5.14)

Подставляя (5.5), (5.6), (5.11) и (5.14) в (5.4) и, сокращая на $(\parallel {\epsilon }_{j+1}{\parallel }_H+\parallel {\epsilon }_j{\parallel }_H),$ получаем

$\frac{1}{2\Delta }(\parallel {\epsilon }_{j+1}{\parallel }_H-\parallel {\epsilon }_j{\parallel }_H)\le 27L_f{2}^{K-2}K\underset{j-M_0\le k\le j}{max}\parallel {\epsilon }_k{\parallel }_H+\frac{3}{4}{C}_5\sqrt{X}({\Delta }^2+h^4)$

или

$\parallel {\epsilon }_{j+1}{\parallel }_H\le \parallel {\epsilon }_j{\parallel }_H+\Delta 27L_f{2}^{K-1}K\underset{j-M_0\le k\le j}{max}\parallel {\epsilon }_k{\parallel }_H+\Delta \frac{3}{2}{C}_5\sqrt{X}({\Delta }^2+h^4).$ (5.15)

Обозначим для краткости константы

$A=27L_f{2}^{K-1}K,B=\frac{3}{2}{C}_5\sqrt{X}$

и запишем следующую из (32) оценку для предыстории погрешности

$\parallel {\left\{{\epsilon }_k\right\}}_{j+1}{\parallel }_{\tau }\le (1+A\Delta )\parallel {\left\{{\epsilon }_k\right\}}_j{\parallel }_{\tau }+B\Delta ({\Delta }^2+h^4).$ (5.16)

Из (5.16) выводится оценка

$\parallel {\left\{{\epsilon }_k\right\}}_j{\parallel }_{\tau }\le \frac{B}{A}{e}^{AT}({\Delta }^2+h^4)$

для всех $j=\mathrm{1,}\dots ,M.$ Откуда следует заключение теоремы.

6. Численные эксперименты

Пример 6.1. Рассмотрим тестовое уравнение с постоянным сосредоточенным запаздыванием:

$\frac{\partial u(x,t)}{\partial t}=d\frac{{\partial }^{\alpha }u}{\partial |x{|}^{\alpha }}+3t^2{x}^4{(1-x)}^4+u(x,t-\tau (t\left)\right)-{(t-0.5)}^3{x}^4{(1-x)}^4$

$+\frac{t^{3}d}{2\cos \frac{\pi \alpha }{2}}\left[\frac{\Gamma \left(9\right)}{\Gamma (9-\alpha )}\right(x^{8-\alpha }+{(1-x)}^{8-\alpha })-\frac{4\Gamma \left(8\right)}{\Gamma (8-\alpha )}(x^{7-\alpha }+{(1-x)}^{7-\alpha })$

$+\frac{6\Gamma \left(7\right)}{\Gamma (7-\alpha )}(x^{6-\alpha }+{(1-x)}^{6-\alpha })-\frac{4\Gamma \left(6\right)}{\Gamma (6-\alpha )}(x^{5-\alpha }+{(1-x)}^{5-\alpha })$

$+\frac{\Gamma \left(5\right)}{\Gamma (5-\alpha )}(x^{4-\alpha }+{(1-x)}^{4-\alpha })],$

где $x\in \left[\mathrm{0,1}\right],$ $t\in \left[\mathrm{0,1}\right],$ $\tau \left(t\right)=\mathrm{0.5,}$ $d=100.$

Заданы начальные

$u(x,t)=t^3{x}^4{(1-x)}^4,x\in \left[\mathrm{0,1}\right],t\in [-\mathrm{0.5,0}],$

и граничные условия

$u\left(\mathrm{0,}t\right)=u\left(\mathrm{1,}t\right)=\mathrm{0,}t\in \left[\mathrm{0,1}\right].$

Точным решением является функция

$u(x,t)=t^3{x}^4{(1-x)}^4.$

Введем следующие обозначения:

$E_H(\Delta ,h)=\underset{0\le k\le M}{max}\parallel u_k-U_k{\parallel }_H,$

$u_k=\left(u\right(x_1,t_k),u(x_2,t_k),\dots ,u(x_{N-1},t_k\left)\right),U_k=(U_k^1,U_k^2,\dots ,U_k^{N-1}).$

Метод (3.1) тестировался для различных шагов $\Delta$ и h. Были получены вычислительные порядки сходимости $orde{r}_{\Delta }={\log }_2\left(\frac{E_H(2\Delta ,h)}{E_H(\Delta ,h)}\right),$ $orde{r}_h={\log }_2\left(\frac{E_H\left(\Delta \mathrm{,2}h\right)}{E_H(\Delta ,h)}\right)$ по $\Delta$ и h соответственно.

В таблице 1 показан вычислительный порядок сходимости по h.

Таблица 1

Зависимость погрешностей и вычислительных порядков сходимости от шагов h

Таблица 2 показывает вычислительный порядок сходимости по $\Delta$.

Таблица 2

Зависимость погрешностей и вычислительных порядков сходимости от шагов $\mathit{\Delta }$

Пример 6.2. Рассмотрим тестовое уравнение с переменным запаздыванием:

$\frac{\partial u(x,t)}{\partial t}=d\frac{{\partial }^{\alpha }u}{\partial |x{|}^{\alpha }}+\alpha {(t+1)}^{\alpha -1}{x}^4{(1-x)}^4+u(x,t-\tau (t\left)\right)-{\left(\frac{2(t+1)}{3}\right)}^{\alpha }{x}^4{(1-x)}^4$

$+\frac{{(t+1)}^{\alpha }d}{2\cos \frac{\pi \alpha }{2}}\left[\frac{\Gamma \left(9\right)}{\Gamma (9-\alpha )}\right(x^{8-\alpha }+{(1-x)}^{8-\alpha })-\frac{4\Gamma \left(8\right)}{\Gamma (8-\alpha )}(x^{7-\alpha }+{(1-x)}^{7-\alpha })$

$+\frac{6\Gamma \left(7\right)}{\Gamma (7-\alpha )}(x^{6-\alpha }+{(1-x)}^{6-\alpha })-\frac{4\Gamma \left(6\right)}{\Gamma (6-\alpha )}(x^{5-\alpha }+{(1-x)}^{5-\alpha })$

$+\frac{\Gamma \left(5\right)}{\Gamma (5-\alpha )}(x^{4-\alpha }+{(1-x)}^{4-\alpha })],$