Программирование, 2023, № 2, стр. 5-12

РЕШЕНИЕ ЗАДАЧИ КОШИ ДЛЯ ТРЕХМЕРНОГО РАЗНОСТНОГО УРАВНЕНИЯ В ПАРАЛЛЕЛЕПИПЕДЕ

М. С. Апанович a*, А. П. Ляпин bc**, К. В. Шадрин a***

a Красноярский государственный медицинский университет имени профессора В.Ф. Войно-Ясенецкого
660022 Красноярск, ул. Партизана Железняка, 1, Россия

b Сибирский федеральный университет
660041 Красноярск, пр. Свободный, 79, Россия

c Фэрмонтский государственный университет
660041 Западная Вирджиния, Фэрмонт, Локуст стрит, 1201, США

* E-mail: marina.apanovich@list.ru
** E-mail: aplyapin@sfu-kras.ru
*** E-mail: kvsh_buffon@mail.ru

Поступила в редакцию 17.07.2022
После доработки 10.08.2022
Принята к публикации 21.10.2022

Полный текст (PDF)

Аннотация

Представлен алгоритм вычисления решения задачи Коши для трехмерного разностного уравнения с постоянными коэффициентами в точке по коэффициентам разностного уравнения и начальным данным Коши, заданным в параллелепипеде. Наш алгоритм использует теоремы Апанович М.С. и Лейнартаса Е.К., характеризующие корректность задачи Коши, и методы компьютерной алгебры для достижения вычислительной эффективности. Алгоритм реализован в среде Matlab.

1. ВВЕДЕНИЕ

Метод конечно-разностной дискретизации, примененный к дифференциальным и интегральным уравнениям, приводит к разностным уравнениям (см. [1, 2]). Комбинация линейных разностных уравнений и метода производящих функций эффективно применяется к многочисленным задачам перечислительного комбинаторного анализа (см., например, [35]). Для уравнений Коши–Римана это дает дискретную теорию аналитических функций (см. [6, 7]), используемую для изучения Римановых поверхностей (см., например, [8, 9]). Также разностные уравнения используются в динамических моделях с дискретным временем (см. [10, 11]).

Дополнительные условия (“начальные”, “граничные”, “данные Коши”) позволяют выбрать единственное решение из бесконечного множества решений, а возникающая при этом задача называется задачей Коши для многомерного разностного уравнения. Одномерная задача Коши проста, так как ее исходные данные конечны. Многомерный случай отличается тем, что исходные данные заданы на бесконечном множестве и имеют сложную структуру. Алгоритмы решения многомерного разностного уравнения с постоянными, полиномиальными или рациональными коэффициентами функции изучались в [1215]. Связь производящей функции решения линейной задачи Коши и производящей функции ее данных Коши исследована в работах [1619]. Алгоритм вычисления решения двумерной задачи Коши в точке через коэффициенты разностного уравнения и данные Коши были разработаны и реализованы в [20, 21].

Отметим, что многомерные случаи должным образом не исследованы, однако они имеют большое значение, например, трехмерный случай имеет большое значение в задачах термодинамики (теплопроводность, см. [2], анизотропная диффузия, см. [22]), математической биологии (распределение морфогена в крыльях насекомых [22], распределение плотности популяции видов [23]). В данной статье рассматривается решение трехмерной задачи Коши для линейного разностного уравнения с постоянными коэффициентами в параллелепипеде в точке через коэффициенты линейного разностного уравнения и данные Коши. Мы определяем задачу Коши в параллелепипеде и описываем функцию начальных данных таким образом, который позволяет нам разработать и реализовать компьютерный алгоритм.

2. ЗАДАЧА КОШИ ДЛЯ ПОЛИНОМИАЛЬНОГО РАЗНОСТНОГО ОПЕРАТОРА В “ПАРАЛЛЕЛЕПИПЕДЕ”

Пусть функция f целочисленных аргументов $f:\mathbb{Z}_{ + }^{n} \to \mathbb{C}$ и ${{\delta }_{j}}$ – оператор сдвига по $j$-й переменной, т.е. δjf(t1, ..., tn) = f(t1, ..., tj – 1, ${{t}_{j}} + 1,{{t}_{{j + 1}}},...,{{t}_{n}})$, $j = 1,2,...,n$. Некоторые полезные свойства оператора сдвига были исследованы в [24]. Если α = $({{\alpha }_{1}},...,{{\alpha }_{n}})$ – мультииндекс, то $\left| \alpha \right| = {{\alpha }_{1}} + ... + {{\alpha }_{n}}$, ${{\delta }^{\alpha }} = \delta _{1}^{{{{\alpha }_{1}}}}...\delta _{n}^{{{{\alpha }_{n}}}}$. Для двух мультииндексов $\alpha = \left( {{{\alpha }_{1}},...,{{\alpha }_{n}}} \right)$, β = = $({{\beta }_{1}},...,{{\beta }_{n}})$ неравенство $\alpha \leqslant \beta $ означает, что ${{\alpha }_{j}} \leqslant {{\beta }_{j}}$, $j = 1,2,...,n$.

Рассмотрим разностный полиномиальный оператор вида

(1)
$P\left( \delta \right) = \sum\limits_{0 \leqslant \alpha \leqslant m} {{c}_{\alpha }}{{\delta }^{\alpha }},$
где ${{c}_{\alpha }}$ – постоянные коэффициенты оператора $P\left( \delta \right)$. Если $\alpha = \left( {{{\alpha }_{1}},...,{{\alpha }_{n}}} \right)$, то будем обозначать $'\alpha = \left( {{{\alpha }_{1}},...,{{\alpha }_{{n - 1}}}} \right)$ и тогда оператор $P\left( \delta \right)$ можно записать в виде

$P\left( \delta \right) = \sum\limits_{{{\alpha }_{n}} = 0}^{{{m}_{n}}} \left( {\sum\limits_{'0 \leqslant '\alpha \leqslant 'm} {{c}_{{'\alpha ,{{\alpha }_{n}}}}}{{\delta }^{{'\alpha }}}} \right)\delta _{n}^{{{{\alpha }_{n}}}}.$

Многочлен

$\begin{gathered} P\left( {z,{{z}_{n}}} \right) = \mathop {\mathop {\sum }\limits^{{{m}_{n}}} }\limits_{{{\alpha }_{n}} = 0} \left( {\sum\limits_{'0 \leqslant '\alpha \leqslant 'm} {c{{{\kern 1pt} }_{{'\alpha ,{{\alpha }_{n}}}}}'{\kern 1pt} {{z}^{{'\alpha }}}} } \right)z_{n}^{{{{\alpha }_{n}}}} = \\ = \mathop {\mathop {\sum }\limits^{{{m}_{n}}} }\limits_{{{\alpha }_{n}} = 0} {{P}_{j}}\left( z \right)z_{n}^{{{{\alpha }_{n}}}} \\ \end{gathered} $
будем называть характеристическим многочленом для разностного оператора $P\left( \delta \right)$, а его степень по переменной ${{z}_{n}}$ – порядком разностного оператора.

В целочисленной решетке $\mathbb{Z}_{ + }^{n}$ выберем точку $x\, = \,({{x}_{1}},...,{{x}_{{n - 1}}},0)$ и обозначим Π'x = $\{ t \in \mathbb{R}_{ + }^{n}$ : $'0 \leqslant 't \leqslant 'x\} $ $\left( {n - 1} \right)$-мерный параллелепипед в гиперплоскости ${{t}_{n}} = 0$.

Зафиксируем точку $\beta = \left( {{{\beta }_{1}},...,{{\beta }_{{n - 1}}},{{m}_{n}}} \right) \in \mathbb{Z}_{ + }^{n}$, $'0 \leqslant '\beta \leqslant 'm$ такую, что коэффициент ${{c}_{\beta }} \ne 0$, и рассмотрим множество ${{\Pi }_{{'\beta ,'m}}}$ = $\{ 't \in {{\Pi }_{{'x}}}$ : $'\beta \, \leqslant \,'t\, \leqslant \,'x\, - \,'m$ + + 'β}. Обозначим L = $\left( {{{\Pi }_{{'x}}}{{\backslash }}{{\Pi }_{{'\beta ,'m}}}} \right) \times [0,{{x}_{n}}]$ множество, на котором будем задавать начальные данные и сформулируем задачу:

найти решение разностного уравнения

(2)
$P\left( \delta \right)f\left( x \right) = g\left( x \right),\quad x \in \Pi = {{\Pi }_{{'x}}} \times \left[ {0,{{x}_{n}}} \right],$

удовлетворяющее условию

(3)
$f(x) = \varphi (x),\quad x \in L,$
где $g\left( x \right)$ и $\varphi \left( x \right)$ – заданные функции целочисленных аргументов.

Однородная двухслойная линейная разностная схема с постоянными коэффициентами рассматривалась в [25], где исследована ее устойчивость при ${{m}_{n}} = 1$. Разрешимость двумерной задачи Коши (2), (3) исследована в [26]. Такие задачи известны как многослойные неявные разностные схемы. В [27] исследована корректность задачи (2), (3) и доказано достаточное условие корректности, а достаточное условие разрешимости трехмерной задачи (2), (3) доказано в [28].

Рассмотрим задачу Коши (2), (3) для $n = 3$. ${{\mathbb{Z}}^{3}} = \mathbb{Z} \times \mathbb{Z} \times \mathbb{Z}$ – целочисленная решетка, а $\mathbb{Z}_{ + }^{3}$ – подмножество этой решетки, состоящее из точек с целыми неотрицательными координатами. Пусть ${{\delta }_{1}}$, ${{\delta }_{2}}$, ${{\delta }_{3}}$ – операторы сдвига по переменной $x$, y, $z$ соответственно, т.е. δ1f(x, y, z) = = $f(x + 1,y,z)$, ${{\delta }_{2}}f(x,y,z)\, = \,f(x,y + 1,z)$, δ3f(x, y, z) = = $f(x,y,z + 1)$. Зададим параллелепипед Π = = $\{ (x,y,z) \in \mathbb{Z}_{ + }^{3},$ $0 \leqslant x \leqslant {{B}_{x}},$ $0 \leqslant y \leqslant {{B}_{y}},z > 0\} $ в положительном октанте целочисленной решетки, число ${{B}_{x}} + 1$ будем называть шириной Π, а ${{B}_{y}} + 1$ – длиной Π. Разностный полиномиальный оператор (1) будет иметь вид

(4)
$\begin{array}{*{20}{c}} {P({{\delta }_{1}},{{\delta }_{2}},{{\delta }_{3}}) = \mathop {\mathop {\sum }\limits_{j = 0} }\limits^m \mathop {\mathop {\sum }\limits_{{{i}_{y}} = 0} }\limits^{{{b}_{y}}} \mathop {\mathop {\sum }\limits_{{{i}_{x}}} }\limits^{{{b}_{x}}} {{c}_{{{{i}_{x}}{{i}_{y}}j}}}\delta _{1}^{{{{i}_{x}}}}\delta _{2}^{{{{i}_{y}}}}\delta _{3}^{j} = } \\ { = \mathop {\mathop {\sum }\limits_{j = 0} }\limits^m {{P}_{j}}\left( {{{\delta }_{1}},{{\delta }_{2}}} \right)\delta _{3}^{j},} \end{array}$
где ${{b}_{x}}$, ${{b}_{y}}$ и $m$ определяют размер разностной схемы, а

${{P}_{j}}({{\delta }_{1}},{{\delta }_{2}}) = \sum\limits_{{{i}_{x}} = 0}^{{{b}_{x}}} \sum\limits_{{{i}_{y}} = 0}^{{{b}_{y}}} {{c}_{{{{i}_{x}}{{i}_{y}}}}}\delta _{1}^{{{{i}_{x}}}}\delta _{2}^{{{{i}_{y}}}},\quad j = 0,1,...,m.$

Характеристический многочлен будет иметь вид

$P(s,w,{v}) = \sum\limits_{j = 0}^m \sum\limits_{{{i}_{x}} = 0}^{{{b}_{x}}} \sum\limits_{{{i}_{y}} = 0}^{{{b}_{y}}} {{c}_{{{{i}_{x}}{{i}_{y}}j}}}{{s}^{{{{i}_{x}}}}}{{w}^{{{{i}_{y}}}}}{{{v}}^{j}},$
где m – порядок разностного оператора $P({{\delta }_{1}},{{\delta }_{2}}$, δ3), ${{b}_{x}} < {{B}_{x}}$, ${{b}_{y}} < {{B}_{y}}$.

Зафиксируем $\beta = \left( {{{x}_{\beta }},{{y}_{\beta }},m} \right)$ такое, что ${{c}_{{{{x}_{\beta }}{{y}_{\beta }}m}}} \ne 0$ и рассмотрим множество ${{\Pi }_{\beta }}\, = \,\{ (x,y,z) \in \mathbb{Z}_{ + }^{3}:$ $0 \leqslant x - {{x}_{\beta }} \leqslant {{B}_{x}}$bx, $0 \leqslant y - {{y}_{\beta }} \leqslant {{B}_{y}} - {{b}_{y}},$ z > m – 1}, тогда $L = \Pi {{\backslash }}{{\Pi }_{\beta }}$ и задача Коши будет иметь вид:

найти решение разностного уравнения

(5)
$P({{\delta }_{1}},{{\delta }_{2}},{{\delta }_{3}})f(x,y,z) = g(x,y,z),\quad (x,y,z) \in \Pi $
удовлетворяющее условию
(6)
$f(x,y,z) = \varphi (x,y,z),\quad (x,y,z) \in L,$
где $g(x,y,z)$ и $\varphi (x,y,z)$ – заданные функции целочисленных аргументов.

Задачу (5), (6) назовем задачей Коши для полиномиального разностного оператора (4).

В [28] доказано, что задача (5), (6) однозначно разрешима, если выполняется условие

(7)
$\left| {{{c}_{{{{x}_{\beta }}{{y}_{\beta }}m}}}} \right| > \sum\limits_{\left( {{{\alpha }_{1}},{{\alpha }_{2}}} \right) \ne \left( {{{x}_{\beta }},{{y}_{\beta }}} \right),{{\alpha }_{3}} = m} \left| {{{c}_{{{{\alpha }_{1}}{{\alpha }_{2}}{{\alpha }_{3}}}}}} \right|.$

Поставим задачу вычислить значение функции $f\left( {x,y,z} \right)$ в точке $A$ с координатами $\left( {{{x}_{1}},{{y}_{1}},{{z}_{1}}} \right)$.

3. ОПИСАНИЕ ВХОДНЫХ ДАННЫХ

Решение задачи Коши (5), (6) для трехмерного разностного уравнения с постоянными коэффициентами в точке A с координатами $({{x}_{1}},{{y}_{1}},{{z}_{1}})$ представляет собой значение функции $f(x,y,z)$ в точке A. Алгоритм вычисления значения функции $f(x,y,z)$ в заданной точке c координатами $({{x}_{1}},{{y}_{1}},{{z}_{1}})$ имеет рекурсивный характер и сводится к вычислению значений функции $f(x,y,z)$ на конечном подмножестве точек $\left( {x,y,z} \right)$ из множества $L = \Pi {{\backslash }}{{\Pi }_{\beta }}$.

Начальные данные (6) задаются трехмерной матрицей $F$, содержащей конечное подмножество значений начальных данных задачи Коши. Коэффициенты трехмерного разностного уравнения задаются трехмерной матрицей C. Для технической реализации алгоритма матрицы коэффициентов C и начальных данных $F$ задаются по слоям, начиная с самого нижнего. Проиллюстрируем процедуру задания матриц $F$ и $C$.

Для разностного уравнения

$\begin{gathered} {{c}_{{000}}}f\left( {x,y,z} \right) + {{c}_{{100}}}f\left( {x + 1,y,z} \right) + \\ \, + {{c}_{{200}}}f\left( {x + 2,y,z} \right) + {{c}_{{010}}}f\left( {x,y + 1,z} \right) + \\ \end{gathered} $
(8)
$\begin{gathered} \, + {{c}_{{110}}}f\left( {x + 1,y + 1,z} \right) + {{c}_{{210}}}f\left( {x + 2,y + 1,z} \right) + \\ \, + {{c}_{{001}}}f\left( {x,y,z + 1} \right) + {{c}_{{101}}}f\left( {x + 1,y,z + 1} \right) + \\ \, + {{c}_{{201}}}f\left( {x + 2,y,z + 1} \right) + {{c}_{{011}}}f\left( {x,y + 1,z + 1} \right) + \\ \end{gathered} $
$\begin{gathered} \, + {{c}_{{111}}}f\left( {x + 1,y + 1,z + 1} \right) + \\ \, + {{c}_{{211}}}f\left( {x + 2,y + 1,z + 1} \right) = 0 \\ \end{gathered} $
трехмерная матрица коэффициентов C на первом слое будет содержать матрицу ${{C}_{0}}$ вида
${{C}_{0}} = \left( {\begin{array}{*{20}{c}} {{{c}_{{000}}}}&{{{c}_{{100}}}}&{{{c}_{{200}}}} \\ {{{c}_{{010}}}}&{{{c}_{{110}}}}&{{{c}_{{210}}}} \end{array}} \right),$
а на втором – матрицу ${{C}_{1}}$ вида

${{C}_{1}} = \left( {\begin{array}{*{20}{c}} {{{c}_{{001}}}}&{{{c}_{{101}}}}&{{{c}_{{201}}}} \\ {{{c}_{{011}}}}&{{{c}_{{111}}}}&{{{c}_{{211}}}} \end{array}} \right).$

Запишем трехмерную матрицу коэффициентов $C$ следующим образом

$C = \left( {\begin{array}{*{20}{c}} {{{C}_{0}}} \end{array}\left| {\begin{array}{*{20}{c}} {{{C}_{1}}} \end{array}} \right.} \right) = \left( {\begin{array}{*{20}{c}} {{{c}_{{000}}}}&{{{c}_{{100}}}}&{{{c}_{{200}}}} \\ {{{c}_{{010}}}}&{{{c}_{{110}}}}&{{{c}_{{210}}}} \end{array}\left| {\begin{array}{*{20}{c}} {{{c}_{{001}}}}&{{{c}_{{101}}}}&{{{c}_{{201}}}} \\ {{{c}_{{011}}}}&{{{c}_{{111}}}}&{{{c}_{{211}}}} \end{array}} \right.} \right).$

Для разностного уравнения (8) ${{b}_{x}} = 2$, ${{b}_{y}} = 1$, $m = 1$. Пусть $\beta = \left( {1,1,1} \right)$, ${{B}_{x}} = 4$, ${{B}_{y}} = 2$, тогда трехмерная матрица начальных данных $F$ будет иметь вид $F = \left( {\begin{array}{*{20}{c}} {{{F}_{0}}} \end{array}\left| {\begin{array}{*{20}{c}} {{{F}_{1}}} \end{array}} \right.} \right)$, где

${{F}_{0}} = \left( {\begin{array}{*{20}{c}} {\varphi \left( {0,0,0} \right)}&{\varphi \left( {1,0,0} \right)}&{\varphi \left( {2,0,0} \right)}&{\varphi \left( {3,0,0} \right)}&{\varphi \left( {4,0,0} \right)} \\ {\varphi \left( {0,1,0} \right)}&{\varphi \left( {1,1,0} \right)}&{\varphi \left( {2,1,0} \right)}&{\varphi \left( {3,1,0} \right)}&{\varphi \left( {4,1,0} \right)} \\ {\varphi \left( {0,2,0} \right)}&{\varphi \left( {1,2,0} \right)}&{\varphi \left( {2,2,0} \right)}&{\varphi \left( {3,2,0} \right)}&{\varphi \left( {4,2,0} \right)} \end{array}} \right),$
${{F}_{1}} = \left( {\begin{array}{*{20}{c}} {\varphi \left( {0,0,1} \right)}&{\varphi \left( {1,0,1} \right)}&{\varphi \left( {2,0,1} \right)}&{\varphi \left( {3,0,1} \right)}&{\varphi \left( {4,0,1} \right)} \\ {\varphi \left( {0,1,1} \right)}&*&*&*&{\varphi \left( {4,1,1} \right)} \\ {\varphi \left( {0,2,1} \right)}&*&*&*&{\varphi \left( {4,2,1} \right)} \end{array}} \right).$

Здесь элементы, обозначенные *, вычисляются при выполнении алгоритма. При этом, вычислить элемент $\varphi (1,2,1)$ не представляется возможным без вычисления элементов $\varphi \left( {2,2,1} \right)$, $\varphi (1,1,1)$, $\varphi (2,1,1)$, $\varphi \left( {3,1,1} \right)$, $\varphi \left( {3,2,1} \right)$. Таким образом, для нахождения неизвестных элементов необходимо решить систему линейных разностных уравнений вида (8) с использованием начальных данных $\varphi (i,j,k)$, где i = 0, ..., 4, $j = 0,1,2$, $k = 0$, $\varphi (i,j,k)$, i = 0, ..., 4, $j = 0$, $k = 1$ и $\varphi \left( {0,1,1} \right)$, $\varphi \left( {0,2,1} \right)$, φ(4, 1, 1), $\varphi \left( {4,2,1} \right)$.

Итак, входные данные конечны и имеют вид:

1. Трехмерная матрица коэффициентов C = = $({{c}_{{{{\alpha }_{1}}{{\alpha }_{2}}{{\alpha }_{3}}}}})$, ${{\alpha }_{1}} = 0, \ldots ,{{b}_{x}}$, ${{\alpha }_{2}} = 0, \ldots ,{{b}_{y}}$, ${{\alpha }_{3}} = 0,...,m$ размера $({{b}_{x}} + 1) \times ({{b}_{y}} + 1) \times (m + 1)$ из коэффициентов ${{c}_{{{{\alpha }_{1}}{{\alpha }_{2}}{{\alpha }_{3}}}}}$ трехмерного разностного уравнения;

2. Точка ${{\beta }_{0}} = ({{x}_{\beta }},{{y}_{\beta }},m)$;

3. Точка A с координатами $({{x}_{1}},{{y}_{1}},{{z}_{1}})$, определяющая координаты искомого значения функции $f(x,y,z)$ и количество слоев в трехмерной матрице начальных данных F;

4. Трехмерная матрица начальных данных F = = $(\varphi (x,y,z))$, размера $({{B}_{x}} + 1) \times ({{B}_{y}} + 1) \times ({{z}_{1}} + 1)$, где $(x,y,z) \in L$, для всех остальных значений $(x,y,z)$ значения $\varphi (x,y,z) = 0$.

Поскольку координаты элементов трехмерной матрицы коэффициентов разностного оператора и трехмерной матрицы начальных данных в декартовой системе координат $(X,Y,Z)$ не совпадают с их координатами в матрице (строка × столбец × слой), то следует перейти из декартовой системы координат $(D({{d}_{1}},{{d}_{2}},{{d}_{3}}))$ в “матричную” $(M({{m}_{1}},{{m}_{2}},{{m}_{3}}))$ по правилу: D(d1, d2, ${{d}_{3}}) \to M({{m}_{1}}$, m2, m3), где ${{m}_{1}} = {{d}_{2}} + 1,$ ${{m}_{2}} = {{d}_{1}} + 1$, ${{m}_{3}} = {{d}_{3}} + 1$.

Далее необходимо будет проверить задачу Коши (5), (6) на разрешимость, т.е. проверить, выполняется ли условие (7) для коэффициентов разностного оператора (4).

4. ПРИМЕР

Будем рассматривать полиномиальный разностный оператор

$\begin{gathered} P({{\delta }_{1}},{{\delta }_{2}},{{\delta }_{3}}) = {{c}_{{000}}} + {{c}_{{100}}}{{\delta }_{1}} + {{c}_{{200}}}\delta _{1}^{2} + {{c}_{{010}}}{{\delta }_{2}} + \\ \, + {{c}_{{110}}}{{\delta }_{1}}{{\delta }_{2}} + {{c}_{{210}}}\delta _{1}^{2}{{\delta }_{2}} + {{c}_{{020}}}\delta _{2}^{2} + {{c}_{{120}}}{{\delta }_{1}}\delta _{2}^{2} + \\ \, + {{c}_{{220}}}\delta _{1}^{2}\delta _{2}^{2} + {{c}_{{001}}}{{\delta }_{3}} + {{c}_{{101}}}{{\delta }_{1}}{{\delta }_{3}} + {{c}_{{201}}}\delta _{1}^{2}{{\delta }_{3}} + \\ \end{gathered} $
(9)
$\begin{gathered} \, + {{c}_{{011}}}{{\delta }_{2}}{{\delta }_{3}} + {{c}_{{111}}}{{\delta }_{1}}{{\delta }_{2}}{{\delta }_{3}} + {{c}_{{211}}}\delta _{1}^{2}{{\delta }_{2}}{{\delta }_{3}} + {{c}_{{021}}}\delta _{2}^{2}{{\delta }_{3}} + \\ \, + {{c}_{{121}}}{{\delta }_{1}}\delta _{2}^{2}{{\delta }_{3}} + {{c}_{{221}}}\delta _{1}^{2}\delta _{2}^{2}{{\delta }_{3}} = 1 + 2{{\delta }_{1}} + 3\delta _{1}^{2} + 4{{\delta }_{2}} + \\ \end{gathered} $
$\begin{gathered} \, + 5{{\delta }_{1}}{{\delta }_{2}} + 6\delta _{1}^{2}{{\delta }_{2}} + 1\delta _{2}^{2} + 2{{\delta }_{1}}\delta _{2}^{2} + 3\delta _{1}^{2}\delta _{2}^{2} + 3{{\delta }_{3}} + \\ \, + 4{{\delta }_{1}}{{\delta }_{3}} + 5\delta _{1}^{2}{{\delta }_{3}} + 6{{\delta }_{2}}{{\delta }_{3}} + 80{{\delta }_{1}}{{\delta }_{2}}{{\delta }_{3}} + 8\delta _{1}^{2}{{\delta }_{2}}{{\delta }_{3}} + \\ \, + 3\delta _{2}^{2}{{\delta }_{3}} + 4{{\delta }_{1}}\delta _{2}^{2}{{\delta }_{3}} + 5\delta _{1}^{2}\delta _{2}^{2}{{\delta }_{3}}, \\ \end{gathered} $
где ${{b}_{x}} = 2$, ${{b}_{y}} = 2$, m = 1. Фиксируем β0 = = $({{x}_{\beta }},{{y}_{\beta }},m) = (1,1,1)$, ${{B}_{x}} = 4$, ${{B}_{y}} = 3$. Множество начальных данных будет иметь вид:
$L = \Pi {{\backslash }}{{\Pi }_{{\left( {1,1} \right)}}},$
где $\Pi = \{ (x,y,z) \in {{\mathbb{Z}}^{3}},$ $0 \leqslant x \leqslant 4,$ $0 \leqslant y \leqslant 3,z \geqslant 0\} $, ${{\Pi }_{{\left( {1,1} \right)}}} = \{ (x,y,z) \in \mathbb{Z}_{ + }^{3}:$ $1 \leqslant x \leqslant 3,$ $1 \leqslant y \leqslant 2,z \geqslant 1\} $. Трехмерная матрица коэффициентов C полиномиального разностного оператора (9) задается по слоям, начиная с самого нижнего:

$\begin{gathered} C = \left( {\begin{array}{*{20}{c}} {{{C}_{0}}} \end{array}\left| {\begin{array}{*{20}{c}} {{{C}_{1}}} \end{array}} \right.} \right) = \left( {\begin{array}{*{20}{c}} {{{c}_{{000}}}}&{{{c}_{{100}}}}&{{{c}_{{200}}}} \\ {{{c}_{{010}}}}&{{{c}_{{110}}}}&{{{c}_{{210}}}} \\ {{{c}_{{020}}}}&{{{c}_{{120}}}}&{{{c}_{{220}}}} \end{array}\left| {\begin{array}{*{20}{c}} {{{c}_{{001}}}}&{{{c}_{{101}}}}&{{{c}_{{201}}}} \\ {{{c}_{{011}}}}&{{{c}_{{111}}}}&{{{c}_{{211}}}} \\ {{{c}_{{021}}}}&{{{c}_{{121}}}}&{{{c}_{{221}}}} \end{array}} \right.} \right) = \\ \, = \left( {\begin{array}{*{20}{c}} 1&2&3 \\ 4&5&6 \\ 1&2&3 \end{array}\left| {\begin{array}{*{20}{c}} 3&4&5 \\ 6&{80}&8 \\ 3&4&5 \end{array}} \right.} \right). \\ \end{gathered} $

Расположение элементов матрицы коэффициентов C в декартовой системе координат представлено на рис. 1

Рис. 1.

Расположение элементов C в декартовой системе координат.

Поставим задачу: найти значение функции $f\left( {x,y,z} \right)$ в точке $A$ c координатами (2, 1, 4).

Трехмерная матрица данных Коши задается по слоям, начиная с самого нижнего:

$F = \left( {\begin{array}{*{20}{c}} {{{F}_{0}}} \end{array}\left| {\begin{array}{*{20}{c}} {{{F}_{1}}} \end{array}} \right.\left| {\begin{array}{*{20}{c}} {{{F}_{2}}} \end{array}} \right.\left| {\begin{array}{*{20}{c}} {{{F}_{3}}} \end{array}} \right.\left| {\begin{array}{*{20}{c}} {{{F}_{4}}} \end{array}} \right.} \right) = \left( {\begin{array}{*{20}{c}} 1&2&3&4&5 \\ 4&5&6&1&3 \\ 1&2&3&8&4 \\ 1&2&3&4&6 \end{array}\left| {\begin{array}{*{20}{c}} 3&4&5&7&6 \\ 6&0&0&0&8 \\ 3&0&0&0&4 \\ 6&5&3&2&3 \end{array}} \right.\left| {\begin{array}{*{20}{c}} 3&4&5&7&6 \\ 5&0&0&0&2 \\ 3&0&0&0&4 \\ 6&8&1&4&3 \end{array}} \right.} \right.\left. {\left| {\begin{array}{*{20}{c}} 3&4&5&7&6 \\ 6&0&0&0&1 \\ 3&0&0&0&4 \\ 6&6&7&1&4 \end{array}} \right.\left| {\begin{array}{*{20}{c}} 3&4&5&7&6 \\ 3&0&0&0&8 \\ 3&0&0&0&4 \\ 6&2&6&7&3 \end{array}} \right.} \right).$

Расположение элементов матрицы начальных данных $F$ в декартовой системе координат представлено на рис. 2.

Рис. 2.

Расположение элементов F в декартовой системе координат.

Произведем переход из декартовой системы координат в “матричную”: ${{\beta }_{0}} = [1,1,1] \to {{\beta }_{1}}$ = [2, 2, 2], $A(2,1,4) \to {{A}_{1}}(2,3,5)$.

Проверим задачу на разрешимость. Так как

$80 = {{c}_{{111}}} > \sum\limits_{\left( {{{\alpha }_{1}},{{\alpha }_{2}}} \right) \ne \left( {1,1} \right),{{\alpha }_{3}} = 1} \left| {{{c}_{{{{\alpha }_{1}}{{\alpha }_{2}}{{\alpha }_{3}}}}}} \right| = 38,$
то задача разрешима.

Для нахождения неизвестных (нулевых) значений в матрице $F$ на втором слое необходимо решить систему разностных уравнений

(10)
$\left( {\begin{array}{*{20}{l}} {P({{\delta }_{1}},{{\delta }_{2}},{{\delta }_{3}})f(x,y,z) = 0,\quad (x,y,z) \in \Pi ,} \\ {{\text{для}}\;{\text{которых}}\;\;z = 0,1,} \\ {f(x,y,z) = \varphi (x,y,z),\quad (x,y,z) \in L,} \\ {{\text{для}}\;{\text{которых}}\;\;z = 0,1.} \end{array}} \right.$

Матрица T системы разностных уравнений (10) будет блочной теплицевой матрицей (подробнее см. [29, 30]) вида

$T = \left( {\begin{array}{*{20}{c}} {{{T}_{0}}}&{{{T}_{1}}} \\ {{{T}_{{ - 1}}}}&{{{T}_{0}}} \end{array}} \right),$
где

${{T}_{{ - 1}}} = \left( {\begin{array}{*{20}{c}} 4&5&0 \\ 3&4&5 \\ 0&3&4 \end{array}} \right),\quad {{T}_{0}} = \left( {\begin{array}{*{20}{c}} {80}&8&0 \\ 6&{80}&8 \\ 0&6&{80} \end{array}} \right),$
${{T}_{1}} = \left( {\begin{array}{*{20}{c}} 4&5&0 \\ 3&4&5 \\ 0&3&4 \end{array}} \right).$

Решив систему уравнений (10), получим

${{F}_{1}} = \left( {\begin{array}{*{20}{c}} {3.00}&{4.00}&{5.00}&{7.00}&{6.00} \\ {6.00}&{ - 2.13}&{ - 1.86}&{ - 2.64}&{8.00} \\ {3.00}&{ - 1.40}&{ - 1.39}&{ - 2.48}&{4.00} \\ {6.00}&{5.00}&{3.00}&{2.00}&{3.00} \end{array}} \right).$

Для нахождения неизвестных (нулевых) значений в матрице $F$ на третьем слое необходимо решить систему разностных уравнений

(11)
$\left( {\begin{array}{*{20}{l}} {P({{\delta }_{1}},{{\delta }_{2}},{{\delta }_{3}})f(x,y,z) = 0,\quad (x,y,z) \in \Pi ,{\kern 1pt} } \\ {{\text{для}}\;{\text{которых}}\;\;z = 1,2,} \\ {f(x,y,z) = \varphi (x,y,z),\quad (x,y,z) \in L,{\kern 1pt} } \\ {{\text{для}}\;{\text{которых}}\;\;z = 1,2.} \end{array}} \right.$

Матрица системы уравнений (11) будет равна матрице системы уравнений (10). Решив систему уравнений (11), получим

${{F}_{2}} = \left( {\begin{array}{*{20}{c}} {3.00}&{4.00}&{5.00}&{7.00}&{6.00} \\ {5.00}&{ - 1.19}&{ - 0.50}&{ - 1.99}&{2.00} \\ {3.00}&{ - 1.01}&{0.02}&{ - 1.22}&{4.00} \\ {6.00}&{8.00}&{1.00}&{4.00}&{3.00} \end{array}} \right).$

Продолжая процесс, найдем

${{F}_{3}} = \left( {\begin{array}{*{20}{c}} {3.00}&{4.00}&{5.00}&{7.00}&{6.00} \\ {6.00}&{ - 1.40}&{ - 0.75}&{ - 1.52}&{1.00} \\ {3.00}&{ - 1.48}&{ - 0.36}&{ - 1.29}&{4.00} \\ {6.00}&{6.00}&{7.00}&{1.00}&{4.00} \end{array}} \right),$
${{F}_{4}} = \left( {\begin{array}{*{20}{c}} {3.00}&{4.00}&{5.00}&{7.00}&{6.00} \\ {3.00}&{ - 1.16}&{ - 0.71}&{ - 2.11}&{8.00} \\ {3.00}&{ - 1.17}&{ - 0.47}&{ - 1.88}&{4.00} \\ {6.00}&{2.00}&{6.00}&{7.00}&{3.00} \end{array}} \right).$

Таким образом, искомое значение $f(2,1,4)$ = = ‒0.71.

5. ОПИСАНИЕ АЛГОРИТМА

Алгоритм был реализован в среде MatLab2014 32bit. Вычисления производились на машине Intel(R) Core(TM) i5-3330S CPU 2.70 GHz, 32bit, ОЗУ 4.00 Гб под управлением Windows 7 Корпоративная SP1. Время счета для приведенного примера составило менее 1 секунды.

Алгоритм 1. Схема алгоритма.

Input: $C$, ${{\beta }_{0}}$, $A$, $F$.
Output: Значение функции $f(x,y,z)$ в точке A.
begin
    ${{\beta }_{1}}: = $ координаты точки ${{\beta }_{0}}$ в “матричной” системе координат;
    if$\left| {C\left( {{{\beta }_{1}}\left( 1 \right),{{\beta }_{1}}\left( 2 \right),{{\beta }_{1}}\left( 3 \right)} \right)} \right| < \sum \left( {\sum \left| {C\left( {:,:,{{\beta }_{1}}\left( 3 \right)} \right)} \right|} \right) - \left| {C\left( {{{\beta }_{1}}\left( 1 \right),{{\beta }_{1}}\left( 2 \right),{{\beta }_{1}}\left( 3 \right)} \right)} \right|$  
    then  
    return Ошибка ввода матрицы C;
    A1 := координаты точки A в “матричной” системе координат; $Q: = C(:,:,end)$;  
    $p = size\left( {Q,2} \right)$; $NAD = {{\beta }_{1}}\left( 1 \right) - 1$;  
    $POD = size\left( {Q,1} \right) - {{\beta }_{1}}\left( 1 \right)$;  
    $sTblock = size\left( {F,1} \right)*size\left( {F,2} \right) - size\left( {F,1} \right)*\left( {size\left( {Q,2} \right) - 1} \right) - $  
    $ - \left( {size\left( {F,2} \right) - \left( {size\left( {Q,2} \right) - 1} \right)} \right)*\left( {size\left( {Q,1} \right) - 1} \right)$;  
    $T: = zeros\left( {sTblock,sTblock} \right)$;  
    $a = Q\left( {{{\beta }_{1}}\left( 1 \right),:} \right)$; $e = {{\beta }_{1}}\left( 2 \right)$;  
    $colum{{n}_{c}} = zeros\left( {p,1} \right)$; $col = 1$;  
    while$e \geqslant 1$do  
    $colum{{n}_{c}}\left( {col} \right): = a\left( e \right)$; $col: = col + 1$; $e: = e - 1$;
    $e: = {{\beta }_{1}}(2)$; $ro{{w}_{c}}: = zeros\left( {1,p} \right)$; $r: = 1$;  
    for $m$ from $e$ to $length\left( a \right)$do  
     $ro{{w}_{c}}(r): = a(e)$; $r: = r + 1$; $e: = e + 1$;
    ${{T}_{0}}: = toeplitz(colum{{n}_{c}},ro{{w}_{c}})$;  
    fori from 0 to $sTblock{\text{/}}(size(Q,2)) - 1$do  
    $T(((i*p) + 1):p*(i + 1),((i*p) + 1):p*(i + 1)): = {{T}_{0}}$;
   forw from 1 to N AD do  
      $a: = Q({{\beta }_{1}}(1) - w,:)$; $e: = {{\beta }_{1}}(2)$;  
      $colum{{n}_{c}}: = zeros(p,1)$; $col: = 1$;  
      while $e \geqslant 1$do  
      $colum{{n}_{c}}(col): = a(e)$; $col: = col + 1$; e := e + 1;
     $e: = {{\beta }_{1}}(2)$; $ro{{w}_{c}} = zeros(1,p)$; $r: = 1$;  
     for $m$from$e$to$length\left( a \right)$do  
      $ro{{w}_{c}}(r): = a(e)$; $r: = r + 1$; $e: = e + 1$;
     ${{T}_{{NAD}}}: = toeplitz\left( {colom{{n}_{c}},ro{{w}_{c}}} \right)$;  
     for $i$ from 0 to $sTblock{\text{/}}(size(Q,2)) - 1 - w$do  
      $T((((i*p) + 1)):p*(i + 1),((i*p) + 1) + p*w:p*(i + 1) + p*w): = {{T}_{{NAD}}}$;
   forw from 1 to PODdo  
     $a: = Q({{\beta }_{1}}(1) + w,:)$; $e: = {{\beta }_{1}}(2)$;  
     $colum{{n}_{c}} = zeros(p,1)$; $col = 1$;  
     while $e \leqslant 1$do  
      $colum{{n}_{c}}(col) = a(e)$; $col = col + 1$; e = e – 1;
     $e = {{\beta }_{1}}(2)$; $ro{{w}_{c}} = zeros(1,p)$; $r = 1$;  
     for $m$ from $e$ to $length\left( a \right)$do  
      $ro{{w}_{c}}(r) = a(e)$; $r = r + 1$; e = e + 1;
     ${{T}_{{POD}}} = toeplitz(colum{{n}_{c}},ro{{w}_{c}})$;  
     for $i$ from 0 to $sTblock{\text{/}}(size(Q,2)) - 1 - w$do  
      $T(((i*p) + 1) + p*w:p*(i + 1) + p*w,((i*p) + 1):p*(i + 1)): = {{T}_{{POD}}}$
   $b: = zeros\left( {size\left( {T,1} \right),1} \right)$;  
   fork from β1(3) to A(3) do  
     $t: = 1$;  
     fori from 0 to $(size(F,2) - size(C,2))$do;  
       forj from 0 to $(size(F,1) - size(C,1))$do  
        do
            ${{f}_{{iter}}}: = F(j + 1:size(C,1) + j,$
            $(i + 1):size(C,2) + i,$
            $k - size(C,3) + 1:k)$;
            ${{f}_{{new}}}: = {{f}_{{iter}}}.*C$;
            $b\left( t \right) = - \left( {\sum \left( {\sum \left( {\sum \left( {{{f}_{{new}}}} \right)} \right)} \right)} \right)$;
            $t: = t + 1$
     $elementy: = linsolve\left( {T,b} \right)$; $h: = 1$;
     fori from 0 to $(size(F,2) - size(C,2))$do
        for$j$ from 0 to $(size(F,1) - size(C,1))$do
          $F\left( {{{\beta }_{1}}\left( 1 \right) + j,{{\beta }_{1}}\left( 2 \right) + i,k} \right): = elementy\left( h \right);$
          $h: = h + 1$;
  returnf(A)

6. ЗАКЛЮЧЕНИЕ

Классические методы решения дифференциальных уравнений (методы Рунге–Кутта четвертого порядка, Эйлера, Ньютона и т.д.) создают определенные трудности при их применении в многомерном случае. В каждой конкретной задаче необходимо индивидуально подбирать способ решения. Это делает разработку универсального подхода решения многомерных дифференциальных уравнений с использованием указанных методов трудоемкой. В проведенном исследовании мы расширили возможность использования стандартных симметричных разностных схем для аппроксимации двумерных дифференциальных уравнений с постоянными коэффициентами, где рассматривается задача “на плоскости”, в которой, как правило, в качестве независимых переменных x,y используют временную и одну пространственную переменные. Предложенная нами схема, во-первых, имеет произвольное число точек, а во-вторых, может применяться для аппроксимации дифференциальных уравнений в трехмерном случае. Для уравнения теплопроводности это означало бы, что при описании теплопереноса появляется возможность учитывать анизотропность этого процесса, через добавление дополнительной пространственной переменной.

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

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

  1. Годунов С.К., Рябенький В.С. Разностные схемы. Введение в теорию. Изд. 2, перераб. и доп. М.: Наука. 1977. 440 с.

  2. Самарский А.А. Введение в теорию разностных схем. М.: Наука. 1971. 553 с.

  3. Bousquet-Melou M., Petkovsek M. Linear Recurrences with Constant Coefficients: the Multivariate Case // Discrete Mathematics. 2000. V. 225. № 1–3. P. 51–75.

  4. Стенли Р. Перечислительная комбинаторика: Пер. с англ. М.: Мир. 1990. 440 с.

  5. Стенли Р. Перечислительная комбинаторика. Деревья, производящие функции и симметрические функции: Пер. с англ. М.: Мир. 2005. 767 с.

  6. Duffin R.J. Basic Properties of Discrete Analytic Functions // Duke Math. J. 1956. V. 23. № 2. P. 335–363.

  7. Duffin R.J. Potential Theory on Rhombic Lattice // J. Combinatorial Theory. 1968. V. 5. № 3. P. 258–272.

  8. Данилов О.A. Интерполяционная формула Лагранжа для дискретной аналитической функции // Вестник НГУ. Серия: Математика, механика, информатика. 2008. Т. 8. № 4. С. 33–39.

  9. Данилов О.A., Медных А.Д. Дискретные аналитические функции многих переменных и формула Тейлора // Вестник НГУ. Серия: Математика, механика, информатика. 2009. Т. 9. № 2. С. 38–46.

  10. Даджион Д., Мерсеро Р. Цифровая обработка многомерных сигналов. М.: Мир. 1988. 448 с.

  11. Изерман Р. Цифровые системы управления. М.: Мир. 1984. 541 с.

  12. Abramov S.A., Barkatou M.A., van Hoeij M., Petkovsek M. Subanalytic Solutions of Linear Difference Equationce and Multidimensional Hypergeometric Sequences // Journal of Symbolic Computation. 2011. V. 46. № 11. P. 1205–1228.

  13. Abramov S.A., Gheffar A., Khmelnov D.E. Rational Solutions of Linear Difference Equations: Universal denominators and denominator bounds // Programming and Computer Software. 2011. V. 37. № 2. P. 78–86.

  14. Abramov S.A., Petkovsek M., Ryabenko A.A. Hypergeometric Solutions of First-order Linear Difference Dystems with Rational-function Coefficients. Lecture Notes in Computer Science 9301. 2015. P. 1–14.

  15. Abramov S.A. Search of Rational Solutions to Differential and Difference Systems by Means of Formal Series // Programming and Computer Software. 2015. V. 41. № 2. P. 65–73.

  16. Kytmanov A.A., Lyapin A.P., Sadykov T.M. Evaluating the Rational Generating Function for the Solution of the Cauchy Problem for a Two-dimensional Difference Equation with Constant Coefficients // Programming and computer software. 2017. V. 43. № 2. P. 105–111.

  17. Leinartas E.K., Lyapin A.P. On the Rationality of Multidimentional Recusive Series // Journal of Siberian Federal University. Mathematics & Physics. 2009. V. 2. № 4. P. 449–455.

  18. Lyapin A.P., Chandragiri S. Generating Functions For Vector Partition Functions And A Basic Recurrence Relation // Journal of Difference Equations & Applications. 2019. V. 25. № 7. P. 1052–1061.

  19. Lyapin A.P. Riordan’s Arrays and Two-dimensional Difference Equations // Journal of Siberian Federal University. Mathematics & Physics. 2009. V. 2. № 2. P. 210–220.

  20. Apanovich M.S., Lyapin A.P., Shadrin K.V. Solving the Cauchy Problem for a Two-Dimensional Difference Equation at a Point Using Computer Algebra Methods // Programming and Computer Software. 2021. V. 47. № 1. P. 1–5.

  21. Апанович М.С., Ляпин А.П., Шадрин К.В. Алгоритм вычисления решения задачи Коши для двумерного разностного уравнения с начальными данными, заданными в “полосе”. Программирование. 2022. № 4. С. 50–56.

  22. Murray J.D. Mathematical Biology. II: Spatial Models and Biomedical Applications. Berlin [Germany] : Springer-Verlag. 2003. 811 p.

  23. Братусь А.С., Новожилов А.С., Платонов А.П. Динамические системы и модели биологии. М.: ФИЗМАТЛИТ. 2010. 400 с.

  24. Ляпин А. П., Кучта Т. Сечения производящего ряда решения разностного уравнения в симплициальном конусе // Известия Иркутского государственного университета. Серия Математика. 2022. Т. 42. C. 75–89.

  25. Федорюк М.В. Асимптотика: интегралы и ряды. М.: Наука. 1987. 544 с.

  26. Рогозина М.С. Разрешимость разностной задачи Коши для многослойных неявных разностных схем // Вестник СибГАУ. Математика, механика, информатика. 2014. Т. 55. № 3. С. 126–130.

  27. Apanovich M.S. Correctness of a Two-dimensional Cauchy Problem for a Polynomial Difference Operator with Constant Coefficients // Journal of Siberian Federal University. Mathematics & Physics. 2017. V. 10. № 2. P. 199–205.

  28. Апанович М.С. О разрешимости трехмерной разностной задачи Коши // ИТНОУ: Информационные технологии в науке, образовании и управлении. 2018. № 5(9). С. 9–13.

  29. Апанович М.С., Ляпин А.П., Шадрин К.В. Вычисление последовательности главных миноров теплицевой ленточной матрицы // Прикладная математика & Физика. 2020. Т. 52. № 1. С. 5–10.

  30. Иохвидов И.С. Ганкелевы и теплицевы матрицы. М.: Наука. 1974. 264 с.

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