Внешние барицентрические координаты для произвольных многоугольников и приближенный метод их вычисления

Полный текст

Введение

Теоретическое изучение колебательно-волновых процессов неизбежно связано с исследованием соответствующих краевых и начально краевых задач математической физики [1–4]. Одним из вычислительно эффективных методов их численного решения является барицентрический метод (БМ) [3]. С учетом выделенных в работах [4–9] алгоритмических особенностей реализаций вычислительная эффективность БМ основывается на формировании глобальной системы базисных функций для заданной области анализа  граница  которой параметризуется в кусочно-линейном представлении. Глобальные для  базисные функции составляются [6] с применением классических интерполяционных методов [5] в вводимой для  барицентрической системы координат [10–13]. Относительно простое аналитическое соотношение, позволяющее с заданной точностью составлять для  барицентрическую систему координат, получено в [13].

В целом текущие математические представления БМ [3–13] определяют его вычислительно эффективную применимость относительно численного решения внутренних краевых и начально краевых задач математической физики. Одно из направлений развития БМ состоит в формировании теоретических решений, унифицирующих его относительно исследования внешних краевых и начально краевых задач математической физики. Первичный этап в получении подобных решения состоит в введении понятия внешних барицентрических координат для  задаваемой произвольной многоугольной областью, и формировании простого аналитического соотношения, позволяющего с заданной точностью составлять для   барицентрическую систему координат.

Получение указанных теоретических решений составляет цель настоящей статьи. В основу их формирования положим результаты [13–15].

1. Постановка задачи

Пусть $\Omega \subset {ℝ}^2:0\in \Omega$ – односвязная область, ограниченная замкнутой ломаной линией без самопересечений при

$\partial \Omega =\underset{i=0}{\overset{N-1}{\cup }}{\Gamma }_i,$

где

${\Gamma }_i=\left\{x_i\left(t\right)=e_{i}t+P_i,t\in \left[\mathrm{0,1}\right]\right\};$

$e_i=P_{i+1\mathrm{mod}N}-P_i,$

$\left\{P_0,P_1\mathrm{,...,}{P}_{N-1}\right\}$ – множество неповторяющихся вершин $\Omega ,$ нумерация которых определена в порядке положительного обхода  [16].

Определение 1. Внешними барицентрическими координатами ${\overline{\zeta }}_i$ для $\Omega$ назовем набор $\overrightarrow{\overline{\zeta }}={\left({\overline{\zeta }}_i\right)}_N$ функций ${\overline{\zeta }}_i\left(x\right)\in \left[\mathrm{0,1}\right]$ $\left(x\in {ℝ}^2\backslash \overline{\Omega }\right),$ которые удовлетворяют условиям:

$\Delta {\overline{\zeta }}_i\left(x\right),$ $x\in {ℝ}^2\backslash \overline{\Omega };$ (1)

${\overline{\zeta }}_i\left(x\right)=t,$ $x\in {\Gamma }_{i-1};$

${\overline{\zeta }}_i\left(x\right)=1-t,$ $x\in {\Gamma }_i;$

${\overline{\zeta }}_i\left(x\right)=\mathrm{0,}$ $x\in \partial \Omega \backslash \left\{{\Gamma }_{i-1},{\Gamma }_i\right\}.$

Решение внешней задачи Дирихле (1) по аналогии с [13] выполним известным методом Фредгольма при представлении функции ${\overline{\zeta }}_i\left(x\right)$ в виде логарифмического потенциала двойного слоя [15, с. 93]:

${\overline{\zeta }}_i\left(x\right)=\underset{\partial \Omega }{\int }{\Phi }_i\left(y\right)\left[1-\frac{1}{2\pi }\frac{\partial \ln \left|x-y\right|}{\partial {\nu }_y}\right]d{l}_y,$ (2)

где $x\in {ℝ}^2\backslash \overline{\Omega };$ $\partial /\partial {\nu }_y$ – частная производная по внутренней нормали ${\nu }_y$ к $\partial \Omega$ в точке $y\in \partial \Omega ;$ $d{l}_y$ – дифференциал кривой $\partial \Omega ;$ ${\Phi }_i\left(y\right)$ – неизвестная плотность на границе $\partial \Omega$ области $\Omega ,$ однозначно определяемая из интегрального уравнения Фредгольма II рода [13, с. 93]:

${\Phi }_i\left(x\right)+2\underset{\partial \Omega }{\int }{\Phi }_i\left(y\right)\left[1-\frac{1}{2\pi }\frac{\partial \ln \left|x-y\right|}{\partial {\nu }_y}\right]d{l}_y=2U_i\left(x\right)$ (3)

$x\in \partial \Omega ,$ $i=\overline{\mathrm{0,}N-1}.$

В выражении (3) через $U_i\left(x\right)$ обозначены заданные в (1) значения ${\overline{\zeta }}_i\left(x\right)$ на $\partial \Omega .$

Следуя результатам [13; 14], для удобства представления решения задачи (2), (3) построим $\Omega \subset ℂ:0\in \Omega$ на комплексной плоскости $ℂ.$ При этом, с учетом параметризации $\partial \Omega$ и граничных условий из (1), интегральное уравнение (3) запишем в виде

${\phi }_j^i\left(t\right)+\sum _{k=0}^{N-1}\underset{0}{\overset{1}{\int }}{\phi }_k^i\left(s\right)K_{jk}\left(t,s\right)ds=u_j^i\left(t\right),$ (4)

где

$j=\overline{\mathrm{0,}N-1};$ ${\phi }_j^i\left(t\right)={\Phi }_i\left(x_j\left(t\right)\right);$

$u_j^i\left(t\right)=2U_i\left(x_j\left(t\right)\right)$;

$K_{jk}\left(t,s\right)$ – ядро интегрального уравнения (4), которое с учетом [13; 14] при $t,s\in \left[\mathrm{0,1}\right]$ определяется соотношением

$K_{jk}\left(t,s\right)=\left\{\begin{array}{l}2\left|e_k\right|-{\pi }^{-1}\times \\ \times \mathrm{Im}\left[e_k/\left(e_{k}s+P_k-e_{j}t-P_j\right)\right],j\ne k;\\ 2\left|e_k\right|,j=k.\end{array}\right.$ (5)

Решение интегрального уравнения (4) относительно ${\phi }_j^i\left(t\right)$ позволяет задать ${\overline{\zeta }}_i\left(x\right)$ при вычислении интеграла:

${\overline{\zeta }}_i\left(x\right)=\sum _{j=0}^{N-1}\underset{0}{\overset{1}{\int }}{\phi }_j^i\left(t\right)H_j\left(t,x\right)dt,$ (6)

где

$H_j\left(t,x\right)=\left|e_j\right|-\mathrm{Im}\left[e_j/\left(e_{j}t+P_j-x\right)\right]/2\pi .$ (7)

2. Приближенно аналитическое определение внешних барицентрических координат

Решение интегрального уравнения (4) по аналогии с [13; 14] предполагается выполнить с применением формулы Гейне [17, с. 169] при разложении ядра $K_{jk}\left(t,s\right)$ в виде

$K_{jk}\left(t,s\right)={\lambda }_0^{jk}\left(t\right)+\sum _{n=1}^{\infty }\left(2n+1\right){\lambda }_n^{jk}\left(t\right)L_n\left(2s-1\right);$ (8)

${\lambda }_0^{jk}\left(t\right)=2\left\{\begin{array}{l}\left|e_k\right|+{\pi }^{-1}\mathrm{Im}\left[Q_0\left(2\left(e_{j}t+R_{jk}\right)/e_k-1\right)\right],j\ne k,\\ \left|e_k\right|,j=k;\end{array}\right.$

${\lambda }_n^{jk}\left(t\right)=\left\{\begin{array}{l}2{\pi }^{-1}\mathrm{Im}\left[Q_n\left(2\left(e_{j}t+R_{jk}\right)/e_k-1\right)\right],j\ne k,\\ \mathrm{0,}j=k;\end{array}\right.$

$n>\mathrm{0,}$

где $R_{jk}=P_j-P_k;$ $L_n\left(\tau \right)$ и $Q_n\left(z\right)$ – многочлены Лежандра первого и второго рода соответственно, задаваемые с учетом известных [14; 17; 18] рекуррентных соотношений.

Разложение (8) ядра (5) интегрального уравнения (4) на две системы линейно независимых интегрируемых с квадратом функций $\left\{{\lambda }_n^{jk}\right(t\left)\right\},$ $\left\{L_n\right(2s-1\left)\right\}$ $(n=\mathrm{0,1,...})$ и результаты лемм [13] позволяют задать приближенно аналитическое решение внешней задачи Дирихле (1) при введении следующих представлений.

Неизвестную функцию плотности ${\phi }_j^i\left(t\right)$ в интегральном уравнении (4) формализуем выражением

${\phi }_j^i\left(t\right)=u_j^i\left(t\right)-\sum _{k=0}^{N-1}\sum _{n=0}^{\infty }{S}_{kn}^i\sqrt{2n+1}{\lambda }_n^{jk}\left(t\right),$ (9)

где

$S_{kn}^i=\sqrt{2n+1}\underset{0}{\overset{1}{\int }}{\phi }_k^i\left(s\right)L_n\left(2s-1\right)ds.$ (10)

Определение $S_{kn}^i$ выполним, применив метод PG-ядер [19], в котором при подстановке (9) в (4) задачу нахождения коэффициентов $S_{kn}^i$ сведем к решению системы линейных алгебраических уравнений (СЛАУ):

$S^i\left(E+T\right)=U^i,$ (11)

где E – единичная матрица размера

$\left[N\left(M+1\right)\right]\times \left[N\left(M+1\right)\right];$

T – блочная матрица, составленная из элементов

$T_{nm}^{jk}=\sqrt{2n+1}\sqrt{2m+1}\underset{0}{\overset{1}{\int }}{\lambda }_m^{jk}\left(t\right)L_n\left(2t-1\right)dt$

$(n,m=\mathrm{0,1,...});$

$S^i$ – блочный вектор искомых коэффициентов разложения $S_{kn}^i;$ $U^i$ – блочный вектор, сформированный элементами [13]:

$U_{kn}^i=\sqrt{2n+1}\underset{0}{\overset{1}{\int }}{u}_k^i\left(t\right)L_n\left(2t-1\right)dt=$ (12)

$\left\{\begin{array}{l}\mathrm{1,}\left(i=k\vee i=k-1\right)\wedge n=\mathrm{0,}\\ -1/\sqrt{3},i=k\wedge n=\mathrm{1,}\\ 1/\sqrt{3},i=k-1\wedge n=\mathrm{1,}\\ \mathrm{0,}\text{ иначе.}\end{array}\right.$

Вычисление

$T_{nm}^{jk}=\sqrt{2n+1}\sqrt{2m+1}\underset{0}{\overset{1}{\int }}{\lambda }_m^{jk}\left(t\right)L_n\left(2t-1\right)dt$

может быть выполнено аналитически с применением результатов лемм [14] или численно по квадратурному методу Гаусса – Лежандра [3].

Перепишем уравнение (4) в операторной форме:

${\phi }^i+\mathcal{K}{\phi }^i=u^i,$ (13)

где ${\phi }^i={\left({\phi }_j^i\right)}_N$; $u^i={\left(u_j^i\right)}_N$; $\mathcal{K}={\left({\mathcal{K}}_{jk}\right)}_{N\times N}$ – матричный оператор;

$\mathcal{K}{\phi }^i=\left(\sum _{k=0}^{N-1}{\mathcal{K}}_{0k}{\phi }_k^i;\mathrm{...};\sum _{k=0}^{N-1}{\mathcal{K}}_{N-1k}{\phi }_k^i\right);$

$\left({\mathcal{K}}_{jk}{\phi }_k^i\right)\left(t\right)\equiv \underset{0}{\overset{1}{\int }}{K}_{jk}\left(t,s\right){\phi }_k^i\left(s\right)ds$

– линейный ограниченный оператор на пространстве функций из $C\left(\left[\mathrm{0,1}\right]\right)$ [13; 14].

Следуя [13; 14], введем в рассмотрение линейный ограниченный оператор:

$\left({\mathcal{K}}_{jk}^M{\phi }_k^i\right)\left(t\right)\equiv$ (14)

$\sum _{n=0}^M\left(2n+1\right){\lambda }_n^{jk}\left(t\right)\underset{0}{\overset{1}{\int }}{\phi }_k^i\left(s\right)L_n\left(2s-1\right)ds,$

и определим уравнение:

${\phi }_M^i+{\mathcal{K}}^M{\phi }_M^i=u^i,$ (15)

где ${\phi }_M^i={\left({\phi }_{jM}^i\right)}_N$ обозначает приближение функции плотности ${\phi }_j^i$ выражением (9) при замене бесконечной суммы по индексу  конечной с ограничением числа слагаемых до M.

C применением формулы Гейне [17, с. 169] запишем разложение (7) в виде

$H_j\left(t,x\right)=\left|e_j\right|-\frac{1}{2\pi }\mathrm{Im}\left[e_j{\left(e_{j}t+P_j-x\right)}^{-1}\right]=$ (16)

$\sum _{n=0}^{\infty }\left(2n+1\right){\overline{\lambda }}_n^j\left(x\right)L_n\left(2t-1\right);$

${\overline{\lambda }}_n^j\left(x\right)=\left\{\begin{array}{l}\left|e_j\right|+{\pi }^{-1}\mathrm{Im}\left[Q_n\left(2\left(x-P_j\right)/e_j-1\right)\right],n=0;\\ {\pi }^{-1}\mathrm{Im}\left[Q_n\left(2\left(x-P_j\right)/e_j-1\right)\right],n>\mathrm{0,}\end{array}\right.$

для которого введем в рассмотрение операторы [14]:

$\left(\mathcal{H}{\phi }^i\right)\left(x\right)\equiv \sum _{j=0}^{N-1}\underset{0}{\overset{1}{\int }}{H}_j\left(t,x\right){\phi }_j^i\left(t\right)dt;$ (17)

$\begin{array}{l}\left({\mathcal{H}}^M{\phi }^i\right)\left(x\right)\equiv \\ \sum _{j=0}^{N-1}\left\{\frac{1}{\pi }\sum _{n=0}^{\infty }\left(2n+1\right){\overline{\lambda }}_n^j\left(x\right)\underset{0}{\overset{1}{\int }}{L}_n\left(2t-1\right){\phi }_j^i\left(t\right)dt\right\}.\end{array}$

Обозначим единичный оператор символом $\mathcal{I}.$

Уточним определение норм некоторой векторной функции $\phi ={\left({\phi }_j\right)}_N$ в пространствах $C\left(\left[\mathrm{0,1}\right]\right)$ и $L_2\left(\left[\mathrm{0,1}\right]\right)$ [19]:

${||\phi ||}_{С}=\underset{\begin{array}{l}t\in \left[\mathrm{0,1}\right]\\ j\in \left\{\overline{\mathrm{0,}N-1}\right\}\end{array}}{\max }\left|{\phi }_j\left(t\right)\right|;$ ${||\phi ||}_{L_2}=\sqrt{\sum _{j=0}^{N-1}\underset{0}{\overset{1}{\int }}{\left[{\phi }_j\left(t\right)\right]}^{2}dt}.$  (18)

Для введенных представлений сформулируем основной результат настоящей статьи следующим утверждением.

Теорема 1. $\exists \tilde{M}\in \mathbb{N}:\forall M\ge \tilde{M}$ решение

${\overline{\zeta }}_i^M\left(x\right)=\sum _{j=0}^{N-1}\sum _{n=0}^M\sqrt{2n+1}{\overline{\lambda }}_n^j\left(x\right)S_{kn}^i$ (19)

внешней задачи Дирихле (1) существует и единственно, при этом справедливы оценки:

${||{\overline{\zeta }}_i^M-{\overline{\zeta }}_i||}_C\le \mathrm{const}\underset{j\in \left\{\overline{\mathrm{0,}N-1}\right\}}{\max }\left\{2{\varpi }_j^{-1}+\left|e_j\right|\right\}\times$ (20)

$\times \left[{\left(M+\mathrm{0,5}\right)}^{-1}+{\left(M+\mathrm{1,5}\right)}^{-1}\right];$

${||{\overline{\zeta }}_i^M-{\overline{\zeta }}_i||}_{L_2}=\mathrm{const}M{2}^{-M}{\left(2M+1\right)}^{-1}{\left(2M+3\right)}^{-1/2},$ (21)

где const – положительна и не зависит от M.

В выражениях (19)–(21) приняты обозначения [13; 14]:

${\varpi }_j=\min \left\{{\Theta }_j,\pi {\left(\pi +\left|\pi -{\alpha }_j\right|\right)}^{-1}\right\}$;

${\Theta }_j=\left\{\begin{array}{l}\left|\sin {\alpha }_j\right|,{\alpha }_j\in \left(\mathrm{0,}\pi /2\right)\cup \left(3\pi /\mathrm{2,2}\pi \right),\\ \mathrm{1,}{\alpha }_j\in \left[\pi /2,\pi \right)\cup \left(\pi ,3\pi /2\right],\end{array}\right.$ [20];

${\alpha }_j$ – внутренний угол $\Omega$ при вершине $P_j$.

Доказательство. Подставив разложение (16) в (6), получим

${\overline{\zeta }}_i\left(x\right)=\sum _{n=0}^{\infty }\left(2n+1\right)\sum _{j=0}^{N-1}{\overline{\lambda }}_n^j\left(x\right)\underset{0}{\overset{1}{\int }}{\phi }_j^i\left(t\right)L_n\left(2t-1\right)dt.$ (22)

Принимая во внимание (10), перепишем (22) в виде

${\overline{\zeta }}_i\left(x\right)=\sum _{j=0}^{N-1}\sum _{n=0}^{\infty }\sqrt{2n+1}{\overline{\lambda }}_n^j\left(x\right)S_{kn}^i.$ (23)

Из [21, с. 76] известно, что для интегрального оператора $\mathcal{K}$ существует обратный оператор ${\left(\mathcal{I}+\mathcal{K}\right)}^{-1},$ ввиду справедливости разложения (8) его ядра. Указанное обеспечивает разрешимость уравнения (13) в виде ${\phi }^i={\left(\mathcal{I}+\mathcal{K}\right)}^{-1}{u}^i$ и с учетом справедливости (23) определяет существование и единственность решения (19) внешней задачи Дирихле (1) для $M\to \infty .$ Принимая во внимание результаты доказательства теоремы 2 в [14], определим следующую оценку в $C\left(\left[\mathrm{0,1}\right]\right)\text{:}$

${||{\overline{\zeta }}_i^M-{\overline{\zeta }}_i||}_C\le {||\rho \mathcal{H}||}_C{||{\xi }_M^i-{\xi }_i||}_C+$ (24)

$+{||{\xi }_M^i||}_C{||\rho \left({\mathcal{H}}^M-\mathcal{H}\right)||}_C,$

где ${\phi }^i={\xi }_i\rho ,$ ${\phi }_M^i={\xi }_M^i\rho ,$ $\rho \equiv \rho \left(t\right)=2\sqrt{t-t^2}$ – весовая функция [13].

Учитывая представления (8), (9), (19), результаты [13; 14] для заданной параметризации $\partial \Omega$ и ${\overline{\zeta }}_i\left(x\right)\le 1$ при $x\in {ℝ}^2\backslash \overline{\Omega },$ установим справедливость следующих неравенств:

${||\rho \mathcal{H}||}_C=$ (25)

$\begin{array}{l}=\underset{\begin{array}{l}t\in \left[\mathrm{0,1}\right]\\ j\in \left\{\overline{\mathrm{0,}N-1}\right\}\\ x\in {ℝ}^2\backslash \overline{\Omega }\end{array}}{\max }\left|2\rho \left(t\right)\left(\left|e_j\right|-\frac{1}{2\pi }\mathrm{Im}\left[\frac{e_j}{e_{j}t+P_j-x}\right]\right)\right|\le \\ \le \underset{j\in \left\{\overline{\mathrm{0,}N-1}\right\}}{\max }\left|e_j\right|+\frac{1}{\pi }\underset{x\in {ℝ}^2\backslash \overline{\Omega }}{\max }\sum _{j=0}^{N-1}\underset{0}{\overset{1}{\int }}\left|\frac{e_j\sqrt{t-t^2}}{e_{j}t+P_j-x}\right|dt\le \\ \le \underset{j\in \left\{\overline{\mathrm{0,}N-1}\right\}}{\max }\left|e_j\right|+\underset{j\in \left\{\overline{\mathrm{0,}N-1}\right\}}{\min }\left(\frac{{\alpha }_j}{2\pi }\right)+1;\end{array}$

${||{\xi }_M^i||}_C={||{{\left[\rho \left(\mathcal{I}+{\mathcal{K}}^M\right)\right]}^{-1}}^i||}_C\le$ (26)

$\le 2{||{\left[\rho \left(\mathcal{I}+{\mathcal{K}}^M\right)\right]}^{-1}||}_C=\mathrm{const}\underset{j\in \left\{\overline{\mathrm{0,}N-1}\right\}}{\max }\left\{2{\varpi }_j^{-1}+\left|e_j\right|\right\};$

${||\rho \left({\mathcal{H}}^M-\mathcal{H}\right)||}_C\le$ (27)

$\begin{array}{l}\le \underset{x\in {ℝ}^2\backslash \overline{\Omega }}{\max }\sum _{j=0}^{N-1}\underset{0}{\overset{1}{\int }}\left|\frac{\rho \left(t\right)}{2\pi }\left(\sum _{n=0}^{\infty }\left(2n+1\right)Q_n\left(2\frac{x-P_j}{e_j}-1\right)\right.\right.\times \\ \times \left.\left.L_n\left(2t-1\right)-\frac{e_j}{x-e_{j}t+P_j}\right)\right|\le \\ \le \frac{\mathrm{const}}{2\pi }\left[{\left(M+\mathrm{0,5}\right)}^{-1}+{\left(M+\mathrm{1,5}\right)}^{-1}\right].\end{array}$

${||{\xi }_M^i-{\xi }_i||}_C\le \mathrm{const}\underset{j\in \left\{\overline{\mathrm{0,}N-1}\right\}}{\max }\left\{2{\varpi }_j^{-1}+\left|e_j\right|\right\}\times$ (28)

$\times \frac{{||\rho \left(\mathcal{K}-{\mathcal{K}}^M\right)||}_C{||\rho {\left(\mathcal{I}+\mathcal{K}\right)}^{-1}||}_C^2}{1-{||\rho {\left(\mathcal{I}+\mathcal{K}\right)}^{-1}||}_C{||\rho \left(\mathcal{K}-{\mathcal{K}}^M\right)||}_C}.$

Принимая во внимание результаты [13; 14] для (24)–(28), окончательно определим справедливость оценки (20), что и требовалось доказать относительно полиномиальной сходимости в $C\left(\left[\mathrm{0,1}\right]\right)$ решения (19) внешней задачи Дирихле (1).

С учетом результатов [13; 14] и полученных соотношений (25)–(28) зададим следующие оценки в $L_2\left(\left[\mathrm{0,1}\right]\right)\text{:}$

${||{\overline{\zeta }}_i^M-{\overline{\zeta }}_i||}_{L_2}\le {||\rho \mathcal{H}||}_{L_2}{||{\xi }_M^i-{\xi }_i||}_{L_2}+$ (29)

$+{||{\xi }_M^i||}_{L_2}{||\rho \left({\mathcal{H}}^M-\mathcal{H}\right)||}_{L_2};$

${||\rho \mathcal{H}||}_{L_2}\le \mathrm{const};$ (30)

${||{\xi }_M^i||}_C\le \mathrm{const};$ (31)

${||\rho \left({\mathcal{H}}^M-\mathcal{H}\right)||}_C\le {\frac{\mathrm{const}M{2}^{-M}}{\left(2M+1\right)\sqrt{2M+3}}}^{-1};$ (32)

${||{\xi }_M^i-{\xi }_i||}_{L_2}\le$ (33)

$\begin{array}{l}\le \mathrm{const}\frac{{||\rho \left(\mathcal{K}-{\mathcal{K}}^M\right)||}_{L_2}{||\rho {\left(\mathcal{I}+\mathcal{K}\right)}^{-1}||}_{L_2}^2}{1-{||\rho {\left(\mathcal{I}+\mathcal{K}\right)}^{-1}||}_{L_2}{||\rho \left(\mathcal{K}-{\mathcal{K}}^M\right)||}_{L_2}}\le \\ \le {\frac{\mathrm{const}M{2}^{-M}}{\left(2M+1\right)\sqrt{2M+3}}}^{-1}.\end{array}$

Принимая во внимание результаты [13; 14] для (29)–(33), окончательно определим справедливость оценки (21), что и требовалось доказать относительно экспоненциальной сходимости в $L_2\left(\left[\mathrm{0,1}\right]\right)$ решения (19) внешней задачи Дирихле (1). Теорема доказана.

3. Алгоритмическая реализация вычисления внешних барицентрических координат и тестовые примеры

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

Псевдокоды программы отразим поэтапным решением, взаимосвязанным с материалом пп. 1 и 2.

1. Задание многоугольной областью:

– число вершин и индексы: $N:=\mathrm{32,}$ $i:=\mathrm{0..}N-\mathrm{1,}$ $j:=\mathrm{0..}N-1;$ 

– положение вершин с нумерацией в порядке положительного обхода Ω:

  $P_0:=\left(\begin{array}{c}10\\ 1\end{array}\right),$ $P_1:=\left(\begin{array}{c}10\\ 2.5\end{array}\right),$ $P_2:=\left(\begin{array}{c}8\\ 2\end{array}\right),$ $P_3:=\left(\begin{array}{c}8\\ 5\end{array}\right),$

  $P_4:=\left(\begin{array}{c}6.5\\ 8.5\end{array}\right), P_5:=\left(\begin{array}{c}10.5\\ 9\end{array}\right), P_6:=\left(\begin{array}{c}8\\ 6\end{array}\right), P_7:=\left(\begin{array}{c}9\\ 5\end{array}\right),$

$P_8:=\left(\begin{array}{c}11.5\\ 9.5\end{array}\right), P_9:=\left(\begin{array}{c}6.5\\ 10.5\end{array}\right), P_{10}:=\left(\begin{array}{c}6.5\\ 12\end{array}\right), P_{11}:=\left(\begin{array}{c}8\\ 13\end{array}\right),$

$P_{12}:=\left(\begin{array}{c}9.5\\ 15.5\end{array}\right), P_{13}:=\left(\begin{array}{c}9\\ 16\end{array}\right), P_{14}:=\left(\begin{array}{c}8.5\\ 17\end{array}\right), P_{15}:=\left(\begin{array}{c}8\\ 18\end{array}\right),$

$P_{16}:=\left(\begin{array}{c}6.5\\ 18.5\end{array}\right), P_{17}:=\left(\begin{array}{c}5\\ 18\end{array}\right), P_{18}:=\left(\begin{array}{c}4.5\\ 17\end{array}\right), P_{19}:=\left(\begin{array}{c}4\\ 16\end{array}\right),$

$P_{20}:=\left(\begin{array}{c}3.5\\ 15.5\end{array}\right), P_{21}:=\left(\begin{array}{c}4.5\\ 13\end{array}\right), P_{22}:=\left(\begin{array}{c}5.5\\ 12\end{array}\right),$

$P_{23}:=\left(\begin{array}{c}5.5\\ 10.5\end{array}\right), P_{24}:=\left(\begin{array}{c}0.5\\ 7.5\end{array}\right), P_{25}:=\left(\begin{array}{c}1.5\\ 6\end{array}\right),$

$P_{26}:=\left(\begin{array}{c}5.5\\ 9.5\end{array}\right),$ $P_{27}:=\left(\begin{array}{c}3.5\\ 4.5\end{array}\right),$ $P_{28}:=\left(\begin{array}{c}3.5\\ 1\end{array}\right),$

$P_{29}:=\left(\begin{array}{c}4.5\\ 1\end{array}\right), P_{30}:=\left(\begin{array}{c}6\\ 4\end{array}\right), P_{31}:=\left(\begin{array}{c}7\\ 1\end{array}\right);$

– для выполнения условия $\Omega \subset {ℝ}^2:0\in \Omega$ осуществляется центрирование точек:

$\tilde{P}:=\frac{1}{N}\sum _{i=0}^{N-1}{P}_i=\left(\begin{array}{c}6.515625\\ 9.546875\end{array}\right),$ $P_j:=P_j-\tilde{P};$

– построение Ω на $ℂ$ и определение элементов, применяемых для параметризации $\partial \Omega \text{:}$

$Z_j:={\left(P_j\right)}_0+i{\left(P_j\right)}_1,$ $e_j:=z_{\mathrm{mod}\left(j+\mathrm{1,}N\right)}-z_j,$ $R_{ij}:=z_i-z_j;$

2. Формирование вспомогательных вычислительных функций:

– рекуррентное вычисление многочленов Лежандра первого и второго рода:

$L\left(x,n\right):=$if $n=-1$

then return 0

else if $n=0$

then return 1

else if $n=1$

then return x

else

$p0\leftarrow 1$

$p1\leftarrow x$

for $i\leftarrow 1$ to $n-1$ do

$p\leftarrow \frac{\left(2i+1\right)\cdot x\cdot p1-i\cdot p0}{i+1}$

$p0\leftarrow p1$

$p1\leftarrow p$

return $p$

$Q\left(x,n\right):=$if $n=-1$

then return $1-0.5\ln \left(1-x^2\right)$

else if $n=0$

then return $\mathrm{arccoth}\left(x\right)$

else if $n=1$

then return $x\mathrm{arccoth}\left(x\right)-1$

else

$p0\leftarrow \mathrm{arccoth}\left(x\right)$

$p1\leftarrow x\mathrm{arccoth}\left(x\right)-1$

for $i\leftarrow 1$ to $n-1$ do

$p\leftarrow \frac{\left(2i+1\right)\cdot x\cdot p1-i\cdot p0}{i+1}$

$p0\leftarrow p1$

$p1\leftarrow p$

return p

– ядро (5) интегрального уравнения (4):

$K\left(t,s,j,k\right):=\left\{\begin{array}{l}2\left|e_k\right|-\frac{1}{\pi }\mathrm{Im}\left[\frac{e_k}{e_{k}s+Z_k-e_{j}t-Z_j}\right],j\ne k,\\ 2\left|e_k\right|,j=k;\end{array}\right.$

– разложение (8) ядра (5) для заданного M (для примера положим $M:=4)\text{:}$

$\lambda \left(t,j,k,n\right):=$if $j=k$

then if $n=0$

then return $2\left|e_k\right|$

else return 0

else

$r\leftarrow \frac{1}{\pi }\mathrm{Im}\left[Q\left(2\frac{e_{j}t+R_{jk}}{e_k}-\mathrm{1,}n\right)\right]$

if $n=0$

then $r\leftarrow r+\left|e_k\right|$

return $2r$

$\tilde{K}\left(t,s,j,k\right):=\sum _{n=0}^M\left[\left(2n+1\right)\lambda \left(t,i,j,n\right)L\left(2s-\mathrm{1,}n\right)\right],$

$\tilde{K}\left(\mathrm{0.2,0.8,1,3}\right)=\mathrm{7.630364353568128,}$

$K\left(\mathrm{0.2,0.8,1,3}\right)=7.630398154688569;$

– функция $H_j\left(t,x\right)$ в (7):

$H\left(j,t,x\right):=\left|e_j\right|-\frac{1}{2\pi }\mathrm{Im}\left[\frac{e_j}{e_{j}t+Z_j-x}\right];$

– разложение (16) функции $H_j\left(t,x\right)\text{:}$

$\overline{\lambda }\left(j,n,x\right):=$$r\leftarrow \frac{1}{\pi }\mathrm{Im}\left[Q\left(2\frac{x-Z_j}{e_j}-\mathrm{1,}n\right)\right]$

if $n=0$

then $r\leftarrow r+\left|e_j\right|$

return r

$\tilde{H}\left(j,t,x\right):=\sum _{n=0}^M\left[\left(2n+1\right)\overline{\lambda }\left(j,n,x\right)L\left(2t-\mathrm{1,}n\right)\right],$

$\tilde{H}\left(\mathrm{1,0.2,0.1}+0.2i\right)=\mathrm{2.1025188715423075,}$

$H\left(\mathrm{1,0.2,0.1}+0.2i\right)=2.102518880606145;$

– элементы $T_{nm}^{jk}$ блочной матрицы T в составляемой СЛАУ:

$T\left(j,k,n,m\right)=\sqrt{2n+1}\sqrt{2m+1}\underset{0}{\overset{1}{\int }}\lambda \left(t,j,k,m\right)L\left(2t-\mathrm{1,}n\right)dt;$

– вспомогательная функция определения индекса в диапазоне от $\left[0;N-1\right]$ и элементы $U_{kn}^i$ блочного вектора $U^i,$ вычисляемые по правилу (12):

$U\left(k,i,n\right):=$if $k=i$

then if $n=0$

then return 1

else if $n=1$ then return $-1/\sqrt{3}$

else return 0

else if $k=\mathrm{sp}\left(i-\mathrm{1,}N\right)$

then if $n=0$

then return 1

else if $n=1$ then return $1/\sqrt{3}$

else return 0

else return 0

$\mathrm{sp}\left(i,n\right)=\left\{\begin{array}{l}\mathrm{mod}\left(i,n\right),i>\mathrm{0,}\\ \mathrm{mod}\left(n+i,n\right),i\le 0.\end{array}\right.$

3. Составление СЛАУ (11) и ее решение при $S_{kn}^i\text{:}$

– число элементов в блочных векторах и матрицах составляемой СЛАУ:

$NM:=N\left(M+1\right)=160;$

– вычисление T и $U^i:$

$T:=$for $q\leftarrow 0$ to $NM-1$ do

$j\leftarrow \mathrm{ceil}\left(\frac{q+1}{M+1}\right)-1$

$n\leftarrow \mathrm{mod}\left(q,M+1\right)$

for $g\leftarrow 0$ to $NM-1$ do

$k\leftarrow \mathrm{ceil}\left(\frac{g+1}{M+1}\right)-1$

$m\leftarrow \mathrm{mod}\left(g,M+1\right)$

$T_{qg}=T\left(j,k,n,m\right)$

return T

$U^i:=$for $q\leftarrow 0$ to $NM-1$ do

$j\leftarrow \mathrm{ceil}\left(\frac{q+1}{M+1}\right)-1$

$n\leftarrow \mathrm{mod}\left(q,M+1\right)$

$f_q=U\left(k,i,n\right)$

return f

– формирование единичной матрицы E и решение СЛАУ (11) при определении $S_{kn}^i\text{:}$

$E:=diag\left(1\right);$ ${\tilde{S}}^i:={\left(E+T\right)}^{-1}{U}^i;$

$S^i:=$for $q\leftarrow 0$ to $NM-1$ do

$j\leftarrow \mathrm{ceil}\left(\frac{q+1}{M+1}\right)-1$

$n\leftarrow \mathrm{mod}\left(q,M+1\right)$

$S_{jn}={\tilde{S}}_q^i$

return S

4. Вычисленная функция плотности ${\phi }_j^i\left(t\right)$ по правилу (9):

$u\left(i,j,t\right):=\left\{\begin{array}{l}2\left(1-t\right),i=j,\\ 2t,\text{ sp}\left(i-\mathrm{1,}N\right)=j,\\ \mathrm{0,}\text{ otherwise;}\end{array}\right.$

$\phi \left(i,j,t\right):=u\left(i,j,t\right)-\sum _{k=0}^{N-1}\sum _{n=0}^M\left[\sqrt{2n+1}\lambda \left(t,j,k,n\right)S_{kn}^i\right];$

5. Определение внешних барицентрических координат по правилу (19):

$\overline{\zeta }\left(i,x\right):=\sum _{j=0}^{N-1}\sum _{n=0}^M\left[\sqrt{2n+1}\overline{\lambda }\left(j,n,x\right)S_{jn}^i\right].$

На рис. 1 приведены результаты расчета линий уровня внешней барицентрической координаты (БК) ${\overline{\zeta }}_1\left(x\right)$ для рассматриваемого в приведенном листинге $\Omega$ при различных M.

Рис. 1. Результаты расчета ${\overline{\mathbf{\zeta }}}_{\mathbf{1}}\mathbf{\left(}\mathbf{x}\mathbf{\right)}$ для $\mathit{\Omega }$

Fig. 1. Calculation results of ${\overline{\mathbf{\zeta }}}_{\mathbf{1}}\mathbf{\left(}\mathbf{x}\mathbf{\right)}$ for $\mathit{\Omega }$

На рис. 2 показаны результаты расчета линий уровня различных внешних БК для рассматриваемого в приведенном листинге $\Omega .$

Рис. 2. Результаты расчета ${\overline{\mathbf{\zeta }}}_{\mathbf{i}}\left(x\right)$ для $\mathit{\Omega }$ при $\mathit{M}\mathbf{=}\mathbf{7}$

Fig. 2. Calculation results of ${\overline{\mathbf{\zeta }}}_{\mathbf{i}}\left(x\right)$ for $\mathit{\Omega }$ at $\mathit{M}\mathbf{=}\mathbf{7}$

Заключение

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

Главный результат настоящей статьи состоит в формировании теоретического решения, составляющего исходную основу для задания глобальных базисных функций в  при численном решении внешних краевых и начально краевых задач математической физики в приближении барицентрического метода. Автор статьи считает, что приведенные подробные результаты алгоритмической реализации вычисления внешних барицентрических координат вызовут интерес и сделают материал публикации доступнее широкому кругу читателей, позволят просто интерпретировать текущий результат в практическую реализацию решений [13; 14] и приведут к развитию барицентрического метода в направлениях, указанных в настоящей публикации, а также в работах [3; 22].

Об авторах

Иван Сергеевич Полянский

Автор, ответственный за переписку.
Email: van341@mail.ru
ORCID iD: 0000-0002-1282-1522
SPIN-код: 3238-3737
ResearcherId: R-1427-2017

доктор физико-математических наук, доцент, сотрудник Академии Федеральной службы охраны Российской Федерации

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

Дополнительные файлы

Доп. файлы

Действие

1. JATS XML

Скачать

2. Рис. 1. Результаты расчета для

Скачать (582KB)

Метаданные

3. Рис. 2. Результаты расчета для при

Скачать (465KB)

Метаданные