Журнал вычислительной математики и математической физики, 2020, T. 60, № 9, стр. 1503-1512

Аналитико-численное исследование вопроса о разрушении за конечное время решения начально-краевой задачи для нелинейного уравнения Клейна–Гордона

М. О. Корпусов^1, *, А. Н. Левашов^1, **, Д. В. Лукьяненко^1, ***

¹ МГУ им. М.В. Ломоносова, физический факультет, кафедра математики
119992 Москва, Ленинские горы, Россия

^* E-mail: korpusov@gmail.com
^** E-mail: levashov.an16@physics.msu.ru
^*** E-mail: lukyanenko@physics.msu.ru

Поступила в редакцию 11.06.2018
После доработки 10.12.2019
Принята к публикации 09.04.2020

DOI: 10.31857/S0044466920090100

Аннотация

В данной работе рассматривается задача об аналитико-численном исследовании вопроса о разрушении за конечное время решения начально-краевой задачи для нелинейного уравнения Клейна–Гордона. Из аналитических результатов получена оценка сверху на время разрушения решения с произвольной положительной начальной энергией. Данная априорная информация о времени разрушения решения используется при численном уточнении процесса разрушения. Показано, что численная диагностика разрушения решения позволяет как уточнить аналитическую оценку времени разрушения решения, так и выявить факт локального разрушения решения по пространственной переменной. Библ. 17. Фиг. 5.

Ключевые слова: дифференциальные уравнения в частных производных, численная диагностика разрушения решения.

ВВЕДЕНИЕ

В работе рассматривается начально-краевая задача для модельного нелинейного уравнения Клейна–Гордона

(1)

$\begin{gathered} {{u}_{{tt}}} - {{u}_{{xx}}} + u = {{u}^{2}} + {{u}^{3}},\quad x \in (0,l),\quad t \in (0,T], \hfill \\ u(0,t) = u(l,t) = 0,\quad t \in (0,T], \hfill \\ u(x,0) = {{u}_{0}}(x),\quad {{u}_{t}}(x,0) = {{u}_{1}}(x),\quad x \in [0,l], \hfill \\ \end{gathered} $

которое описывает механику кристаллов [1].

Структура данной работы следующая. В разд. 1 получено достаточное условие на начальные функции ${{u}_{0}}(x)$ и ${{u}_{1}}(x)$ в начально-краевой задаче (1), при которых время разрушения конечно. Для получения достаточных условий разрушения решения соответствующей первой смешанной краевой задачи используется модифицированный метод Х. Левина, изложенный в работе [2]. Также в этом разделе получена аналитическая оценка сверху на время существования решения. В разд. 2 подробно описаны методы, которые позволяют, используя аналитически полученную априорную информацию о времени разрушения решения, численно диагностировать факт разрушения решения и уточнить его локализацию как во времени, так и в пространстве.

1. АНАЛИТИЧЕСКАЯ ДИАГНОСТИКА РАЗРУШЕНИЯ РЕШЕНИЯ

Пусть $u(x,t) \in {{\mathbb{C}}^{{(2)}}}([0,T];{{\mathbb{C}}^{{(2)}}}[0,l])$ – это классическое решение задачи (1). Умножим обе части первого уравнения из (1) на решение $u(x,t)$ и проинтегрируем обе части равенства по $x \in [0,l]$. Тогда с учетом граничных условий из задачи (1) получим

(2)

$\frac{1}{2}\Phi {\text{''}} - J + \left\| {{{u}_{x}}} \right\|_{2}^{2} - \left\| u \right\|_{4}^{4} = \int\limits_0^l {{{u}^{3}}} (x,t)dx,$

где

$\Phi (t) = \int\limits_0^l {{{u}^{2}}} (x,t)dx,\quad J(t) = \int\limits_0^l {u_{t}^{2}} (x,t)dx,\quad {{\left\| u \right\|}_{p}} = \mathop {\left( {\int\limits_0^l {{{{\left| u \right|}}^{p}}} dx} \right)}\nolimits^{1/p} .$

Теперь умножим обе части уравнения из (1) на функцию ${{u}_{t}}(x,t)$ и проинтегрируем по $x \in [0,l]$. В результате получим

(3)

$\frac{d}{{dt}}\left( {\frac{1}{2}J(t) + \frac{1}{2}\left\| {{{u}_{x}}} \right\|_{2}^{2} + \frac{1}{2}\left\| u \right\|_{2}^{2} - \frac{1}{4}\left\| u \right\|_{4}^{4}} \right) = \frac{1}{3}\frac{d}{{dt}}\int\limits_0^l {{{u}^{3}}} (x,t)dx.$

Далее проинтегрируем обе части равенства (3) по времени и с учетом начальных условий из (1) получим равенство

(4)

$\frac{1}{2}J(t) + \frac{1}{2}\left\| {{{u}_{x}}} \right\|_{2}^{2} + \frac{1}{2}\left\| u \right\|_{2}^{2} - \frac{1}{4}\left\| u \right\|_{4}^{4} - E(0) = \frac{1}{3}\int\limits_0^l {{{u}^{3}}} (x,t)dx,$

где

(5)

$E(0) = \frac{1}{2}\int\limits_0^l {u_{1}^{2}} (x)dx + \frac{1}{2}\int\limits_0^l {[{{{({{u}_{{{{0}_{x}}}}})}}^{2}}(x) + u_{0}^{2}(x)]dx} - \frac{1}{4}\int\limits_0^l {u_{0}^{4}} (x)dx - \frac{1}{3}\int\limits_0^l {u_{0}^{3}} (x)dx.$

Из равенств (2) и (4) вытекает равенство

$\frac{1}{2}\Phi {\text{''}}(t) - \frac{5}{2}J(t) - \frac{1}{2}\left[ {\left\| {{{u}_{x}}} \right\|_{2}^{2} + \left\| u \right\|_{2}^{2}} \right] - \frac{1}{4}\left\| u \right\|_{4}^{4} + 3E(0) = 0,$

из которого приходим к неравенству

(6)

$\frac{1}{2}\Phi {\text{''}}(t) + 3E(0) \geqslant \frac{5}{2}J(t).$

Учтем, что из неравенства Коши–Буняковского следует, что

(7)

${{(\Phi {\text{'}}(t))}^{2}} \leqslant 4J(t)\Phi (t).$

И, далее, из неравенств (6) и (7) вытекает неравенство

$\Phi (t)\Phi {\text{''}}(t) - \frac{5}{4}{{(\Phi {\text{'}}(t))}^{2}} + 6E(0)\Phi (t) \geqslant 0.$

Математический интерес представляет случай $E(0) > 0$ (см. формулу (5)).

Далее заметим, что

$\Phi (0) = \int\limits_0^l {u_{0}^{2}} (x)dx,\quad \Phi {\text{'}}(0) = 2\int\limits_0^l {{{u}_{0}}} (x){{u}_{1}}(x)dx.$

Тогда при выполнении условий

${{(\Phi {\text{'}}(0))}^{2}} > 8E(0)\Phi (0) > 0$

из результатов работы [2] на интервале $t \in [0,T]$ выполняется следующее неравенство:

(8)

$\Phi (t) \geqslant {{(\Phi {{(0)}^{{ - 1/4}}} - At)}^{{ - 1/4}}}.$

Имеем

(9)

$T \leqslant {{T}_{0}} = \Phi {{(0)}^{{ - 1/4}}}{{A}^{{ - 1}}},$

где

${{A}^{2}} = \frac{1}{{16}}\Phi {{(0)}^{{ - 5/2}}}[{{(\Phi {\text{'}}(0))}^{2}} - 8E(0)\Phi (0)].$

Таким образом, доказана

Теорема 1. Пусть $u(x,t) \in {{\mathbb{C}}^{{(2)}}}([0,T];{{\mathbb{C}}^{{(2)}}}[0,l])$ – классическое решение задачи Коши (1). Тогда для любых начальных функций ${{u}_{0}}(x),{{u}_{1}}(x) \in {{\mathbb{C}}^{{(2)}}}[0,l]$ таких, что выполнены неравенства

$\Phi {\text{'}}(0) = 2\int\limits_0^l {{{u}_{1}}} (x){{u}_{0}}(x)dx > 0,\quad \Phi (0) = \int\limits_0^l {u_{0}^{2}} (x)dx > 0,\quad E(0) > 0,$

где величина $E(0)$ определена равенством (5), справедлива оценка снизу (8). Более того, справедлива оценка сверху (9). Следовательно, $T < + \infty $, и поэтому классическое решение задачи (1) не существует глобально во времени.

2. ЧИСЛЕННАЯ ДИАГНОСТИКА РАЗРУШЕНИЯ РЕШЕНИЯ

Детально обсудим методы, которые помогут нам численно диагностировать факт разрушения решения и уточнить его локализацию как во времени, так и в пространстве. Напомним, что априорная информация (9), полученная аналитически в теоретической части работы, дает нам верхнюю оценку ${{T}_{0}}$ времени разрушения решения, однако не дает точного времени разрушения ${{T}_{{{\text{bl}}}}}$ (индекс “bl” – сокращение от “blow-up”, “разрушение”) и детального описания процесса разрушения. Численный же подход, использующий аналитически полученную априорную информацию, может помочь детализировать процесс разрушения и уточнить момент разрушения.

Будем использовать обозначения ${{u}_{{{\text{init}}\;{\text{0}}}}}(x)$ и ${{u}_{{{\text{init}}\;1}}}(x)$ вместо ${{u}_{0}}(x)$ и ${{u}_{1}}(x)$ в связи с тем, что мы используем подобные индексы для определения сеточных значений функций.

Для начала сведем исходную начально-краевую задачу (1) к системе первого порядка по времени, который необходим для применения эффективных численных методов, которые будут описаны ниже:

(10)

$\begin{gathered} \frac{{\partial u}}{{\partial t}} = {v},\quad x \in (0,l),\quad t \in (0,T], \hfill \\ \frac{{\partial {v}}}{{\partial t}} - \frac{{{{\partial }^{2}}u}}{{\partial {{x}^{2}}}} + u = {{u}^{2}} + {{u}^{3}},\quad x \in (0,l),\quad t \in (0,T], \hfill \\ u(0,t) = u(l,t) = 0,\quad t \in (0,T], \hfill \\ u(x,0) = {{u}_{{{\text{init}}\;{\text{0}}}}}(x),\quad {{u}_{{{\text{init}}}}}(x,0) = {{u}_{{{\text{init}}\;{\text{1}}}}}(x),\quad x \in [0,l]. \hfill \\ \end{gathered} $

Замечание. Мы ставим задачу численного нахождения решения до момента времени $T = {{T}_{0}}$ включительно, хотя знаем, что решение в этот момент времени не существует и, более того, даже может разрушиться ранее. В связи с этим может сложиться впечатление, что мы ставим заведомо неразрешимую задачу: построить несуществующее решение. Однако это не так. Поясним основную идею работы. Мы имеем аналитическую оценку времени разрушения ${{T}_{0}}$ сверху (9). Реальное время разрушения решения ${{T}_{{{\text{bl}}}}}$ мы не знаем. Поэтому мы вынуждены проводить расчеты до момента времени ${{T}_{0}}$ включительно. Численное решение, благодаря используемой ниже численной схеме, будет существовать при всех $t \in [0,{{T}_{0}}]$, хотя в какой-то момент времени ${{T}_{{{\text{bl}}}}}$ реальное решение разрушается. Возникает вопрос, при каких значениях $t$ можно доверять результату численных расчетов, а при каких нет. Этому и посвящена дальнейшая часть работы.

2.1. Поиск численного решения

C целью численного решения системы (10) мы применим жесткий метод прямых (SMOL) [3], [4] для того, чтобы свести исходную систему уравнений с двумя уравнениями в частных производных к системе обыкновенных дифференциальных уравнений.

Сначала мы введем равномерную сетку ${{X}_{N}}$ только по пространственной переменной $x$ с шагом $h = (l - 0){\text{/}}N$ и содержащую $N + 1$ узлов (что соответствует $N$ интервалам): ${{X}_{N}} = \left\{ {{{x}_{n}},\;0 \leqslant n \leqslant N:{{x}_{n}} = 0 + nh} \right\}$. Таким образом, после конечно-разностной аппроксимации пространственных производных со вторым порядком точности мы получим следующую дифференциально-алгебраическую систему, из которой требуется определить $N - 1$ неизвестную функцию ${{u}_{n}} \equiv {{u}_{n}}(t) \equiv u({{x}_{n}},t)$ ($n = \overline {1,N - 1} $, ${{u}_{0}}$ и ${{u}_{N}}$ определяются из соотношений, заданных граничными условиями) и $N - 1$ вспомогательную функцию ${{{v}}_{n}} \equiv {{{v}}_{n}}(t) \equiv {v}({{x}_{n}},t)$ ($n = \overline {1,N - 1} $, ${{{v}}_{0}}$ и ${{{v}}_{N}}$ не входят в систему):

$\begin{gathered} \frac{{d{{u}_{n}}}}{{dt}} = {{{v}}_{n}},\quad n = \overline {1,N - 1} ,\quad t \in (0,T], \hfill \\ \frac{{d{{{v}}_{n}}}}{{dt}} = \frac{{{{u}_{{n + 1}}} - 2{{u}_{n}} + {{u}_{{n - 1}}}}}{{{{h}^{2}}}} - {{u}_{n}} + u_{n}^{2} + u_{n}^{3},\quad n = \overline {1,N - 1} ,\quad t \in (0,T], \hfill \\ {{u}_{0}} = {{u}_{N}} = 0,\quad t \in (0,T], \hfill \\ {{u}_{n}}(0) = {{u}_{{{\text{init}}\;0}}}({{x}_{n}}),\quad {{{v}}_{n}}(0) = {{u}_{{{\text{init}}\;{\text{1}}}}}({{x}_{n}}),\quad n = \overline {0,N} . \hfill \\ \end{gathered} $

Эта система может быть сведена к чисто дифференциальной и переписана в следующем виде:

(11)

$\begin{gathered} \frac{{d{\mathbf{y}}}}{{dt}} = {\mathbf{f}}({\mathbf{y}}),\quad t \in (0,T], \\ {\mathbf{y}}(0) = {{{\mathbf{y}}}_{{{\text{init}}}}}, \\ \end{gathered} $

где ${\mathbf{y}} = {{({{u}_{1}}{{u}_{2}}\; \ldots \;{{u}_{{N - 1}}}{{{v}}_{1}}{{{v}}_{2}}\; \ldots \;{{{v}}_{{N - 1}}})}^{{\text{т}}}}$, ${\mathbf{f}} = {{({{f}_{1}}{{f}_{2}}\; \ldots \;{{f}_{{2N - 2}}})}^{{\text{т}}}}$ и ${{{\mathbf{y}}}_{{{\text{init}}}}} = {{({{u}_{{{\text{init}}\;{\text{0}}}}}({{x}_{1}})\mathop {{{u}_{{{\text{init}}}}}}\nolimits_0 ({{x}_{2}})\; \ldots \;{{u}_{{{\text{init}}\;{\text{0}}}}}({{x}_{{N - 1}}}){{u}_{{{\text{init}}\;1}}}({{x}_{1}}){{u}_{{{\text{init}}\;1}}}({{x}_{2}})\; \ldots \;{{u}_{{{\text{init}}\;1}}}({{x}_{{N - 1}}}))}^{{\text{т}}}}$ ${{{\mathbf{y}}}_{{{\text{init}}}}} = {{({{u}_{{{\text{init}}\;{\text{0}}}}}({{x}_{1}})\mathop {{{u}_{{{\text{init}}}}}}\nolimits_0 ({{x}_{2}})\; \ldots \;{{u}_{{{\text{init}}\;{\text{0}}}}}({{x}_{{N - 1}}}){{u}_{{{\text{init}}\;1}}}({{x}_{1}}){{u}_{{{\text{init}}\;1}}}({{x}_{2}})\; \ldots \;{{u}_{{{\text{init}}\;1}}}({{x}_{{N - 1}}}))}^{{\text{т}}}}$. Здесь вектор-функция ${\mathbf{f}}$ имеет следующую структуру:

${{f}_{n}} = \left\{ \begin{gathered} {{y}_{{n + N - 1}}},\quad \quad \quad \quad \quad \,\;\,\quad \quad \quad \quad \,\;\,\quad \quad \quad \quad \,\;\,\;\;\,\,{\text{если}}\quad n = \overline {1,\;N - 1} , \hfill \\ \tfrac{1}{{{{h}^{2}}}}({{y}_{2}} - 2{{y}_{1}} + 0) - {{y}_{1}} + y_{1}^{2} + y_{1}^{3},\quad \quad \quad \quad \quad \,\;\,\quad \quad \;\,\,\,{\text{если}}\quad n = N, \hfill \\ \tfrac{1}{{{{h}^{2}}}}({{y}_{{n - N + 2}}} - 2{{y}_{{n - N + 1}}} + {{y}_{{n - N}}}) - {{y}_{{n - N + 1}}} + y_{{n - N + 1}}^{2} + y_{{n - N + 1}}^{3},\quad {\text{если}}\quad n = \overline {N + 1,\;2N - 3} , \hfill \\ \tfrac{1}{{{{h}^{2}}}}(0 - 2{{y}_{{N - 1}}} + {{y}_{{N - 2}}}) - {{y}_{{N - 1}}} + y_{{N - 1}}^{2} + y_{{N - 1}}^{3},\quad \quad \quad \quad \,\;\,\,\,{\text{если}}\quad n = 2N - 2. \hfill \\ \end{gathered} \right.$

Для численного решения системы (11) мы будем использовать одностадийную схему Розенброка с комплексным коэффициентом CROS1 [5], [6], которая является наилучшим выбором для решения задач такого рода по той причине, что она не приводит к переполнению даже в том случае, если решение задачи устремляется к бесконечности [7], [8]. Для ее применения мы введем равномерную сетку ${{T}_{M}}$ по времени $t$ с шагом $\tau = (T - 0){\text{/}}M$, которая имеет $M + 1$ узлов (т.е. $M$ интервалов): ${{T}_{M}} = \left\{ {{{t}_{m}},\;0 \leqslant m \leqslant M:{{t}_{m}} = 0 + m\tau } \right\}$. После этого запишем схему CROS1 для решения системы (11):

(12)

$\begin{gathered} {\mathbf{y}}({{t}_{{m + 1}}}) = {\mathbf{y}}({{t}_{m}}) + ({{t}_{{m + 1}}} - {{t}_{m}})\operatorname{Re} {\mathbf{w}}, \hfill \\ {\text{где}}\;{\mathbf{w}}\;{\text{является}}\;{\text{решением}}\;{\text{СЛАУ,}} \hfill \\ \left[ {{\mathbf{E}} - \frac{{1 + i}}{2}({{t}_{{m + 1}}} - {{t}_{m}}){{{\mathbf{f}}}_{{\mathbf{y}}}}({\mathbf{y}}({{t}_{m}}))} \right]{\mathbf{w}} = {\mathbf{f}}({\mathbf{y}}({{t}_{m}})). \hfill \\ \end{gathered} $

Здесь ${\mathbf{E}}$ – единичная матрица, а ${{{\mathbf{f}}}_{{\mathbf{y}}}}$ – матрица Якоби.

2.2. Численная диагностика разрушения решения

При численных расчетах важно не только получить приближенный численный результат, но также и выполнить некоторую оценку его точности. Метод вычисления апостериорной асимптотически точной оценки погрешности [7] позволяет это сделать. Но этот метод также может помочь и диагностировать факт разрушения точного решения [8]. Основные формулы и утверждения этого пункта впервые были представлены в работах [7]–[9]. Мы развили эти идеи в работах [10]–[17].

Для начала введем базовую сетку ${{X}_{N}} \times {{T}_{M}}$: ${\text{\{ }}{{x}_{n}},{{t}_{m}}{\text{\} }}$, $0 \leqslant n \leqslant N$, $0 \leqslant m \leqslant M$. После этого произведем последовательное сгущение сеток, начиная с базовой, и вычислим набор сеточных значений решения

${{u}_{{(s)}}}(x,t) \equiv {{u}^{{(r_{x}^{{s - 1}}N,r_{t}^{{s - 1}}M)}}}(x,t)$

на полученных сетках ${{X}_{{r_{x}^{{s - 1}}N}}} \times {{T}_{{r_{t}^{{s - 1}}M}}}$ ($s$ – номер сетки из набора $S$ сеток, $s = \overline {1,S} $).

В связи с тем, что мы аппроксимировали все пространственные производные в (10) с точностью $O({{h}^{2}})$, а при численном интегрировании системы (11) используем схему CROS1 (12), которая имеет точность $O({{\tau }^{2}})$, построенный метод решения системы (10) имеет точность $O({{\tau }^{2}} + {{h}^{2}})$ или, другими словами, теоретический порядок точности по времени ${{p}_{t}} = p_{t}^{{{\text{theor}}}} \equiv 2$, а по пространству ${{p}_{x}} = p_{x}^{{{\text{theor}}}} \equiv 2$. Мы выполняем последовательное сгущение сеток по времени в целое число раз ${{r}_{t}}$ и сгущение пространственной сетки в целое число раз ${{r}_{x}}$ так, чтобы выполнялось условие $r_{t}^{{{{p}_{t}}}} = r_{x}^{{{{p}_{x}}}}$ (подробности см. в [9]), т.е. в нашем случае $r_{t}^{2} = r_{x}^{2}$. Наиболее удобно для счета выбрать ${{r}_{t}} = {{r}_{x}} \equiv 2$. В этом случае каждая последующая сетка ${{X}_{{r_{x}^{{s - 1}}N}}} \times {{T}_{{r_{t}^{{s - 1}}M}}}$ имеет узлы, совпадающие с узлами $({{x}_{n}},{{t}_{m}})$ базовой сетки. В этих узлах $(x,t)$ мы можем выполнить апостериорную асимптотически точную оценку погрешности [8], [9]

(13)

${{\Delta }_{{(s)}}}({{x}_{n}},{{t}_{m}}) = \frac{{{{u}_{{(s)}}}({{x}_{n}},{{t}_{m}}) - {{u}_{{(s - 1)}}}({{x}_{n}},{{t}_{m}})}}{{r_{t}^{2} - 1}},$

и оценить эффективный порядок точности [8], [9] на всем промежутке времени $t \in [0,{{T}_{0}}]$

$\mathop p\nolimits_{{{t}_{{(s)}}}}^{{\text{eff}}} = lo{{g}_{{{{r}_{t}}}}}\frac{{\sqrt {\sum\limits_{n = 0}^N {\sum\limits_{m = 0}^M {{{{\left( {{{u}_{{(s - 1)}}}({{x}_{n}},{{t}_{m}}) - {{u}_{{(s - 2)}}}({{x}_{n}},{{t}_{m}})} \right)}}^{2}}} } } }}{{\sqrt {\sum\limits_{n = 0}^N {\sum\limits_{m = 0}^M {{{{\left( {{{u}_{{(s)}}}({{x}_{n}},{{t}_{m}}) - {{u}_{{(s - 1)}}}({{x}_{n}},{{t}_{m}})} \right)}}^{2}}} } } }}$

(здесь выбрана наиболее удобная из возможных норм – евклидова). В случае, если на всем промежутке времени $t \in [0,{{T}_{0}}]$ решение задачи имеет непрерывные производные порядка ${{p}_{x}}$ по $x$ и ${{p}_{t}}$ по $t$, имеет место сходимость [8]

$\mathop p\nolimits_{{{t}_{{(s)}}}}^{{\text{eff}}} \;\xrightarrow[{s \to \infty }]{}\;p_{t}^{{{\text{theor}}}} \equiv {{p}_{t}}.$

Нарушение этой сходимости говорит о потере гладкости точного решения на промежутке времени $t \in [t,{{T}_{0}}]$.

С целью локализации конкретного момента времени, в который произошло разрушение решения, можно оценить эффективный порядок точности и поточечно в каждом узле ${{t}_{m}} \in {{T}_{M}}$, $0 < m \leqslant M$,

(14)

$\mathop p\nolimits_{{{t}_{{(s)}}}}^{{\text{eff}}} ({{t}_{m}}) = lo{{g}_{{{{r}_{t}}}}}\frac{{\sqrt {\sum\limits_{n = 0}^N {{{{\left( {{{u}_{{(s - 1)}}}({{x}_{n}},{{t}_{m}}) - {{u}_{{(s - 2)}}}({{x}_{n}},{{t}_{m}})} \right)}}^{2}}} } }}{{\sqrt {\sum\limits_{n = 0}^N {{{{\left( {{{u}_{{(s)}}}({{x}_{n}},{{t}_{m}}) - {{u}_{{(s - 1)}}}({{x}_{n}},{{t}_{m}})} \right)}}^{2}}} } }}.$

В точках $t$, в которых решение исходной задачи имеет непрерывные производные порядка ${{p}_{t}}$ по времени и порядка ${{p}_{x}}$ по пространству, имеет место сходимость

$\mathop p\nolimits_{{{t}_{{(s)}}}}^{{\text{eff}}} (t)\;\xrightarrow[{s \to \infty }]{}\;p_{t}^{{{\text{theor}}}} \equiv {{p}_{t}}$

и соответствующая оценка погрешности (13) является асимптотически точной при $s \to \infty $ (или, что то же самое, $N,M \to \infty $). Нарушение этой сходимости говорит о потере гладкости точного решения [7]–[9]. В частности, в случае степенной “сингулярности” $u(x,t) \sim {{({{T}_{{{\text{bl}}}}} - t)}^{{ - \beta }}}$, где ${{T}_{{{\text{bl}}}}}$ – время разрушения (напомним, что индекс “bl” – сокращение от “blow-up”, “разрушение”), для любого $t > {{T}_{{{\text{bl}}}}}$ эффективный порядок точности $\mathop p\nolimits_{{{t}_{{(s}}}_{)}}^{{\text{eff}}} (t)\;\xrightarrow[{s \to \infty }]{}\; - {\kern 1pt} \beta $. Это позволяет нам найти соответствующую степень $\beta $. Если $\mathop p\nolimits_{{{t}_{{(s)}}}}^{{\text{eff}}} (t)\;\xrightarrow[{s \to \infty }]{}\; - {\kern 1pt} \infty $ для любого $t > {{T}_{{{\text{bl}}}}}$, мы можем утверждать, что решение экспоненциально возрастает, т.е. $u(x,t) = \infty $; если $\mathop p\nolimits_{{{t}_{{(s)}}}}^{{\text{eff}}} (t)\;\xrightarrow[{s \to \infty }]{}\;0$ для любого $t > {{T}_{{{\text{bl}}}}}$, то рост решения в окрестности “сингулярности” является логарифмическим: $u(x,t) \sim ln({{T}_{{{\text{bl}}}}} - t)$. Вывод соответствующих утверждений можно найти в [7], [8]. Момент разрушения решения ${{T}_{{{\text{bl}}}}}$ может быть найден с точностью до величины интервала базовой сетки по времени ${{T}_{M}}$.

После локализации момента времени разрушения $t$ можно локализовать и пространственные точки разрушения решения по переменной $x$ для каждого конкретного момента времени $t$ в каждом узле ${{t}_{m}} \in {{T}_{M}}$, $0 \leqslant m \leqslant M$,

$\mathop p\nolimits_{x{{t}_{{(s)}}}}^{{\text{eff}}} ({{x}_{n}},{{t}_{m}}) = lo{{g}_{{{{r}_{t}}}}}\frac{{\left| {{{u}_{{(s - 1)}}}({{x}_{n}},{{t}_{m}}) - {{u}_{{(s - 2)}}}({{x}_{n}},{{t}_{m}})} \right|}}{{\left| {{{u}_{{(s)}}}({{x}_{n}},{{t}_{m}}) - {{u}_{{(s - 1)}}}({{x}_{n}},{{t}_{m}})} \right|}}.$

Если нарушение гладкости решения возникает во всей области по пространственной переменной одновременно, то отклонение сходимости $\mathop p\nolimits_{x{{t}_{{(s)}}}}^{{\text{eff}}} (x,t)$ от 2 возникает во всех точках сетки ${\text{\{ }}{{x}_{n}}{\text{\} }}$ с первого временного слоя ${{t}_{m}} \geqslant {{T}_{{{\text{bl}}}}}$ (см., например, пример 1a в [11]). Если разрушение решения возникает в одной-единственной точке $x{\kern 1pt} {\text{*}}$, то описанный метод позволяет оценить локализацию этой точки (см., пример 1 ниже). Такая диагностика процесса разрушения решения возможна в связи с тем, что схема CROS1 не приводит к переполнению даже в том случае, если решение задачи устремляется к бесконечности [7], [8].

Пример 1. Для начала рассмотрим пример для следующего набора входных данных:

(16)

${{u}_{{{\text{init}}\;{\text{0}}}}}(x) \equiv asin2x,\quad {{u}_{{{\text{init}}\,{\text{1}}}}}(x) \equiv bsin2x,\quad l = \pi ,$

(17)

$a = 7,\quad b = 10.$

К сожалению, точное решение уравнения (1) не может быть получено аналитически, однако в теоретической части работы была получена априорная оценка сверху ${{T}_{0}}$ (9) для времени разрушения решения ${{T}_{{{\text{bl}}}}}$. В случае набора параметров (16) мы получим

(18)

${{T}_{{{\text{bl}}}}} \leqslant {{T}_{0}} = {{T}_{0}}(a,b) \equiv \frac{{2a}}{{\sqrt {\tfrac{{3{{a}^{4}}}}{8} - 5{{a}^{2}}} }},$

что дает для (17): ${{T}_{0}} \simeq 0.58$.

Применим численный алгоритм диагностики момента и места разрушения решения. Для численного решения задачи (1) с входными данными (16), (17) мы возьмем следующий набор параметров: ${{T}_{0}} = 0.581$, $N = 50$, $M = 50$, ${{r}_{x}} = 2$, ${{r}_{t}} = 2$, $S = 7$ (число последовательно используемых для вычислений сеток, включая начальную). Численное решение ${{u}^{{(r_{x}^{{S - 1}}N,r_{t}^{{S - 1}}M)}}}$ представлено на фиг. 1 (отмечены только узлы, совпадающие с узлами базовой сетки ${{X}_{N}} \times {{T}_{M}}$).

Фиг. 1.

Пример 1: численное решение задачи в разные моменты времени. Красным цветом отмечено численное решение в те моменты времени, при которых ему нельзя доверять.

Получив приближенное численное решение на разных сетках, мы можем проверить сходимость эффективного порядка точности к теоретическому для каждого временного слоя по формуле (14). После вычислений на $S$ вложенных сетках эффективный порядок точности $\mathop p\nolimits_{{{t}_{{(s)}}}}^{{\text{eff}}} $ для каждого временного слоя ${{t}_{m}}$ сходится к ${{p}^{{{\text{theor}}}}} = 2$ (см. фиг. 2) за исключением временных слоев, соответствующих моментам времени ${{t}_{m}}$, $m \geqslant 21$: ${{p}^{{{\text{eff}}}}}({{t}_{{21}}}) \to - 1$, что означает, что ${{T}_{{{\text{bl}}}}} \in ({{t}_{{20}}},{{t}_{{21}}}] \equiv (0.232,\;0.242]$ является временем разрушения решения и что в точке ${{T}_{{{\text{bl}}}}}$ решение имеет особенность типа полюса $u(x,t) \sim {{({{T}_{{{\text{bl}}}}} - t)}^{{ - 1}}}$. Таким образом, численный алгоритм позволил уточнить аналитическую оценку времени разрушения.

Фиг. 2.

Пример 1: эффективный порядок точности численной схемы в каждый момент времени. Разрушение решения диагностировано в момент времени ${{T}_{{{\text{bl}}}}} \in ({{t}_{{20}}},{{t}_{{21}}}] \equiv (0.232,\;0.242]$.

Мы также можем оценить эффективный порядок точности для каждой пространственной точки выбранного временного слоя, соответствующего моменту времени ${{t}_{m}}$, по формуле (14). Мы можем использовать эту формулу, например, для временного слоя, соответствующего моменту времени ${{t}_{{21}}}$, в котором было диагностировано разрушение решения изначально, с целью определить, возникло ли разрушение на всем временном слое или только в отдельных точках по пространственной переменной. Исходя из фиг. 3 мы можем предположить, что решение разрушилось сначала только в некоторой подобласти рассматриваемой области $x \in [0,\pi ]$.

Фиг. 3.

Пример 1: эффективный порядок точности в каждой точке по пространству в моменты времени ${{t}_{{20}}}$ и ${{t}_{{21}}}$. Разрушение решения диагностировано на некоторой подобласти области $x \in [0,\pi ]$.

Также лишний раз отметим, что данная методика позволяет сделать вывод о том, при каких ${{t}_{m}}$, $0 \leqslant m \leqslant M$, мы можем доверять численному решению, а при каких нет (см. фиг. 1).

Пример 2. В данном численном эксперименте мы рассматриваем серию задач из примера 1 с набором входных данных (16) для различных значений параметров $a$ и $b$. В этом случае по-прежнему верна аналитическая оценка времени разрушения (18) (отметим, что она не зависит от параметра $b$), из которой следует, что условия применимости теоремы 1 выполнены при $\tfrac{{3{{a}^{4}}}}{8} - 5{{a}^{2}} > 0$, т.е. при $a > \sqrt {\tfrac{{40}}{3}} \simeq 3.65$. При меньших значениях параметра $a$ аналитический ответ на вопрос о наличии факта разрушения решения остается открытым. Однако численные методы могут помочь дополнить информацию, полученную аналитически. На фиг. 4 представлены графики зависимости времени разрушения ${{T}_{{{\text{bl}}}}}$ (оцененного численно) решения задачи (1) с набором входных данных (16) от амплитуды $b$ начальной функции ${{u}_{{{\text{init}}\;{\text{1}}}}}(x)$ при фиксированных значениях параметра $a$ ($a = 0.25$ и $a = 0.5$).

Фиг. 4.

Пример 3: график зависимости времени разрушения ${{T}_{{{\text{bl}}}}}$ решения задачи (1) с набором входных данных (16) от амплитуды $b$ начальной функции ${{u}_{{{\text{init}}\;{\text{1}}}}}(x)$ при фиксированных значениях параметра $a$.

Из характера соответствующих графиков можно предположить, что они имеют вертикальные асимптоты, в окрестностях которых реальное время разрушения ${{T}_{{{\text{bl}}}}}$ уходит на бесконечность. Это означает, что можно оценить предполагаемую область глобальной разрешимости. Более того, на фиг. 5 представлена зависимость времени разрушения решения задачи (1) с набором входных данных (16) от амплитуд $a$ и $b$, $(a,b) \in [0,3] \times [0,5]$, начальных функций. В результате мы можем оценить область определения параметров $a$ и $b$, в которой поставленная задача предположительно является глобально разрешимой по времени. Отметим, что этот результат не удалось получить аналитическими методами.

Фиг. 5.

Пример 2: график зависимости времени разрушения ${{T}_{{{\text{bl}}}}}$ решения задачи (1) с набором входных данных (16) от амплитуд $a$ и $b$ начальных функций (16).

ЗАКЛЮЧЕНИЕ

В данной работе продемонстрирована возможность эффективного комбинирования аналитического и численного подходов при диагностике факта разрушения численного решения нелинейного уравнения Клейна–Города, которое описывает механику кристаллов. Комбинация этих подходов позволяет получить качественную информацию о факте и характере разрушения приближенного численного решения.

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

Li Kaitaia, Zhang Quande. Existence and nonexistence of global solutions for the equation of dislocation // J. Differential Equations. 1998. V. 10. № 1. P. 5–21.
Корпусов М.O. Разрушение в нелинейных волновых уравнениях с положительной энергией // URSS, 2012, 256 с.
Hairer E., Wanner G. Solving of Ordinary Differential Equations. Stiff and Differential-Algebraic Problems. Springer, 2002.
Kalitkin N.N. Numerical methods for solving stiff systems [in Russian] // Matematicheskoye modelirovanie. 1995. V. 7. № 5. P. 8–11.
Rosenbrock H.H. Some general implicit processes for the numerical solution of differential equations // Computer Journal. 1963. V. 5. № 4. P. 329–330.
Al’shinn A.B., Al’shina E.A., Kalitkin N.N., Koryagina A.B. Rosenbrock schemes with complex coefficients for stiff and differential algebraic systems // Comp. Math. and Math. Phys. 2006. V. 46. № 8. P. 1320–1340.
Al’shin A.B., Kalitkin N.N., Koryakin P.V. Diagnostics of singularities of exact solutions in computations with error control // Comp. Math. and Math. Phys. 2005. V. 45. № 19. P. 1769–1779.
Al’shin A.B., Al’shina E.A. Numerical diagnosis of blow-up of solutions of pseudoparabolic equations // J. Math. Sci. 2008. V. 148. № 1. P. 143–162.
Kalitkin N.N., Al’shin A.B., Al’shina E.A., Rogov B.V. Computations on Quasi-Uniform Grids [in Russian]. M.: Fizmatlit. 2005.
Korpusov M.O., Lukyanenko D.V., Panin A.A., Yushkov E.V. Blow-up for one sobolev problem: theoretical approach and numerical analysis // J. Math. Analysis and Applicat. 2016. V. 442. № 52. P. 451–468.
Korpusov M.O., Lukyanenko D.V., Panin A.A., Yushkov E.V. Blow-up phenomena in the model of a space charge stratification in semiconductors: analytical and numerical analysis // Math. Methods in the Appl. Sci. 2017. V. 40. № 7. P. 2336–2346.
Лукьяненко Д.В., Панин А.А. Разрушение решения уравнения стратификации объемного заряда в полупроводниках: численный анализ при сведении исходного уравнения к дифференциально-алгебраической системе // Вычисл. методы и программирование: Новые вычисл. технологии (Электронный научный журнал). 2016. Т. 17. № 1. С. 437–446.
Корпусов М.О., Лукьяненко Д.В., Овсянников Е.А., Панин А.А. Локальная разрешимость и разрушение решения одного уравнения с квадратичной некоэрцитивной нелинейностью // Вест. Южно-Уральского гос. университета. Сер. Матем. моделирование и программирование. 2017. Т. 10. № 2. С. 107–123.
Korpusov M.O., Lukyanenko D.V. Instantaneous blow-up versus local solvability for one problem of propagation of nonlinear waves in semiconductors // J. Math, Analy, and Applicat. 2018. V. 459. № 1. P. 159–181.
Korpusov M.O., Lukyanenko D.V., Panin A.A., Shlyapugin G.I. On the blow-up phenomena for a one-dimensional equation of ion-sound waves in a plasma: analytical and numerical investigation // Math. Methods in the Appl. Sci. 2018. V. 41. № 8. P. 2906–2929.
Korpusov M.O., Lukyanenko D.V., Panin A.A., Yushkov E.V. On the blow-up of solutions of a full non-linear equation that describes ion-sound waves in plasma with non-coercive non-linearities // Izvestiya. Mathematics. 2018. V. 82. № 2. P. 283–317.
Колотов И.И., Панин А.А. О непродолжаемых решениях и разрушении решений псевдопараболических уравнений с коэрцитивной и знакопостоянной нелинейностями // Матем. заметки. 2019. Т. 105. № 5. С. 708–723.

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

Инструменты

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

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

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