Вычислительные технологии
Том 13, № 2, 2008
Прямой метод декомпозиции без наложения подобластей для решения краевых задач на прямоугольных квазиструктурированных сетках*
В. М. Свешников Институт вычислительной математики и математической геофизики СО РАН, Новосибирск, Россия e-mail: victor@lapasrv.sscc.ru
The method of solving of the boundary value problems on rectangular quasistructured grids is proposed. The grids are constructed in two stages: first, the computational domain is decomposed into subdomaines using a uniform rectangular macrogrid, then a uniform rectangular and locally-modified grid is constructed in each macrocell. An essence of the method is a construction and solving of the system of linear algebraic equations for values of the required function in nodes of the grid which covers an interface boundary between subdomains. Properties of the system have been investigated and the existence and uniqueness theorem has been proved. The proposed method is a direct one as no iterations on subdomains is implied.
Введение
Разностная сетка, на которой ищется решение краевой задачи, должна удовлетворять следующим, как правило, противоречивым требованиям: с одной стороны, содержать узлы, сгущающиеся в подобластях физической неоднородности, и учитывать конфигурацию границы, а с другой — иметь простую структуру. Первые требования призваны увеличить точность расчетов, второе — уменьшить трудозатраты на его получение.
В настоящей работе для решения краевых задач в двумерных областях строятся ква-зиструктурированные прямоугольные сетки, которые сочетают простоту прямоугольных сеток с адаптивными свойствами квазиструктурированных сеток [1]. Суть рассматриваемого подхода состоит в двухэтапном покрытии расчетной области О прямоугольными равномерными по каждому направлению сетками. На первом этапе в прямоугольнике, окружающем О, строится прямоугольная равномерная сетка Пя с шагами, по каждому направлению значительно превышающими шаги результирующей сетки П^, на которой проводится дискретизация исходной задачи, т. е. фактически проводится декомпозиция О на подобласти О к, принадлежащие равным прямоугольникам Як (к = 1 , К,К — число подобластей). На втором этапе в Як строятся прямоугольные равномерные подсетки П^,к. В прямоугольниках Як', содержащих границу расчетной
* Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований (грант № 08-01-00526).
© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2008.
области Г, осуществляется локальная модификация П^,к' путем сдвига приграничных узлов на Г. к
Решение краевой задачи проводится на сетке П = и П^,к по алгоритму, содержа-
к= 1
щему три этапа. На первом из них на сетке ш^, покрывающей границу 7 сопряжения подобластей и состоящей из узлов подсеток П^,к, принадлежащих 7, строится система линейных алгебраических уравнений относительно значений искомой функции в узлах ш^. Для этого решается серия вспомогательных задач в подобластях. Доказаны леммы, устанавливающие свойства матрицы данной системы, и теорема о том, что ее решение существует и оно единственно. Второй этап алгоритма — решение данной системы. Здесь необходимо подчеркнуть, что порядок системы уравнений намного меньше числа узлов результирующей сетки П^. На третьем этапе решаются краевые задачи в подобластях, причем в качестве граничных условий берутся значения функции в узлах ш^, рассчитанные на предыдущем этапе. Важным свойством третьего этапа является то, что расчеты проводятся только в подобластях, представляющих интерес для последующего анализа и обработки. Решение сеточных уравнений в подобластях осуществляется быстрыми решателями, включая оптимальный по быстродействию метод циклической редукции.
Предлагаемый алгоритм — это по сути дела прямой метод декомпозиции в том смысле, что в нем отсутствуют итерации по подобластям. Кроме того, он является в большей части параллельным и поэтому применим для решения краевых задач на многопроцессорных суперкомпьютерах.
Оценки эффективности данного алгоритма, численные эксперименты, подтверждающие его работоспособность, и перечень его основных свойств приведены в работе [2].
Настоящий подход имеет схожие черты с методом конечных суперэлементов Р.П. Федоренко [3], но в то же время имеет и ряд существенных отличий, из которых главным является то, что в предлагаемом подходе размер макроэлемента может без снижения точности расчетов значительно превышать шаг результирующей сетки.
1. Постановка задачи. Построение сеток
Пусть в замкнутой области О = О и Г с границей Г требуется решить краевую задачу
А и = д1, 1и|г = 02- (1)
Здесь и = и(Т) — искомая функция; д1 = д1(Т), д2 = д2(Т) — заданные функции (Т = (х,у) — текущая точка, где х,у — декартовы координаты); А — оператор Лапласа; I — оператор граничных условий. В дальнейшем мы будем предполагать, что условия задачи (1) обеспечивают существование и единственность достаточно гладкого решения и.
Опишем вокруг расчетной области прямоугольник Я = {0 < х < Пх, 0 < у < Оу}, где Ох, Оу — заданы (О С Я). Построим в Я прямоугольную равномерную макросетку
Пя = {хI = 1Их, YJ = зву, I = 0Ж, з = оТ^у, Нх = ^, Ну = ^| ,
где Ых,Ыу — заданные целые числа, с шагами Hx,Hy >> к (к — максимальный шаг сетки, на которой аппроксимируется задача (1)). Тем самым мы фактически проводим
декомпозицию расчетной области О на непустые подобласти О к = Оі,з = О П Кіз = 0,
к = 1,К, (1,3) £ М, где Яг,з — открытая прямоугольная подобласть, представляющая собой макроэлемент сетки Пя; К — известное целое число, которое определяется геометрией С; М — множество пар индексов (1,3), для которых выполняется указанное неравенство (цифрой 0 в данном случае обозначена пустая подобласть). Точки Тр = Тгз £ С, р = 1, Р, пересечения координатных линий х = Хг, у = Уз макросетки назовем макроузлами. Граница сопряжения подобластей, которую мы обозначим через 7, состоит из отрезков координатных линий макросетки Пя. Среди С к будем различать
Здесь множество Мг С М содержит (1, 3) такие, что ГП Ягз = 0, где Яг з — замыкание Яг,з. Пусть 7к С 7 — отрезок границы подобласти С к. Если рассматриваемая подобласть
В макроэлементах Кк = Кі з ((I, J) Є М) построим равномерные прямоугольные подсетки
причем будем предполагать без существенного ограничения общности, что пх^ = 2т'х’к, пУ)к = 2тУ’к, где тх,к, ту,к > 0 — целые числа. По аналогии с подобластями будем
целью их адаптации к границе подвергнем локальной модификации [4] путем сдвига приграничных узлов в точки пересечения координатных линий с границей.
Отметим, что при сделанных предположениях Пь,к являются квазисогласованными в том смысле, что на отрезках 7к границы сопряжения все узлы редких сеток совпадают с узлами густых.
Решение исходной задачи (1) будем искать на результирующей сетке
Кроме рассмотренных выше, введем сетку Шн С П н, Шн = {(хг,Уг) £ ч, (Хг,Уг) = Тр, г = 1, п}, где (хг, уг) — узлы Пн, а также подсетки Шн,к, состоящие из узлов Шн, лежащих на
1к (п — известно). Подчеркнем, что в Шн, Шн,к не входят макроузлы Тр.
Наряду с введенными, будем иметь дело с сеткой Пн и подсетками Пн,к, которые совпадают с сетками Пн и Пн,к, но не включают узлы, лежащие на границе сопряжения.
Примем следующие соглашения по нумерации элементов и узлов во введенных сетках. Для макроэлементов и макроузлов в структурированной сетке Пя принимается стандартная нумерация: слева направо, снизу вверх. В соответствии с ней нумеруются
внутренняя (Ск = С^), то вся ее граница совпадает с 7к, а если — граничная (Ск С^), то помимо 7к ее граница включает Гк — отрезок границы Г.
П ]%,к {хік Хі + ^к hx,k) Уік ^3 + Зк hy,k > Ік ° пх,к > І к ° пу,к\
с шагами
различать подсетки двух типов: внутренние П^к и граничные П^к. Подсетки П^к с
к
(2)
к=1
подобласти С к (подчеркнем, что внешние подобласти не входят в нумерацию по индексу к). В структурированных подсетках Пн,к также вводится стандартная нумерация. В сетке Шн узлы нумеруются следующим образом: сначала на горизонтальных отрезках координатных линий Пя, которые перебираются стандартным образом, а затем аналогично — на вертикальных отрезках.
2. Алгоритм решения задачи
Проведем декомпозицию расчетной области С на подобласти С к, принадлежащие макроэлементам Як сетки Пя, как это было сделано в предыдущем пункте. Алгоритм решения задачи (1) относительно функции и строится таким образом, что в первую очередь находится и|7 на границе 7 сопряжения подобластей, а затем, с использованием полученных значений в качестве граничных условий Дирихле, решаются краевые задачи в подобластях, т. е. определяется функция и во всей области. Основной проблемой, алгоритм решения которой предлагается в настоящей статье, является отыскание и|7. Перейдем к его изложению.
Запишем условия сопряжения на 7 в виде
где п — нормаль к 7, а знаки ± означают, что соответствующие величины относятся к подобластям, лежащим по разные стороны от 7. Условия (3), (4) должны выполняться на решении исходной задачи в С. Если же на 7 задана произвольная функция V, то условие (3) выполняется, а условие (4) — нет, причем нормальные производные в этом случае являются функциями, зависящими от V. Можно сказать, что в пространстве V функций, определенных на 7, задан оператор Г : V ^ V, ставящий в соответствие каждой функции V функцию f £ V по формулам
Тогда исходную краевую задачу (1) можно переформулировать следующим образом: в замкнутой области С требуется найти решение операторного уравнения
(3)
(4)
где функция и есть решение задачи
А и = дь I и|г = д2, и|7 = V.
Fv = 0
(5)
совместно с решением краевой задачи
А и = д!, 1и|г = д2, и|7 = V.
Это положение составляет суть предлагаемого подхода.
Формулы (6) определяют краевую задачу, которая редуцируется к решению подзадач
А и(к) = дъ и(к)\ = ^к) (7)
Я1 \1к
во внутренних подобластях С^ и подзадач
А и(к) = дъ 1и(к)\г = д2, и(к)\ = ^к) (8)
'г к ' 1к
в граничных подобластях С^, где и(к — искомая функция в подобласти С к, а v(k') — ее значения на границах сопряжения 7к. Основной проблемой, алгоритм решения которой излагается ниже, является решение операторного уравнения (5), т. е. отыскание функции V на границе 7 сопряжения подобластей, а точнее — функций v(k') на 7к.
Подзадачи (7), (8) решаются приближенно на подсетках Пн,к. Для этого они методом конечных разностей или конечных объемов [5] сводятся к следующим дискретным
подзадачам: \
Аь = ді, = ь1К> (9)
,(к)
1к
п(1)
на подсетках к и к подзадачам
Аьи1к) = ді, ки1к) г = д2, и1к) = (10) Г к
ь
1к
на подсетках М^к. Здесь и^, — приближенные значения функций и(к), ь(к) в узлах
подсеток Мь,к, &н,к, операторы Аь, 1ь — аппроксимации оператора Лапласа и оператора граничных условий на данных подсетках. Задача (6) при этом заменяется следующей приближенной задачей:
АЬ иЬ д1> 1Н иЬ|г д2) иь\1 (11)
к (к) к (к)
где иь = и иь , = У ), которая решается на сетке Мь.
Рассмотрим разность
ди \(+) / ди \(-)
и = 1 ап)г Чоп), ■ г = 1-п- <12>
в узлах сетки Шн. Пусть п — внутренняя нормаль к 7 в подобластях. Тогда соотношение (12) приобретет вид
ди \(+) ( ди \(-)
Л =1 еП)і ЧдГп), • г = ^
Заменим производные здесь конечно-разностными соотношениями
ди
т— ~ ^ин, (13)
дп) ь
к=1 к=1
где йн — линейный оператор, определяющий конкретную схему, которая по порядку точности должна быть согласована с порядком точности решения приближенных
задач (9), (10). При этом величины ^ заменяются их разностными аналогами Д,г, которые определяются как
Ь,г = (4+)ил+)) . + ., г =1,п (14)
где , й-- = й-. По аналогии с функциями f непрерывного аргумента величины fh,г являются функциями, зависящими от , т. е. fh,г = ^,г (юн,1, ю-,2, ■ ■ ■ , ю-,п). Более того, в силу линейности исходной задачи, fh,г - линейные функции, поэтому их можно представить как
п
Ь,г = ^2 Vh,j + Ьг, (15)
j=!
где аг^, Ьг — неизвестные пока коэффициенты. Запишем (15) в матричном виде
Ь = Аю- + Ь, (16)
где А = {аг^} — квадратная матрица, а fh = {fh,г}, V- = {ю-,г}, Ь = {Ьг} — векторы (г,3 = 1,п).
Вычисление компонент вектора Ь проводится по следующему алгоритму:
1) дадим V- пробное значение V- = е(0) = 0 и решим краевую задачу
Аь = д1, 1н = д2, ФГ =° (17)
(0)
относительно функции ф(0);
2) рассчитаем компоненты /ь,і = /^ вектора / = /^° согласно формуле (14), в которой вместо иь фигурирует Фь0). Тогда из (16) следует, что
Ь = /Г-
Расчет элементов матрицы А осуществляется по алгоритму, который также основан на том, что вектор ьь получает ряд пробных значений. Будем присваивать вектору ьь пробные векторы Ьь = Є(3 = |, компоненты которых равны единице в І-м узле
сетки шь и нулю — в остальных, т. е.
,(3) = Я. . Я. . = / 1
°, * = 3
еі = кз, кз = 1 п „• -I ■ , 1,3 = 1,п-
Тогда для каждого І = 1,п из решения краевых задач
,(з) - „ і ,„(з)'
Аь = дь 1ь № = д2> р(н =е(3) (18)
Г 7
относительно функций и вычисления /ь = /(3 согласно (14), можно определить А(3) — І-й столбец матрицы А по формуле
А(3) = /3 - Ь. (19)
Рассмотрим процесс этих вычислений подробнее. В силу линейности задач (18) пред-
(з)
ставим <рь в виде суммы
= фн] + ф(н , 3 = 1п',
7
где ф" 1 — решение задачи (17), а ф" — решения следующих краевых задач:
Дн Фі'
к ф'Л]
фЛ
Лз)
(20)
Тогда для компонент вектора будем иметь
№ = ^Г)г + (4Г' ^
и, следовательно, с учетом (19)
4+Ч3л)) . + (4-)фЛ-)) + ь
А+)фЛ+)) + (4-)фЛ-)
і,3
1, п.
(21)
Здесь, как обычно, верхние индексы ± обозначают операторы и функции по разные стороны от 7.
Задачи (17), (20) также сводятся к решению подзадач
п(1)
на внутренних подсетках \Гк к и подзадач
(Л,к)
н
(Л,к)
(22)
1к
дн фН3к) = *1Л Л, 1н фн
(Л,к)
Г*
Л
Аз > к)
(23)
1к
на граничных подсетках О^к относительно функций ф""’к^ = 0,п, к = 1, К. Здесь е{з,к) — подвекторы векторов е(з), состоящие из компонент, лежащих на 7к, т. е. е(з’к') = е3’к) = е(3 (Тг) , Тг Е 'ук г, а величины г2" определяются как
*1,з
0, з> 0,
9ъ 3 = 0,
*2,
0, 3> 0, 92, З = 0
Для вычисления элементов всех столбцов А(3 (] = 1,п) по формулам (21) необходимо, вообще говоря, решать п задач (20), каждая из которых сводится к решению К вспомогательные подзадач (22), (23), т. е. всего необходимо решить пК вспомогательных подзадач в подобластях. Столь большое количество вспомогательных задач значительно сокращается (до приемлемой величины) благодаря принятому способу построения квазиструктурированных сеток по следующим двум основным причинам.
Пусть отличная от нуля компонента вектора е(3 соответствует узлу, лежащему на отрезке 7; Е 7, который ограничен двумя соседними макроузлами сетки Он. Решение вспомогательных подзадач ф3'1 будет нулевым во всех подобластях за исключением двух подобластей, имеющих отрезок 7; своей общей границей. Отсюда следует, что для определения ненулевых значений ф^’к требуется решение только 2п подзадач — и это первая причина сокращения объема вычислений.
Рассмотрим две внутренние подобласти О^1, Оу). Напомним, что они равны между собой. Пусть Т", Т" — узлы подсеток Ш"’к, ^Н’к', определенных на границах 'ук, 1к‘ данных подобластей, в которых векторы е(з\ е(з 1 имеют компоненты, равные единице. Кроме того, будем считать, что данные узлы имеют равные координаты в локальных
0
0
Г
7
системах координат, введенных в подобластях одинаковым образом (например, центр локальной системы координат лежит в левом нижнем углу рассматриваемой подобласти). Тогда, в силу равенства подобластей, если подсетки в них равны между собой (П^к = П*^/), то равны и решения вспомогательных подзадач (22) (ф*,к') = ф* ,к )). Если же П**^ = ', то тем не менее все узлы редкой подсетки (допустим, это будут узлы
Тг, Є П*1^') являются общими с узлами густой подсетки Тг Є П*1^. Тогда в совпадающих
узлах Тг = Тг/ справедливо приближенное равенство ф^’к) (Тг) ф( ,к ) (Тг>). Отсюда следует, что, решив подзадачу на густой подсетке, мы без потери точности будем иметь решение и на редкой подсетке. Применительно ко всем внутренним подобластям это означает, что для получения решения в них достаточно решить серию следующих базовых вспомогательных задач относительно функции, которую мы обозначим через )ш.
Аа(н) =0, а(^(Тм) = Тм Є 7а, (24)
на самой густой подсетке П*1^ = П*1^. Здесь 70 — граница внутренней подобласти
Є Я0, в которой построена подсетка П*10; V, і — номера узлов подсетки ы*,0 —
проекции П(10 на 70 (без макроузлов). Дальнейшее сокращение объема вычислений связано с тем, что нет необходимости решать вспомогательные подзадачи (24) для всех узлов подсетки ы*,0: из соображений симметрии задачи и ы *,0 это достаточно выполнить для узлов, лежащих на двух половинах смежных сторон Я0. Все эти соображения составляют вторую причину сокращения числа вспомогательных задач.
Если же 7' — общая граница двух граничных подобластей О±±\ то вспомогательные задачи (23) в общем случае необходимо решать для каждой из них (см. по этому
(2)
поводу работу [2]). При этом, так как подсетки П*± являются квазисогласованными,
(?)
может возникнуть ситуация, изображенная на рис. 1. Если компонента єї пробного вектора е(? отлична от нуля в каком-либо из узлов, обозначенных целыми числами, то решение задач на подсетках П(2+, П(2- не вызывает затруднений. Если же е?) =0 в
одном из дробных узлов (например, в узле 1/2 компонента е? = 1), то для решения
(2)
вспомогательной задачи на редкой подсетке П* + предварительно проводится локальная
^ди^карм йанйшпоїасетйки^іс^осїїаййа&ил?дчи5езе5л-ижсат®ліе^ц сотв-ляютевеличины,
п(2)
1 /1/2 2 / ✓ / / / / / / / / / / / 3/2 п(2) іігА.+
Рис. 1. К решению вспомогательных задач в граничных подобластях
вычисленные по формуле (21) из решения ^3) = I% = 1,п| задачи (20), в которой в качестве условий на границе сопряжения подобластей берется вектор с компонентами, равными единице в ^-м узле сетки Шн и нулю — в остальных узлах. Строка с номером % — это величины, рассчитанные по формуле (21) по %-й компоненте решений ] = 1,п.
Лемма 1. Для элементов матрицы А справедливы соотношения
< 0, а^- > 0, = 1,п. (25)
Здесь неравенства понимаются асимптотически при Н 0.
Рассмотрим %-ю строку матрицы А. Диагональный элемент этой строки соответствует %-му узлу сетки Шн. Пусть данный узел лежит на отрезке 7' Е 7, который является общей границей двух для определенности внутренних подобластей С+, С— (см. рис. 2, на котором отрезок 7' обозначен цифрой 2). Ненулевые элементы %-й строки — это величины а^-, рассчитанные по формуле (21) для ] Е Мо, где Мо — множество номеров узлов сетки Шн, лежащих на сторонах подобластей Согласно принципу максимума [6] решения ^Н?,±) в окрестности ]-го узла будут убывающими функциями вдоль нормали п при ] = % и возрастающими — при ] = %. Тогда разностные производные в формуле (21) при достаточно малых Н будут отрицательными при ] = % и положительными при ] = %. Аналогичные рассуждения можно провести и для граничных подобластей &± , не забывая, что в данном случае под Мо понимаются номера узлов, лежащих только на границе сопряжения 7&. Так как а^- = 0 при ] / Мо, то лемма 1 доказана. □
Лемма 2. Матрица А имеет диагональное преобладание по строкам, т. е.
п
|а^| > ^ а^, % = 1,п. (26)
3=1
3 =
Пусть %-й узел лежит на стороне 2, общей для внутренних подобластей С(1) (рис. 2, а). Обозначим через Мг множества, состоящие из номеров узлов подсеток Шн,г на сторонах
о[1) 7 . о[2)
2 ( " 2
0а) б г- О?
1 г 1
а б
Рис. 2. К расчету элементов матрицы А: а —внутренние подобласти; б — граничные подобласти
___ 7
данных подобластей (шн,г С Шн, г =1, 7), а через М0 — их объединение: М0 = У Мг.
Г=1
Ясно, что аг ? = 0 при 3 / М0. Рассмотрим вначале случай, когда данные подобла-
п(1) - о(1)
сти покрываются одинаковыми подсетками ъгн ±, равными самой густой подсетке ъгн 0.
Представим элементы матрицы А в виде
аг,3 = а( + ’ + а( >3-)) 3 Е М0, (27)
где
431 = (^’фГ’) г ■
Здесь -фН>>±) — решения подзадач (22) в узлах подсеток • Так как в качестве ф3’*’
берутся базовые решения подзадач (24), которые можно обозначить как а(з,±’, то последнее равенство приобретает вид
«5 > = К^Г’). ■ (28)
В формуле (27) оба слагаемых отличны от нуля лишь при і Є М2. Если і Є М3, М5, М7, то а(+) = 0, а при і Є Мі, М4, М6 — а(.,-) = 0. Сумму 8, элементов і-й строки можно представить в виде
8, = 4+) +
где
*Т> = £ 3 (29)
3ем (±)
Здесь
М(+> = М1 и М2 и М4 и М6, М(-> = М2 и М3 и М5 и М7. (30)
Так как 4±>, фигурирующие в формуле (28), — линейные операторы, то, подставив (28)
в (29), вынесем их за знак суммы, что дает соотношения
*(±> = (4*44 • (3!)
где
„Н±> = £ 3>. (32)
3ем (±)
Из (32) следует, что аН±> являются решениями краевых подзадач, в которых искомая сеточная функция принимает значения, равные единице во всех узлах, лежащих на границах подобластей С±>. Такие задачи можно представить в виде
ДнСТ1±>(Г)=0, Т Е С±>, стн±>(г.) = !, т, Е шН±>, (33)
где
шН+> = шН,1 и шН,2 и шН,4 и шН,6 э Шн ’ = шН,2 и шН,3 и шН,5 и шН,7- (34)
Из (33) видно, что аН±> = 1. Тогда по определению оператора 4н (13) и на основании
соотношений (31) получим з(±> = 0. Последнее с учетом неравенств (25) из леммы 1
означает равенство в соотношениях (26), т. е. справедливость настоящей леммы в данном случае.
Если же хотя бы одна из подсеток (допустим, это будет ) не совпадает с - п(!)
подсеткой \Гк 0, то возможны следующие два варианта.
1. Все узлы подсеток — это по-прежнему узлы самой густой подсетки ^ ,о, и тогда проведенные рассуждения остаются в силе.
2. Хотя бы одна из подсеток ^,г (допустим, ^д) включает не все узлы густой подсетки ^,0, лежащие на стороне 1. Тогда относительно ) справедливы приведенные
выше рассуждения, а а^+)есть сеточная функция, принимающая на границе подобласти
^(1) (+)
значения, равные единице в совпадающих узлах подсеток ^ и ^,0 и равные нулю
в несовпадающих узлах. Разумеется, внутри подобласти ст^+) удовлетворяет сеточному уравнению Лапласа. Таким образом, данная функция есть решение подзадачи
Д„а<+)(Т)=0, Т 6 С+°, а'+)(Т,) = 1, Т 6 Ш<+), а'+)(Ту) = 0, Т 6 и*,о\4+)- (35)
Отсюда видно, что функция ст^+) < 1 в С+1) принимает максимальное значение, равное единице, в узлах ^+), в частности, в г-м узле. Поэтому в окрестности г-го узла данная функция является убывающей вдоль нормали п, разностные производные в формулах (31) — отрицательными, следовательно, з(+) < 0, что означает справедливость строгого неравенства (26).
Пусть теперь одна из подобластей 0± или обе сразу являются граничными 0± = С±2) (рис. 2, б). В этом случае также справедливы рассуждения, проведенные выше для внутренних подобластей, со следующими отличиями.
1. Множества М(±) (см. (30)) определяются как
М(+) = М1 и М2 и М4, М(-) = М2 и М3 и М5.
2. В формулах (31), (32) вместо функций а;,±), а|1±) фигурируют функции ^^?’±), ^1±). Функции ^^;,±) определены выше, а ^^±) — это решение подзадач
Дь#)(Т) = 0, Т 6 С±2), /^(±)(Т) = 0, Т 6 Г±, ^(Т,) = 1, з 6 М(±).
С учетом данных замечаний и из принципа максимума следует, что з(±) < 0, а это в свою очередь означает справедливость строгого неравенства (26) для граничных подобластей.
Лемма 2 доказана полностью. □
Из лемм 1 и 2 следует, что (—А) — матрица положительного типа и А — невырождена [6].
Возвратимся теперь к задаче отыскания функции V на границе сопряжения подобластей, т. е. к решению операторного уравнения (5). Будем решать ее приближенно относительно функции ^ исходя из того, что должны выполняться условия сопряжения
для разностных аналогов (13) нормальных производных. Данное требование на основании изложенного выше означает, что вектор Д, компоненты которого определяются по формуле (14), должен быть нулевым, т. е. должно выполняться равенство
А^ + Ь = 0. (36)
Выражение (36) представляет собой систему линейных алгебраических уравнений относительно вектора Vh = компоненты которого суть значения Vh в узлах сетки
По определению данная сетка не содержит макроузлов Tp. Для отыскания в них значений искомой функции дополним (36) системой уравнений
(Ahnh)p = gi (Тр) , p =1, p, (37)
где Дh — аппроксимация оператора Лапласа на сеточном шаблоне, включающем узлы сетки ^h и один из узлов Tp, а П, — значения искомой функции. Если оператор Ah, определенный в подобластях, не содержит макроузлов, что соответствует, например, обычной пятиточечной схеме, то система уравнений (36), (37) является распадающейся и тогда вычисления по формулам (37) проводятся после того, как получено решение системы (36).
Из лемм 1 и 2 следует, что верна теорема.
Теорема. Решение системы линейных алгебраических уравнений (36) существует и единственно.
Подведем итоги. Предлагаемый алгоритм решения краевой задачи (1) на квазиструк-турированной сетке Пh (2) состоит из следующих этапов.
1. Формирование системы линейных алгебраических уравнений (36).
1.1. Вычисление вектора b. Для этого решаются краевые подзадачи (22), (23) с нулевыми значениями искомых функций на границе сопряжения, рассчитываются нормальные производные по формуле (13) и компоненты вектора /h° по формуле (14).
1.2. Расчет элементов матрицы A. Для этого выполняются следующие действия.
1.2.1. Решается серия базовых вспомогательных подзадач (24) в одной внутренней подобласти.
1.2.2. Решается серия вспомогательных подзадач (23) в каждой граничной подобласти.
По решениям данных подзадач вычисляются нормальные производные (13) и элементы матрицы A по формуле (21).
2. Вычисление приближенных значений искомой функции Vh на границе сопряжения подобластей. Для этого выполняются следующие действия.
2.1. Решается система линейных алгебраических уравнений (36), в результате чего определяются значения Vh в узлах сетки ^h.
2.2. Проводятся вычисления Vh в макроузлах по формуле (37).
3. Решение краевых подзадач (9), (10), в которых в качестве условий на границе сопряжения задаются значения Vh, в результате чего находится приближенное значение Uh искомой функции u.
Список литературы
[1] Ильин В.П., Свешников В.М., Сынах В.С. О сеточных технологиях для двумерных краевых задач // Сиб. журн. индустр. математики. 2000. Т. 3, № 1. С. 124-136.
[2] Свешников В.М. Оценки эффективности прямого метода декомпозиции // Вычисл. технологии. 2008. (в печати).
[3] Федоренко Р.П. Введение в вычислительную физику. М.: МФТИ, 1996.
[4] МАцокин А.М. Автоматизация триангуляции областей с гладкой границей при решении уравнений эллиптического типа. Новосибирск: ВЦ СО АН СССР, 1975. Препринт. № 15.
[5] Ильин В.П. Методы конечных разностей и конечных объемов для эллиптических уравнений. Новосибирск: ИВМиМГ (ВЦ) СО РАН, 2001.
[6] Самарский А.А. Теория разностных схем. М.: Наука, 1977.
Поступила в редакцию 25 апреля 2007 г., в переработанном виде —11 октября 2007 г.