Журнал вычислительной математики и математической физики, 2023, T. 63, № 9, стр. 1438-1445

Численный метод оценки скорости роста ошибок округления в равномерной метрике

М. И. Зуев¹, С. И. Сердюкова^1, *

¹ ЛИТ ОИЯТ
141980 М.о., Дубна, ул. Жолио-Кюри, 6, Россия

Поступила в редакцию 16.02.2023
После доработки 20.03.2023
Принята к публикации 28.04.2023

EDN: WTZUZO
DOI: 10.31857/S0044466923080173

Аннотация

В настоящей работе разработан численно-аналитический алгоритм оценки ошибок округления в равномерной метрике. Установлена их ограниченность на всем интервале вычисления вольт-амперных характеристик длинных джозефсоновских переходов при использовании предлагаемой схемы второго порядка точности. На примере системы двух разностных уравнений показано, как можно исследовать численно скорость роста ошибок округления в равномерной метрике в случае степенной неустойчивости. Кроме того, получены оценки скорости роста ошибок округления в раномерной метрике для схемы Русанова третьего порядка точности. Расчеты проводились на суперкомпьютере “Говорун” с использованием системы REDUCE. Библ. 9. Фиг. 6.

Ключевые слова: конечно-разностные методы, оценка роста ошибок округления в равномерной метрике, численный метод, система REDUCE, суперкомпьютер “Говорун”.

ВВЕДЕНИЕ

Цель этой работы – показать, как можно оценивать скорость роста ошибок округления в равномерной метрике, используя численные методы. Устойчивые в равномерной метрике разностные схемы пригодны для счета разрывных решений. Равномерная метрика является более естественной и наглядной по сравнению с ${{l}_{2}}$ метрикой. Но для получения оценок в равномерной метрике приходится использовать более сложные методы. Так, основные результаты по исследованию устойчивости разностных схем в равномерной метрике были получены с помощью метода перевала. При этом, кроме определения точек перевала и линий наискорейшего спуска, немалого труда стоило доказывать определяющую роль “асимптотических вычетов”. С появлением суперкомпьютеров возникла реальная возможность оценивать скорость роста ошибок округления путем прямого вычисления норм степеней характеристической матрицы $\left\| {{{A}^{n}}} \right\|$. Элементы ${{A}^{n}}$ являются тригонометрическими многочленами $p({{e}^{{i\phi }}}) = \sum\nolimits_{k = - n}^n \,{{a}_{k}}{{e}^{{ik\phi }}}.$ Норма $p({{e}^{{i\phi }}})$ в метрике ${{l}_{1}}$ есть $\sum\nolimits_k \,{\text{|}}{{a}_{k}}{\text{|}}$. Норма $\left\| {{{A}^{n}}} \right\|$ в метрике ${{l}_{1}}$ есть максимум норм ее элементов. В этой работе получены оценки скорости роста ошибок округления в равномерной метрике для трех конкретных примеров.

В разд. 1 получены оценки скорости роста ошибок округления в равномерной метрике для разностной схемы, которая используется при вычислении вольт-амперных характеристик (ВАХ) длинных джозефсоновских переходов. Установлена ограниченность роста ошибок округления в равномерной метрике на всем интервале вычисления ВАХ.

В разд. 2 получены оценки скорости роста ошибок округления в равномерной метрике для системы двух разностных уравнений путем прямого вычисления норм степеней характеристической матрицы в метрике ${{l}_{1}}$. Полученные численно результаты согласуются с теоретическими: в работе [5] была установлена неустойчивость в $C$ рассматриваемой системы порядка $O({{n}^{{1/8}}})$.

В разд. 3 численно получены оценки скорости роста ошибок округления в равномерной метрике для схемы Русанова [6, 7] при различных значениях параметров схемы. Все расчеты были проведены на суперкомпьютере “Говорун” [8] с использованием системы REDUCE [9].

1. ВЫЧИСЛЕНИЕ ВАХ ДЛИННЫХ ДЖОЗЕФСОНОВСКИХ ПЕРЕХОДОВ. ЧИСЛЕННЫЙ МЕТОД ОЦЕНКИ ОШИБОК ОКРУГЛЕНИЯ В РАВНОМЕРНОЙ МЕТРИКЕ

Вычисление вольт-амперных характеристик систем $n$ длинных джозефсоновских переходов связано (см. [1], (7)) с решением системы $n$ сушественно нелинейных дифференциальных уравнений. Был разработан алгоритм, позволяющий свести задачу к решению одного уравнения ${{u}_{{tt}}} = - \beta {{u}_{t}} + {{u}_{{xx}}} - \sin (u) + I.$ Для решения последнего предлагается использовать разностную аппроксимацию второго порядка точности:

$\frac{{u_{\nu }^{{n + 1}} - 2u_{\nu }^{n} + u_{\nu }^{{n - 1}}}}{{{{\tau }^{2}}}} = - \beta \frac{{u_{\nu }^{{n + 1}} - u_{\nu }^{{n - 1}}}}{{2\tau }} + \frac{{u_{{\nu + 1}}^{n} - 2u_{\nu }^{n} + u_{{\nu - 1}}^{n}}}{{{{h}^{2}}}} - \sin (u_{\nu }^{n}) + I(n),$

$\tau $, $h$ – шаги сетки по $t$, $x$ соответственно, ниже используются обозначения $\gamma = \tau {\text{/}}h$, $\delta = \beta \tau {\text{/}}2$. Линейная часть без $I(n)$ рассматриваемого разностного уравнения эквивалентна системе разностных уравнений:

(1)

$u_{\nu }^{{n + 1}} - u_{\nu }^{n} = {v}_{\nu }^{{n + 1}},\quad \frac{{{v}_{\nu }^{{n + 1}} - {v}_{\nu }^{n}}}{{{{\tau }^{2}}}} = - \beta \frac{{{v}_{\nu }^{{n + 1}} + {v}_{\nu }^{n}}}{{2\tau }} + \frac{{u_{{\nu + 1}}^{n} - 2u_{\nu }^{n} + u_{{\nu - 1}}^{n}}}{{{{h}^{2}}}}.$

Решается задача Коши, заданы начальные данные $u_{\nu }^{0}$, $v_{\nu }^{0}$. После преобразования Фурье

${{U}^{n}} = \sum\limits_{\nu = - \infty }^\infty \,{{e}^{{i\nu \phi }}}[u_{\nu }^{n},v_{\nu }^{n}]{\kern 1pt} *$

(знак $*$ означает преобразование вектора-строки в вектор-столбец) и некоторых алгебраических манипуляций система (1) преобразуется в векторную систему:

${{U}^{{n + 1}}} = A({{e}^{{i\phi }}})\,{{U}^{n}},\quad {{U}^{n}} = {{A}^{n}}({{e}^{{i\phi }}})\,{{U}^{0}}.$

Здесь $A({{e}^{{i\phi }}})$ – характеристическая матрица:

$A({{e}^{{i\phi }}}) = \left[ {\begin{array}{*{20}{c}} {(1 - (4{{\gamma }^{2}}{{{\sin }}^{2}}(\phi {\text{/}}2)){\text{/}}(1 + \delta )}&{(1 - \delta ){\text{/}}(1 + \delta )} \\ { - (4{{\gamma }^{2}}{{{\sin }}^{2}}(\phi {\text{/}}2)){\text{/}}(1 + \delta )}&{(1 - \delta ){\text{/}}(1 + \delta )} \end{array}} \right].$

Функция Грина системы (1) есть

$\Gamma _{j}^{n} = \frac{1}{{2\pi }}\int\limits_{ - \pi }^\pi {{{A}^{n}}({{e}^{{i\phi }}})\,{{e}^{{ij\phi }}}d\phi } .$

Порядок роста функции Грина $\Gamma _{j}^{n}$ в метрике ${{l}_{1}}$ определяет скорость роста ошибок округления $\epsilon _{j}^{n}$ в равномерной метрике:

$\mathop {\sup }\limits_{1 \leqslant p,k \leqslant q} \sum\limits_j \left| {\Gamma _{j}^{n}(p,k)} \right| = \left( {\mathop {\sup }\limits_{1 \leqslant j \leqslant q} \left| {\epsilon _{j}^{n}} \right|} \right){\text{/}}\left( {\mathop {\sup }\limits_{1 \leqslant j \leqslant q} \left| {\epsilon _{j}^{0}} \right|} \right).$

В работах [2, 3] были получены оценки, гарантирующие ограниченность скорости роста ошибок округления в ${{l}_{2}}$ метрике для всех $n$.

Собственные значения характеристической матрицы $A$ удовлетворяют уравнению

${{\lambda }^{2}} - 2\frac{{1 - 2{{\gamma }^{2}}{{{\sin }}^{2}}(\phi {\text{/}}2)}}{{1 + \delta }}\lambda + \frac{{1 - \delta }}{{1 + \delta }} = 0.$

Решая это уравнение, находим

(2)

${{\lambda }_{{1,2}}} = \frac{{1 - 2{{\gamma }^{2}}{{{\sin }}^{2}}(\phi {\text{/}}2) \pm \sqrt {{{{(1 - 2{{\gamma }^{2}}{{{\sin }}^{2}}(\phi {\text{/}}2))}}^{2}} - (1 - {{\delta }^{2}})} }}{{1 + \delta }}.$

В работе [2] доказано, что $\left| {{{\lambda }_{1}}} \right| \leqslant 1$ и $\left| {{{\lambda }_{2}}} \right| \leqslant 1$ для всех $\gamma \leqslant 1$ и всех $0 \leqslant \phi \leqslant 2\pi $. Это необходимые условия устойчивости рассматриваемой задачи Коши в метрике ${{L}_{2}}$, гарантирующие отсутствие неустойчивости экспоненциального типа. Но при наличии кратных собственных значений может возникнуть неустойчивость степенного типа [4]: $\left\| {{{G}^{n}}} \right\| \asymp {{n}^{\theta }}$, $\theta > 0$, $G$ – оператор перехода от слоя к слою. При наличии кратных собственных значений нормальная форма характеристической матрицы имеет блочно-треугольную форму [4]. В [2] проведен детальный спектральный анализ характеристической матрицы системы (1). В рассматриваемом случае нормальной формой характеристической матрицы является треугольная матрица второго порядка. Матрица преобразования подобия

$T = \left[ {\begin{array}{*{20}{c}} {{{\lambda }_{1}}}&1 \\ {{{\lambda }_{1}} - 1}&1 \end{array}} \right], \quad \det (T) = 1,$

приводит матрицу $A$ к треугольному виду:

(3)

${{T}^{{ - 1}}}AT = B = \left[ {\begin{array}{*{20}{c}} {{{\lambda }_{1}}}&1 \\ 0&{{{\lambda }_{2}}} \end{array}} \right],\quad {{B}^{n}} = \left[ {\begin{array}{*{20}{c}} {\lambda _{1}^{n}}&{{{D}_{n}}} \\ 0&{\lambda _{2}^{n}} \end{array}} \right],\quad {{D}_{n}} = \frac{{\lambda _{1}^{n} - \lambda _{2}^{n}}}{{{{\lambda }_{1}} - {{\lambda }_{2}}}}.$

Сделав обратное преобразование, получаем

${{A}^{n}} = \left[ {\begin{array}{*{20}{c}} { - {{\lambda }_{1}}{{\lambda }_{2}}{{D}_{n}} - {{D}_{{n + 1}}}}&{{{\lambda }_{1}}{{\lambda }_{2}}{{D}_{n}}} \\ {({{\lambda }_{1}} + {{\lambda }_{2}} - {{\lambda }_{1}}{{\lambda }_{2}} - 1){{D}_{n}}}&{\;\;{{\lambda }_{1}}{{\lambda }_{2}}({{D}_{n}} - {{D}_{{n - 1}}})} \end{array}} \right].$

Так как ${{\lambda }_{1}} + {{\lambda }_{2}} - {{\lambda }_{1}}{{\lambda }_{2}} - 1$ = $2{{\gamma }^{2}}(\cos (\phi ) - 1){\text{/}}(1 + \delta )$, $\gamma \leqslant 1$, оценка скорости роста ошибок округления в равномерной метрике сводится к оценке роста $\left\| {{{D}_{n}}} \right\|$ в метрике ${{l}_{1}}$:

$\frac{{1 - \delta }}{{1 + \delta }}\,\left\| {{{D}_{n}}} \right\| \leqslant \mathop {\sup }\limits_{1 \leqslant p,k \leqslant 2} \sum\limits_j \left| {\Gamma _{j}^{n}(p,k)} \right| \leqslant \frac{4}{{1 + \delta }}\left\| {{{D}_{n}}} \right\|.$

Замечание. Норма $p({{e}^{{i\phi }}}) = \sum\nolimits_{k = - n}^n \,{{a}_{k}}{{e}^{{ik\phi }}}$ в метрике ${{l}_{1}}$ есть $\sum\nolimits_{k = - n}^n \left| {{{a}_{k}}} \right|$. При умножении $p({{e}^{{i\phi }}})$ на $\cos (\phi )$ величина нормы не возрастает: $\left\| {\cos (\phi )p({{e}^{{i\phi }}})} \right\| \leqslant \left\| {p({{e}^{{i\phi }}})} \right\|.$

Алгоритм прямого вычисления норм D_n в метрике ${{l}_{1}}$.

В случае четных $n$, $n = 2nc$, имеем

${{D}_{n}} = \frac{{{{p}_{n}}(z)}}{{{{{(1 + \delta )}}^{{n - 1}}}}},\quad z = 2(1 - {{\gamma }^{2}} + {{\gamma }^{2}}\cos (\phi )),\quad {{p}_{n}}(z) = z\sum\limits_{j = 1}^{nc} \,{{g}_{j}}{{({{\delta }^{2}} - 1)}^{{nc - j}}}{{z}^{{2(j - 1)}}} = \sum\limits_{j = 1}^{nc} \,{{b}_{j}}{{z}^{{2j - 1}}},$

${{g}_{j}} = \sum\limits_{k = 1}^j \,C_{n}^{{2k - 1}}C_{{nc - k}}^{{j - k}},\quad {{z}^{2}} = {{c}_{1}} + {{c}_{2}}\cos (\phi ) + {{c}_{3}}\cos (2\phi ),\quad {{p}_{n}} = \sum\limits_{j = 1}^n \,{{d}_{j}}\cos ((j - 1)\phi ),\quad \left\| {{{p}_{n}}} \right\| = \sum\limits_{j = 1}^n \left| {{{d}_{j}}} \right|.$

Последовательно вычисляем ${{z}^{{2k - 1}}} = {{a}_{1}} + {{a}_{2}}\cos (\phi ) + ... + {{a}_{{2k}}}\cos ((2k - 1)\phi ),$ $2 \leqslant k \leqslant nc$.

Коэффициенты ${{a}_{j}}$ зависят от $k$ и пересчитываются при переходе от ${{z}^{{2k - 1}}}$ к ${{z}^{{2k + 1}}}$:

${{z}^{{2k + 1}}} = {{z}^{{2k - 1}}}({{c}_{1}} + {{c}_{2}}\cos (\phi ) + {{c}_{3}}\cos (2\phi )) = {{a}_{1}}({{c}_{1}} + {{c}_{2}}\cos (\phi ) + {{c}_{3}}\cos (2\phi )) + {{R}_{k}},$

${{R}_{k}} = \sum\limits_{j = 2}^{2k} \,{{a}_{j}}\left[ {{{c}_{1}}\cos ((j - 1)\phi ) + {{c}_{2}}\frac{{\cos ((j - 2)\phi ) + \cos (j\phi )}}{2} + {{c}_{3}}\frac{{\cos ((j - 3)\phi ) + \cos ((j + 1)\phi )}}{2}} \right].$

В случае нечетных $n$, $n = 2nc - 1$, имеем

${{D}_{n}} = \frac{{{{p}_{n}}(z)}}{{{{{(1 + \delta )}}^{{n - 1}}}}},\quad {{p}_{n}}(z) = \sum\limits_{j = 1}^{nc} \,{{g}_{j}}{{({{\delta }^{2}} - 1)}^{{nc - j}}}{{z}^{{2(j - 1)}}},$

${{g}_{j}} = \sum\limits_{k = 1}^j \,C_{n}^{{2(k - 1)}}C_{{nc - k}}^{{j - k}},\quad {{z}^{2}} = {{c}_{1}} + {{c}_{2}}\cos (\phi ) + {{c}_{3}}\cos (2\phi ),\quad {{p}_{n}} = \sum\limits_{j = 1}^n \,{{d}_{j}}\cos ((j - 1)\phi ),\quad \left\| {{{p}_{n}}} \right\| = \sum\limits_{j = 1}^n \left| {{{d}_{j}}} \right|.$

Последовательно вычисляем ${{z}^{{2(k - 1)}}} = {{a}_{1}} + {{a}_{2}}\cos (\phi ) + ...\;{{a}_{{2k - 1}}}\cos (2(k - 1)\phi )$, $2 \leqslant k \leqslant nc$. Коэффициенты ${{a}_{j}}$ зависят от $k$ и пересчитываются при переходе от ${{z}^{{2(k - 1)}}}$ к ${{z}^{{2k}}}$:

${{z}^{{2k}}} = {{z}^{{2(k - 1)}}}({{c}_{1}} + {{c}_{2}}\cos (\phi ) + {{c}_{3}}\cos (2\phi )) = {{a}_{1}}({{c}_{1}} + {{c}_{2}}\cos (\phi ) + {{c}_{3}}\cos (2\phi )) + {{R}_{k}},$

${{R}_{k}} = \sum\limits_{j = 2}^{2k - 1} \,{{a}_{j}}\left[ {{{c}_{1}}\cos ((j - 1)\phi ) + {{c}_{2}}\frac{{\cos ((j - 2)\phi ) + \cos (j\phi )}}{2} + {{c}_{3}}\frac{{\cos ((j - 3)\phi ) + \cos ((j + 1)\phi )}}{2}} \right].$

На фиг. 1 представлены графики $\left\| {{{D}_{n}}} \right\|$ в метрике ${{l}_{1}}$. Нижний график относится к случаю $\tau = h = 0.1$, средний – к случаю $\tau = h = 0.05$, и верхний – для $h = 0.05$, $\tau = h{\text{/}}5$. Аналогично случаю оценки $\left\| {{{D}_{n}}} \right\|$ в метрике ${{l}_{2}}$ (см. фиг. 3, [3]), самый быстрый рост наблюдается при $\gamma \ne 1$. Нижний график выходит на плато $\left\| {{{D}_{{2000}}}} \right\| = 50.5$. Средний график также выходит на плато $\left\| {{{D}_{{2000}}}} \right\| = 100.999...$. А при $\gamma = 1{\text{/}}5$ имеем $\left\| {{{D}_{{2000}}}} \right\| = 491.333034...$, продолжается рост, но видна тенденция выхода на плато:

$\left\| {{{D}_{{1000}}}} \right\| - \left\| {{{D}_{{500}}}} \right\| = 116.386...,\quad \left\| {{{D}_{{1500}}}} \right\| - \left\| {{{D}_{{1000}}}} \right\| = 42.816...,\quad \left\| {{{D}_{{2000}}}} \right\| - \left\| {{{D}_{{1500}}}} \right\| = 15.731...$

Чтобы разобраться, мы продолжили счет до $n = 5500$. В результате выяснилось, что верхний график также выходит на плато, см. фиг. 2, где точками отмечены значения $\left\| {{{D}_{{n + 500}}}} \right\| - \left\| {{{D}_{n}}} \right\|$. Хорошо видно, что с ростом $n$ значения $\left\| {{{D}_{{n + 500}}}} \right\|$ и $\left\| {{{D}_{n}}} \right\|$ сближаются, достигая $\left\| {{{D}_{{5500}}}} \right\| - \left\| {{{D}_{{5000}}}} \right\| = 0.014...$

Фиг. 1.

Графики $\left\| {{{D}_{n}}} \right\|,$ $n \leqslant 2000$ в метрике ${{l}_{1}}$ для уравнения (1). Нижний график относится к случаю $\tau = h = 0.1$. Средний – для $\tau = h = 0.05,$ и верхний график – для $h = 0.05,$ $\tau = h{\text{/}}5$.

Фиг. 2.

Выход на плато графика $\left\| {{{D}_{n}}} \right\|$ в метрике ${{l}_{1}}$ (система (1), $h = 0.05,$ $\tau = h{\text{/}}5$). Точками отмечены значения $\left\| {{{D}_{{n + 500}}}} \right\| - \left\| {{{D}_{n}}} \right\|$, $n = 500$, 1000, 1500, 2000, …, 5000.

Фиг. 3.

Нижний график представляет $\left\| {{{A}^{n}}} \right\|,$ $500 \leqslant n \leqslant 7000$ в метрике ${{l}_{1}}$ для системы (4), $\gamma = 1{\text{/}}2$. Верхний график представляет функцию $\left\| {{{A}^{{500}}}} \right\| * {{n}^{{1/8}}}$.

В рассматриваемом случае оценка скорости роста ошибок округления в равномерной метрике свелась к оценке $\left\| {{{D}_{n}}} \right\|$ в метрике ${{l}_{1}}$. Но это не всегда так. Рассмотрим еще один пример.

2. УСТОЙЧИВАЯ В ${{L}_{2}}$ СИСТЕМА, НЕУСТОЙЧИВАЯ В $C$

Система разностных уравнений

(4)

$\begin{gathered} u_{j}^{{n + 1}} = u_{j}^{n} - \frac{\tau }{{2h}}(v_{{j + 1}}^{n} - v_{{j - 1}}^{n}) - \frac{{{{\tau }^{2}}}}{{{{h}^{2}}}}(u_{{j + 1}}^{n} - 2u_{J}^{n} + u_{{j - 1}}^{n}) = 0, \\ v_{j}^{{n + 1}} = v_{j}^{n} - \frac{\tau }{{2h}}(u_{{j + 1}}^{n} - u_{{j - 1}}^{n}) - \frac{{{{\tau }^{2}}}}{{{{h}^{2}}}}(v_{{j + 1}}^{n} - 2v_{J}^{n} + v_{{j - 1}}^{n}) = 0, \\ \end{gathered} $

устойчива в ${{L}_{2}}$ при $\tau \leqslant h$ (см. [5]).

Характеристическая матрица

$A = \left[ {\begin{array}{*{20}{c}} a&b \\ b&a \end{array}} \right],\quad a = 1 - {{\gamma }^{2}}(1 - \cos (\phi )),\quad b = i\gamma \sin (\phi ),$

имеет два собственных значения ${{\lambda }_{{1,2}}} = a \pm b$. Матрица преобразования подобия

$T = \left[ {\begin{array}{*{20}{c}} {{{\lambda }_{1}}}&1 \\ {{{\lambda }_{1}} - 1}&1 \end{array}} \right]$

приводит характеристическую матрицу к треугольному виду

$TA{{T}^{{ - 1}}} = \left[ {\begin{array}{*{20}{c}} {{{\lambda }_{2}}}&0 \\ b&{{{\lambda }_{1}}} \end{array}} \right] = B,\quad {{B}^{n}} = \left[ {\begin{array}{*{20}{c}} {\lambda _{2}^{n}}&0 \\ {(\lambda _{1}^{n} - \lambda _{2}^{n}){\text{/}}2}&{\;\;\lambda _{1}^{n}} \end{array}} \right].$

Сделав обратное преобразование, получаем

${{A}^{n}} = 0.5\left[ {\begin{array}{*{20}{c}} {\lambda _{1}^{n} + \lambda _{2}^{n}}&{\;\;\lambda _{1}^{n} - \lambda _{2}^{n}} \\ {\lambda _{1}^{n} - \lambda _{2}^{n}}&{\;\;\lambda _{1}^{n} + \lambda _{2}^{n}} \end{array}} \right].$

Скорость роста ошибок округления в равномерной метрике определяется скоростью роста норм элементов матрицы ${{A}^{n}}$ в метрике ${{l}_{1}}$ [5]. В рассматриваемом случае

${{\lambda }_{1}} = 1 - {{\gamma }^{2}} + \frac{{{{\gamma }^{2}} + \gamma }}{2}{{e}^{{i\phi }}} + \frac{{{{\gamma }^{2}} - \gamma }}{2}{{e}^{{ - i\phi }}} = p({{e}^{{i\phi }}}){{e}^{{ - i\phi }}},\quad p(x) = d1 + d2x + d3{{x}^{2}},$

${{\lambda }_{2}} = q({{e}^{{i\phi }}}){{e}^{{ - i\phi }}},$ $q(x) = d3 + d2x + d1{{x}^{2}}$. Последовательно вычисляем ${{p}^{k}}$, ${{q}^{k}}$, $2 \leqslant k \leqslant n$. Сумма модулей коэффициентов многочленов $({{p}^{n}} + {{q}^{n}}){\text{/}}2$, $({{p}^{n}} - {{q}^{n}}){\text{/}}2$ определяет скорость роста ошибок округления в равномерной метрике.

C использованием системы REDUCE удалось вычислить при $\gamma = 1{\text{/}}2$ на суперкомпьютере “Говорун” $r(n) = \left\| {{{A}^{n}}} \right\|$, $1 \leqslant n \leqslant 7000$, в метрике ${{l}_{1}}$, которые медленно растут от $2.0166...$ до $2.645...$, что обеспечивает ограниченность роста ошибок округления в равномерной метрике на интервале $[0,7000\tau ]$. См. нижний график на фиг. 3. Верхний график представляет функцию $r(500){{n}^{{1/8}}}$. В работе [5] аналитически была установлена неустойчивость рассматриваемой системы в $C$ порядка $O({{n}^{{1/8}}})$. На фиг. 4 точками отмечены значения $\Delta (k) = r(14){\text{/}}r(k) - {{(14{\text{/}}k)}^{{1/8}}}$, $k = 1,2,...,13$, где $r(k) = \left\| {{{A}^{{500k}}}} \right\|$. Хорошо видно, что величина $\Delta (k)$ монотонно убывает, достигая в конце значения $0.00159...$, что согласуется с результатом, полученным теоретически.

Фиг. 4.

На графике точками отмечены значения $r(14){\text{/}}r(k) - {{(14{\text{/}}k)}^{{1/8}}}, $ $k = 1,2,...,13$, где $r(k) = \left\| {{{A}^{{500k}}}} \right\|$ в метрике ${{l}_{1}}$. График демонстрирует скорость роста $\left\| {{{A}^{n}}} \right\|$ порядка $O({{n}^{{1/8}}})$ для системы (4) .

3. ОЦЕНКА ОШИБОК ОКРУГЛЕНИЯ В РАВНОМЕРНОЙ МЕТРИКЕ ДЛЯ СХЕМЫ РУСАНОВА

Был разработан и реализован с использованием системы REDUCE на суперкомпьютере “Говорун” численный метод исследования устойчивости в равномерной метрике схемы Русанова третьего порядка точности, которая используется при решении задач газовой динамики. Для модельного уравнения ${{u}_{t}} = {{u}_{x}}$ рассматривается [6] такая разностная схема:

(5)

$\begin{gathered} u_{\nu }^{{n + 1}} = \left\{ {1 + \frac{\gamma }{{12}}( - {{T}^{2}} + 8T - 8{{T}^{{ - 1}}} + {{T}^{{ - 2}}}) + \frac{{{{\gamma }^{2}}}}{8}({{T}^{2}} - 2I + {{T}^{{ - 2}}})} \right. + \\ \left. {\, + \frac{{{{\gamma }^{3}}}}{{12}}({{T}^{2}} - 2T + 2{{T}^{{ - 1}}} - {{T}^{{ - 2}}}) - \frac{\omega }{{24}}({{T}^{2}} - 4T + 6I - 4{{T}^{{ - 1}}} + {{T}^{{ - 2}}})} \right\}u_{\nu }^{n}, \\ n \geqslant 0,\quad \nu = 0,\, \pm {\kern 1pt} 1,\; \pm {\kern 1pt} 2,...;\quad x = \nu h,\quad t = n\tau ,\quad \gamma = \tau {\text{/}}h,\quad Tu_{\nu }^{n} = u_{{\nu + 1}}^{n}. \\ \end{gathered} $

Характеристическая функция рассматриваемой разностной схемы есть

$f({{e}^{{i\phi }}}) = 1 - \frac{{{{\gamma }^{2}}}}{2}{{\sin }^{2}}(\phi ) - \frac{2}{3}\omega {{\sin }^{4}}\left( {\frac{\phi }{2}} \right) + \gamma i\sin (\phi )\left[ {1 + \frac{2}{3}(1 - {{\gamma }^{2}}){{{\sin }}^{2}}\left( {\frac{\phi }{2}} \right)} \right],$

$\begin{gathered} {{\left| {f({{e}^{{i\phi }}})} \right|}^{2}} = 1 + 4P(z){{z}^{2}}{\text{/}}9,\quad P(z) = [4{{\gamma }^{2}}{{(1 - {{\gamma }^{2}})}^{2}} - {{(\omega - 3{{\gamma }^{2}})}^{2}}]{{z}^{2}} - \\ - \;2[2{{\gamma }^{2}}{{(1 - {{\gamma }^{2}})}^{2}} + 3{{\gamma }^{2}}(\omega - {{\gamma }^{2}} - 2)]z + 3(\omega - 4{{\gamma }^{2}} + {{\gamma }^{4}}),\quad z = {{\sin }^{2}}(\phi {\text{/}}2). \\ \end{gathered} $

Последовательно вычисляем суммы модулей коэффициентов $f$, ${{f}^{2}}$, ${{f}^{3}},...$ (cумма модулей коэффициентов ${{f}^{n}}$ есть норма ${{f}^{n}}$ в метрике ${{l}_{1}}$). На фиг. 5 представлены графики $\left\| {{{f}^{n}}} \right\|$ в метрике ${{l}_{1}}$ для схемы Русанова при $\gamma = 1{\text{/}}2$: верхний график относится к $\omega = 3$, средний – к $\omega = 3{\text{/}}2$, два оставшихся графика (практически неразличимые при $n > 500$) относятся к $\omega = 0.5,\;0.25$. Во всех случаях $\left\| {{{f}^{n}}} \right\|$ ограничены при $n \leqslant 3000$, что гарантирует ограниченность ошибок округления в равномерной метрике при $t \leqslant 3000\tau $.

Фиг. 5.

Графики $\left\| {{{A}^{n}}} \right\|$ в метрике ${{l}_{1}}$ – для схемы Русанова при $\gamma = 0.5$: верхний график относится к $\omega = 3$, средний – к $\omega = 3{\text{/}}2$, два оставшихся (практически неразличимых при $n > 500$) относятся к $\omega = 0.5,\;0.25$.

На фиг. 6 представлены графики $\left\| {{{f}^{n}}} \right\|$ в метрике ${{l}_{1}}$ для схемы Русанова при $\omega = 3$: верхний график относится к $\gamma = 0.75$, средний – для $\gamma = 0.5$, нижний относится к $\gamma = 0.25$. Во всех случаях $\left\| {{{f}^{n}}} \right\|$ ограничены при $n \leqslant 3000$, что гарантирует ограниченность ошибок округления в равномерной метрике при $t \leqslant 3000\tau $. Начиная с $n = 1000$, графики практически неразличимы, выходят на одну и ту же асимптоту. Полученные численно результаты согласуются с теоретическими результатами [7].

Фиг. 6.

Графики $\left\| {{{A}^{n}}} \right\|$ в метрике ${{l}_{1}}$ – для схемы Русанова при $\omega = 3$: верхний график относится к $\gamma = 0.75$, средний – для $\gamma = 0.5$, нижний – для $\gamma = 0.25$. При $n > 1000$ графики практически неразличимы, выходят на одну и ту же асимптоту.

Авторы благодарны Г.М. Кобелькову за полезные критические замечания и Е.Е. Тыртышникову за поддержку проводимого исследования.

Список литературы

Башашин М.В., Земляная Е.В., Рахмонов И.Р., Шукринов Ю.М., Атанасова П.Х., Волохова А.В. Вычислительная схема и параллельная реализация для моделирования системы длинных джозефсоновских переходов // Компьютерные исследования и моделирование. 2016. Т. 8. № 4. С. 593–604.
Сердюкова С.И. Моделирование динамических процессов в длинных джозефсоновских переходах. Проблема вычисления вольт-амперных характеристик.Численный метод оценки скорости роста ошибок округления // Ж. вычисл. матем. и матем. физ. 2020. Т. 60. № 1. С. 3–11.
Зуев М.И., Сердюкова С.И. Моделирование динамических процессов в длинных джозефсоновских переходах. Проблема вычисления вольт-амперных характеристик. Оценки скорости роста ошибок округления для разностной схемы второго порядка точности // Ж. вычисл. матем. и матем. физ. 2022. Т. 62. № 1. С. 159–166.
Урм В.Я. О необходимых и достаточных условиях устойчивости систем разностных уравнений // Докл. АН СССР. 1961. Т. 139. № 1. С. 40–43.
Сердюкова С.И. Об устойчивости в равномерной метрике систем разностных уравнений // Ж. вычисл. матем. и матем. физ. 1967. Т. 7. № 3. С. 497–509.
Rusanov V.V. On difference schemes of third order accuracy for nonlinear hyperbolic systems // Jour. of Comp. Phys. 1970. V. 5. № 3. P. 507–516.
Сердюкова С.И. Схема Русанова. Исследование устойчивости в равномерной метрике. Асимптотика в окрестности изолированного разрыва // Сообщение ОИЯИ. 1977. Р5-10708.
Adam Gh., Bashashin M., Belyakov D., Kirakosyan M., Matveev M., Podgainy D., Sapozhnikova T., Streltsova O., Torosyan Sh., Vala M., Valova L., Vorontsov A., Zaikina T., Zemlyanaya E., Zuev M. IT‑ecosystem of the Hybr-iLIT heterogeneous platform for high‑performance computing and training of IT‑specialists. Selected Papers of the 8th International Conference “Distributed Computing and Grid-technologies in Science and Education” (GRID 2018), Dubna, Russia, September 10–14, 2018, CEUR-WS.org/Vol‑2267”.
REDUCE User’s Guide for Unix Systems. Version 3.8 by Winfried Neun ZIB, 2004.

Дополнительные материалы отсутствуют.

Инструменты

следующая статья выпуска предыдущая статья выпуска содержание выпуска

Журнал вычислительной математики и математической физики

Архивы выпусков Информация о журнале Отправить рукопись в журнал