ВАРИАЦИОННЫЕ МЕТОДЫ ПОСТРОЕНИЯ РАЗНОСТНЫХ УРАВНЕНИИ И ИТЕРАЦИОННЫХ АЛГОРИТМОВ ДЛЯ ПРИБЛИЖЕННОГО РЕШЕНИЯ ЗАДАЧИ ДИРИХЛЕ Текст научной статьи по специальности «Математика»

АПВПМ-2019

ВАРИАЦИОННЫЕ МЕТОДЫ ПОСТРОЕНИЯ РАЗНОСТНЫХ

УРАВНЕНИИ И ИТЕРАЦИОННЫХ АЛГОРИТМОВ ДЛЯ ПРИБЛИЖЕННОГО РЕШЕНИЯ ЗАДАЧИ ДИРИХЛЕ

Н, М, Темирбеков1, 2, Ж, Р, Жакеылыкова3

1 Национальная Инженерная академия РК, Алматы, Казахстан 2Казахстанкий инженерно-технологический университет, Алматы, Казахстан 3Казахский Национальный педагогический университет имени Абая, Алматы, Казахстан

УДК 519.6

Б01: 10.24411/9999-016А-2019-10079

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

Ключевые слова: эллиптическое уравнение второго порядка, задача Дирихле, разностная схема, метод Ритца, метод сопряженных градиентов.

Введение

В данной работе рассматривается эллиптическое дифференциальное уравнение. При решении:

1) для построения разностных уравнений был применен вариационный метод, предложенный в 1908 году немецким математиком В. Ритцем, который называется методом Ритца. Найденное этим методом решение ип(х) при некоторых условиях, стремится к точи ому решению и (х), тогда п ^ то. Вопросы сходимости решений, получаемых методом Ритца, рассматриваются в многочисленных работах и монографиях.

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

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

1 Вариационные методы построения разностных уравнений

Рассмотрим в общем виде задачу в операторной форме Lu = f,u G Ф (L) где Ф (L) — область определения оператора L. Данная задача эквивалентна соответствующей вариационной задаче J (и) = min J (v) где

-УёФ(Ь)

J (v) = (Lu,v) - 2(f,u).

Рассмотрим применение метода Ритца в решении эллиптического дифференциального уравнения вида

2 а

£ i iA-и щ)=f м • (1)

i,j= 1

в^с краевым условием

U\dD =0 (2)

ISBN 978-5-901548-42-4

где В ограниченная область с кусочно-линейной границей дВ, А^^ (х) = А^ (х) ограниченные функции и для произвольного вектора £ = (£1 ,£2/ выполняется неравенство я-о ¿=1 С? < ^ 21=1 (х) && <

хеВ ^

вир ^2 ^.=1 Ац (х) < = с положительными постоянными ^о <

хео

Можно показать, что оператор задач (1) и (2) является симметричным и положительно определенным, а сама задача сводится к нахождению функции, минимизирующей в пространстве IV^ (В) квадратичный функционал

■7(и)= I I (ж) ^ щ

\г,3 = 1

¿В - 2 и]¿В ■¡в

(3)

Для нахождения приближенного решения задачи (3) применим метод Ритца с подпространствами ^ специального вида, удовлетворяющих условию (2)

Далее, построим подпространства РК. Для упрощения изложения проиллюстрируем это па примере кусочно-линейных аппроксимаций, когда область В = {(х1:х2) : 0 <х1,х2 < 1} является квадратом. Покроем эту область равномерной квадратной сеткой с шагом Н = (М — положительное целое число), а затем осуществим триангуляцию области В. Пронумеруем все внутренние вершины треугольников через {Рк,1}^1=1, объединение всех треугольников имеющих точку своей вершиной ри^ обозначим через ВВ1*1

можно представить в виде объединения шести треугольников < В* 1 т> , порядок нумерации которых

I ' ' ) т=1

указан на рис. 1.

Рис. 1: Обозначения треугольников

Каждому внутреннему узлу (хь,у{), к,1 = 1, N сетки поставим в соответствие кусочно-линейную базисную функцию (х, у)- Определяем каждую из этих функций (х, у) Для введенной сетки.

Чтобы задать (х, у) аналитически, достаточно для каждого треугольника, входящего в носитель (х, у), составить уравнение плоскости, проходящей через единицу в @к,1, а в остэшьных его двух вершинах через нуль. Вычисляя базисные функции <^¡-,1 (х, у) в каждом го треугольников \ В* 1 т> , построим

' I ' ' J т=1

систему базисных функций

¥к,1 (х,у) = <

- * (хк - х) - * (У1 - у),

- * (хк - х), (х,у) е

+ * (У1 - у), (х,у)е

В*

(х,У) е

в*

+ * (Хк - х) + * (У1 - у),

+ К (Хк - х) ■ К (У1 - У),

К

(х,у) е (х,у) е

Вк

Вк

ик,1,6■

(х,у) е Щ .А,

(4)

6

т

Для построения приближенного решения иК(х) задачи (1) и (2) воспользуемся методом Ритца с применением

базиса (х)^

к, 1=1

N

К

иН(х) = ак,1фк,1 (х) (5)

к, 1=1

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

Аа = д (6)

где а = (а1,а2, )Т — вектор, составленный из коэффициентов {аN(к-1)+1 = ак,1^ 1 = 1 разложения,

д = (д1, д2,gN2)Т — вектор с компонентами

9N(к-1)+1 = 9к,1 = !<Рк,1 (х)(Ю, к,1 =1,Ы (7)

А

2

ак,1 = ^(к-1)+1 ^^-1)+] = (х)к,1,г,0 = 1,Я (8)

А

(9)

А1,1 А1,2 0. .. 0 0 \

А = А2,1 А2,2 А2,3 . .. 0 0

0 0 0. • АN,N-1 АN,N )

где к,к =к к, к,к+1 =к+1 к' к,1 = 1, N и каждая из матриц к,1 является трехдиагональной матрицей порядка N.

Рассмотрим задачу Дирихле в области Б с границей дБ с переменными коэффициентами р(х, у) и д(х, у)

- ^р(х>у) ё- дъя(х>у) % =1 (х'у) • (10)

с краевым условием

ч1до = 0 (И)

Уравнение (10) умножаем на функцию и(х, у)ш интегрируем по Б учитывая граничное условие и1до = 0

J Р(х, дд^ + ! Я(х> = ! у)

I п J в

Согласно (3) построим функционал

2 л / ~ \ 2

3 (и) = ^ р(х, у)(^¿Б + I ц(х, у) (¿Б - 2 ^ /(х, у) • ьЯБ (12)

10 уиху JD \дУу JD

Для применения метода Ритца и(х, у) заменяем разложением =1 ак,1¥к,1 (х)

3 (иК (а.к,1)) = / р(х,у) ( ^ ак,1^к,1 (х) | ¿Б + [ д(х, у) ( ^ ак,1^к,1 (х) | ¿Б-

^ \дхк7=1 ) ^ \дУк7=1 )

N

-2 1(х, у) У2 ак,№к,1 (х)3,В, к,1 =1 ^ к,=1

Далее, производные д'7(^а("1°'1)) (к, I = 1,N) приравниваем нулю и ак,1 выносим за скобку

I( :^ ) ^™ +1*** ) ^

-! / (х,у)^к,1 (х)(1В = 0, к,1 = 1К. ■¡о

,к, 1=1

(13)

Таким образом, согласно (8) и (9) вводя обозначение

гк I = (к-1) + 1 (г

1)+3 = Р(х, У) (

\к,1 =

+ / Ч(х,У) ( У]

9<Рк,1 (х) ) д^ (х)

дх

к,1 = 1

(х) ) д(х) ду

дх

¿В+

(14)

¿В, = 1,М,

получаем систему линейных алгебраических уравнений (6).

Учитывая симметричность и блочную трехдиагональность матрицы А, нам достаточно определить

к,1 к,1 -1 к-1,1 к-1,1 +1 , , —Т7

ак,1, ак,1 , ак,1 , ак,1 = 1,Н.

1. Согласно (12), где г= к,1 находим акк,\ :

к,1

к,1 %к,1

Р(Х> У){+Ч(Х, У){^)

¿В

(15)

Вычисляем интегралы в каждой из области В* ь учитывая, вид базисной функции (х, у) в рассматриваемой области. Далее, все найденные шесть значения подставляем в (15). При этом, объединяем интегралы с одинаковыми значениями и используем обозначения

Рг

к± 1,1 = Н2

р (х, у) ¿хс!у

я

к,1 ±

1 Н2

ц (х, у) <1х<1у

Значит,

к,1 = гк,1 =

рк-2,1 + Я к,1 - 2 + рк+1,1 + Як,1+2 ■

(16)

(17)

(18)

2. Для того чтобы найти элемент а1 как и в первом случае воспользуемся формулой (14), где г,] = к,1 - 1 и для того чтобы определить место рь,1 -1 вершину рь,1 спускаем па один шаг вниз. По рис.1 видим что, вершины и р^ I-1 лежат в треугольниках В^ и В^ Значит, а^-1 определяем только в этих областях

1

1

а^-1 = 1п 1 ипв (Ф, У)'^ + Ч(х, У) ^ гШу =

= 1п 1 {р(х, У) ^ + ч(х, У)-^ )^У+ (!9)

+ 1пе (р(х, У) ^ ^ + 1(х,у)-^-

Найденные значения интегралов, подставляя в (19) и используя обозначение (17) получаем

аы-1 = -Т2 I Ч^^ у)ЛхЛу - I ц(х,у)д,хду = -■[ ц(х,у)д,хду = -(5к1 -2 (20)

Н ->В1 Н ./Аз Н JDl{JD 6 2

3. Чтобы найти акк11,1 по формуле а£ = ^ р(х, у) д^к-1}х) д<Р%}х) ¿Б + ^ Ч(х, у) ¿В, где 1,з = к - 1,1 передвигавм вершину рк}I на один шаг влево и по рис.1 видим, что рк}I и Рк-1,1 являются вершинами треугольников Б^ш Значит, а'к г 1,1 определяется в треугольниках и Б2-

ак—1,1 = fD1 ^ (Р(х, У) ^ Чх^ + Ж*, У) ^ ) <МУ = = -ц (р(х,у) ^ + Ч(х, У)-^ д-%^)<1х<1у+ (21)

+ fD2 (р(-> У)-^ ^ + ч(х,у) ^

Найденные значения а'к г 1,1 в треугольниках и Б2 подставляем в (21) с учетом (16) и получим

ак—1,1 = -То • [ Р(х,у)3,хд,у - ^г • [ р(х, у)д,хд,у = -^г • [ р(х, у)д,хд,у = -Рк— 1 , (22)

' Н2 Ь2 .¡D2 Ь-2 JDl UD2 2

4. Осталось определить элемент а,к 1 1 Для этого вершину рк,1 нужно будет с начало передвигать па один шаг влево, затем на один шаг поднять вверх. Рассуждая по рис.1 видим, что рк}1 и рк—1,1+1 являются вершинами треугольников и Бз. Но при выбранном сдвиге рк,1 вер шина Рк—1, г+1 не будут принадлежать ни одному из шести треугольников. Поэтому

а

к—1, +1 = к, =

0

А

получаем

(Рк— 1 ,1 + Як,1 —1 + Рк+1,1 + Як,1 + ^ ак,1 - Як,1 — 2 ак,1 — 1 - Як,1 + ^ ак,I +1 - Рк—1,1ак — 1,1 -

(24)

-Рк+\,1ак+1,1; М =

Пусть в (24)

(А1а)к,1 = -Рк—1,1ак—1,1 + (Рк—1,1 + Рк+2,1) ак,1 - Рк+1,1ак+1,1; к,1 = 1^ (25)

(А2а)к,1 = -Як,1—2ак,1—1 + (як,1—1 + Як,1+1) ак,1 - Як,1+1ак,1+1; М = 1^ (26)

тогда, используя А1 + А2 = получим систему (6).

2 Вариационные методы построения итерационных алгоритмов

Для численного решения системы (6) можно применять методы вариационного типа, такие как, метод скорейшего спуска, метод минимальных поправок, метод сопряженных градиентов и т.д.

Метод сопряженных градиентов является наиболее предпочтительным для систем с самосопряженной положительной матрицей А = А* > 0. При плохой обусловленности матрицы этот метод не всегда становится вычислительно устойчивым.

Оператор А является самосопряженным и положительно определенным оператором в пространстве сеточных функций Н0.В Н0 введем скалярное произведение (и, у) = ^норма ||м|| = у/(и, и).

Для решения уравнения (6) воспользуемся методом сопряженных градиентов. Итерационный процесс реализуется в следующем порядке:

0. Подготовка перед итерационным процессом

0 „• вычисляется невязка в0 „• ., , .. .

, , , ,

к

По заданному а0 ^ вычисляется невязка в0^ = г0^ = (Аа0). - дг = 1,...,п1, ] = 1,п2;

1. Вычисляется параметр: ^ = ((^^-т'-т), г^ = 1 - ^кАз

2. Следующее приближение решения вычисляется по формуле: ^ = -1 —

г,3

(гк - ,гк - )

3. Вычисляется параметр: 3к = А-1

4. Вспомогательное значение вычисляется по формуле: вк+1 = гк^ + 3п$к^ где

(А

5 ' ь 3 = ь1

Рг+1 ( — — Р.1 ^ (в ^ — 1,э)

+

1

ь2

(3i+1 ^ 1,3+1 — ^,з) — $ 1,3-2 (вьз — -1) ; ЬЗ = 1,Н

Этот процесс продолжается до тех пор, пока не выполнится критерий остановки итераций |ак+1 — ак II < е.

3 Примеры расчетов

Для иллюстрации предложенного метода рассмотрим пример задачи (10)—(11) в круге.

Пусть в прямоугольной области Р имеется к руг П радиуса И с центром (А1,А2). Требуется найти приближенное решение задачи Дирихле для уравнения Пуассона

д2и д2и

т: 2 + 2 х1 ' х2 дх( дх^

(27)

в круге, при условии и (х1,х2)|ап = 0.

Вспомогательная задача метода фиктивных областей

_д

дх

1 (р (х1,х2 ) + ддх1 (1 (хь х) = —х1 • х2,

ие (хъх2)1до = 0

, . , . ( 1, (х1,х2) е П,

р(х1,х2 ) = Ч(х1,х<21/£К (х1

ге / \ \ —х1 • х2, (х1, х2) е П (х1'хН 0, (х1 ,х2'_

, х ) еР/П. 0, (х1 ,х2) е Р/П.

Были проведены численные эксперименты. Ниже приводится таблица 1 для числа итераций N и погрешности вычисления \\иь — и||с при П1 = п2 = 101 и различных значениях е для метода Ритца.

Таблица 1: Результаты численных экспериментов

П1 п2 1/е N \\иь —и\\С

101 101 102 1221 0.003888

101 101 103 3395 0.003820

101 101 104 7920 0.003819

101 101 105 12737 0.003818

101 101 106 20042 0.003817

101 101 107 26592 0.003816

Полученные результаты показывают, что при увеличении числа узлов сетки погрешность решения уменьшается.

На рис. 2 приведен результат решения задачи по явной разностной схеме при е* = е = 10 -4 па равномерной сетке с размером 101x101. Значит, исходная система линейных алгебраических уравнений имеет 101x101 неизвестных.

Заключение

На сегодняшний день остается актуальной задача разработки и модификации численных методов. Однако процесс развития вычислительной техники смещает упор с создания новых численных методов, на исследование и классификацию старых, с целью выявлению лучших. Теперь для современных мощных компьютеров, такие характеристики как объём требуемой памяти, и число арифметических операции, не обязательно

Рис. 2: Решение задачи Дирихле для уравнения Пуассона в круге

стоят на первом плане. Более предпочтительными становятся те методы, которые отличаются удобством реализации на ЭВМ. и позволяют решать более широкий класс задач.

Особыми достоинствами данного метода являются его простота и низкие затраты памяти, что делает его эффективным при решении задач большой размерности.

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

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

[1] Самарский А.А. Теория разностных схем М.: Наука. 1983. 616 с.

[2] Марчук Г.И. Методы вычислительной математики М.: Наука. 1989. 608 с.

[3] Мухамбетжанов А. Т.. Отелбаев М. О.. Смагулов Ш. С. Об одном новом приближенном методе решение нелинейных краевых задач.// Вычислительные технологии. Том 3. №4. 1998.

[4] Н.М. Темирбеков, Ж.Р. Жаксылыкова. Об итерационном методе фиктивных областей для решения эллиптического уравнения в областях со сложной геометрией // Известия МКТУ им. Х.А.Ясави, Серия математика, физика, информатика. 2018. Т.1. №1, С.128 132.

fo] N.M.Temirbekov, Zh.R.Zhaksylykova. Arilterative Method for Solving Non linear Navier-Stokes Equationsin Complex Domains Taking into Account Boundary Conditions with Uniform Accuracy // AIP Conference Proceedings 2018. V.1997, 020036 (2018): doi: 10.1063/1.5049030.

[6] Temirbekov A.N., Wojcik W. Numerical Implementation of the Fictitious Domain Method for Elliptic Equatioris//Iritcrriatiorial Journal of Electronics and Telecommunications. 2014. Vol.60. Issue 3. P.219-223.

[7] Temirbekov A.N. Numerical implementation of the method of fictitious domains for elliptic equations // 3rd International Conference on Analysis and Applied Mathematics, ICAAM 2016. Vol. 1759, P. 020053-1 020053-6: doi: 10.1063/1.4959667.'

Темирбеков Нурлаи Муханович д.ф.-м.и., профессор вице-президент Национальной Инженерной академии РК;

e-mail: teinirbekovOrwinbler.ru;

Жаксылыкова Жадра Рахме.товна докторант 2-курса по специальности "6D060100 Математики"

Казахского национального педагогического университета имени Абая;

email: zhaksylykova0507Qmail. ru Дата поступления 30 апреля 2019 г.