Участник:Rpu6/Алгоритм построения приближённого решения нерегулярно вырождающегося эллиптического дифференциального уравнения

Алгоритм построения приближённого решения нерегулярно вырождающегося эллиптического дифференциального уравнения
Последовательный алгоритм
Последовательная сложность	$O(KN\mathrm{max}(K, N))$
Объём входных данных	$KN + 2N$
Объём выходных данных	$KN$
Параллельный алгоритм
Высота ярусно-параллельной формы	$?$
Ширина ярусно-параллельной формы	$?$

Автор описания: Д. П. Емельянов.

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

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

Алгоритм решает следующую задачу:

$\begin{cases} y^2u''_{yy}+u''_{xx}+c(y)u'_y-a(y)u=f(x,y),(x,y)\in(0,1)\times(0,b), \\ u(0,y)=u(1,y)=u(x,b)=0,|u(x,0)|\lt +\infty. \end{cases}$

Функции $a(y), c(y)$ полагаем аналитическими в круге $|y|\lt R$, где $R\gt b$. Кроме того, $c(y), a(y) \ge 0$ в $[0,b]$, $c(0)=0$.

1.2 Математическое описание алгоритма

1. Приблизим функцию $f(x,y)$ конечными суммами ряда Фурье по $x$:

$f(x,y)=\sum_{k=1}^{K}f_k(y)\sin \pi kx.$

$f_k(y)=2\int\limits_0^1f(x,y)\sin \pi kx\,dx.$

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

2. Приблизим функции $f_k(y),a(y),c(y)$ полиномами:

$a(y)=\sum_{n=0}^{N}a_ny^n, c(y)=\sum_{n=1}^{N}c_ny^n, f_k(y)=\sum_{n=0}^{N}f_{kn}y^n.$

Разложение можно проводить по аналогии с разложением в ряд Фурье, но с использованием полиномов Лежандра.

3. Введём в рассмотрение числа:

$\lambda_k=\pi^2k^2+a_0, k = \overline{1,K},$

$r_k=\frac{1-c_1 + \sqrt{(c_1-1)^2+4\lambda_k}}{2}, k = \overline{1,K}.$

4. Построим функции $\eta_k(y),\varphi_k(y)$, они будут использованы при построении приближения к решению:

$\eta_k(y)=\sum_{n=0}^{N}\eta_{kn}y^n, \stackrel{\circ}{\varphi}_k(y)=\sum_{n=0}^{N}\varphi_{kn}y^n,$

$\begin{matrix} \eta_{kn}=\frac{f_{kn}-\sum_{l=2}^{n+1}c_l(n+1-l)\eta_{k,n+1-l}+ \sum_{l=1}^{n}a_l\eta_{k,n-l}}{n(n-1)+nc_1-\pi^2k^2-a_0},n = \overline{1,N},\\ \varphi_{kn}=\frac{ \sum_{l=1}^na_l\varphi_{k,n-l}- \sum_{l=2}^{n+1}(n-l+1+r_k)c_l\varphi_{k,n-l+1} }{n(n-1)+2nr_k+nc_1} ,n= \overline{1,N},\\ \eta_{k,0}=\frac{f_{k,0}}{-\pi^2k^2-a_0}, \varphi_{k,0}=1, \varphi_k(y)=-\frac{\eta_k(b)}{\stackrel{\circ}{\varphi}_k(b)}\stackrel{\circ}{\varphi}_k(y). \end{matrix}$

5. Построим приближение к решению:

$u(x,y)=\sum_{k=1}^{K}\left(\left(\frac{y}{b}\right)^{r_k}\varphi_k(y)+\eta_k(y)\right)\sin \pi kx.$

За счёт выбора $K$ и $N$ можно добиться той или иной точности.

1.3 Вычислительное ядро алгоритма

Алгоритм решения задачи можно разделить на три крупные части примерно одинаковой вычислительной сложности:
1. Разложение заданной функции $f(x,y)$ в ряд Фурье, разложение её коэффициентов Фурье в ряды Тейлора, разложение заданных функций $a(y), c(y)$ в ряды Тейлора. На этом этапе известные сеточные функции скалярно умножаются на тригонометрические функции или полиномы Лежандра. В первом случае коэффициент является скалярным произведением, во втором - линейной комбинацией соответствующих коэффициентов полиномов Лежандра, помноженных на отвечающие данному полиному скалярные произведения.
2. Построение коэффициентов ряда Тейлора для функций $\eta_k(y), \varphi_k(y)$ с последующим построением самих функций. Используются формулы этапа 4 пункта 1.2.
3. По построенным функциям $\eta_k(y), \varphi_k(y)$ строится решение задачи $u(x,y)$ (формула этапа 5 пункта 1.2).

1.4 Макроструктура алгоритма

Макроструктура этого алгоритма фактически совпадает с его вычислительным ядром (см. пункт 1.3).

1.5 Схема реализации последовательного алгоритма

Алгоритм принимает на вход сеточные функции $f(x_i, y_j), a(y_j), c(y_j); x_i = \frac{i}{K}, y_j = \frac{j}{N}; i = \overline{0,K-1}, j=\overline{0,N-1}$.

$f_k(y_j) = \frac{2}{K}\sum_{l=1}^K f(x_l, y_j) sin (\pi kx_l), k=\overline{1,K},\\ f_{kn} = \sum_{m=0}^{N-1} (\mathrm{coeff}_n L_m) \frac{2}{N} \sum_{l=0}^{N-1} f_k(y_l) L_m(y_l), n=\overline{0,N-1},\\ a_{n} = \sum_{m=0}^{N-1} (\mathrm{coeff}_n L_m) \frac{2}{N} \sum_{l=0}^{N-1} a(y_l) L_m(y_l),\\ c_{n} = \sum_{m=0}^{N-1} (\mathrm{coeff}_n L_m) \frac{2}{N} \sum_{l=0}^{N-1} c(y_l) L_m(y_l).$

Тут $\mathrm{coeff}_n$ - коэффициент при n-й степени y в полиноме-аргументе (подразумевается, что алгоритм их знает), $L_m$ - нормированный полином Лежандра степени m.
Формулы для вычислений на других шагах алгоритма изложены в пункте 1.2. Достаточно заметить, что при все функции вычисляются в узлах заданной сетки, суммирование можно вести, запоминая $y^{n-1}$ с предыдущего шага, вычисляется лишь $y^{n-1}y$.

1.6 Последовательная сложность алгоритма

Будем оценивать количество операций умножения вещественных чисел. Произведём оценку отдельно для каждого шага алгоритма.
1. $(K + 1)N$ умножение для вычисления $f_k(y_j)$ (N возникает в силу необходимости вычислять функцию в каждой точке сетки по y), $((N + 1)N + 1)K$ операций для вычисления $f_{kn}$ (считаем, что значения гармоник и полиномов Лежандра известны заранее), $2((N + 1)N + 1)$ операций для вычисления $a_n, c_n$. 2. Для вычисления n-го коэффициента требуется $5 + 3(n-1)$ операций, для вычисления всех коэффициентов - $\sum_{n=1}^{N-1} (5 + 3(n-1))=5(N-1) + 3\frac{N-2}{2}(N-1)=(N-1)(5 + 1,5 (N-2))$, суммирование $\stackrel{\circ}{\varphi}_k(b)$ требует $2N$ операций, коррекция $\varphi_k(y_j)$ требует $N + 2$ операции, построение самих функций - $N(2N)=2N^2$ операций. Это всё нужно проделать $K$ раз.
3. Вычисление значения функции в узле сетки стоит $3K$ операций умножения. Это нужно проделать в $KN$ узлах.
Выделим главный член асимптотики по всем шагам: $KN^2 + 3,5KN^2 + 3K^2N=O(KN\mathrm{max}(K,N)).$
Если положить $K=N$, то сложность алгоритма имеет порядок $O(N^3)$.

3 Литература

<references \>