Введение к работе
1.1 Актуальность темы
Существует широкий класс задач, решение которых содержит резкие неоднородности, проявляющиеся на мелких по отношению к размеру области пространственных масштабах. Численное решение таких задач сеточными методами требует специальных сеток для разрешения особенностей. Для этого необходимо использовать либо адаптивные к решению сетки, сгущающиеся в окрестности особенностей, либо достаточно мелкие сетки с шагом h и огромным количеством точек. Первый вариант требует применения специальных алгоритмов, второй - соответствующей памяти ЭВМ. В то же время проявления особенностей зачастую являются локальными, сосредоточенными в мелкомасштабных подобластях. Подтверждением тому является характерный вид функций Грина и типичных решений многих задач математической физики, а также физические эффекты такие, как принцип Сен-Венана в теории упругости и другие. Наличие областей сосредоточения неоднородностей позволяет ввести сетку с характерным размером Я 3> Л, узлы и ребра которой проходят по участкам относительной гладкости решения. При этом сетка размером Я заведомо не позволит разрешить особенности при использовании обычных численных методов, но зато число ее узлов достаточно мало.
Для решения подобных задач на сетках размером Н в работах Л.Г. Страховской и Р.П. Федоренко1 был предложен метод конечных суперэлементов (МКСЭ). Метод появился более 25 лет назад и использовался при решении ряда сложных задач диффузии, теории упругости, кинетики ядерных реакторов и других.
Метод конечных элементов (МКЭ) основан на представлении решения задачи в виде разложения по системе базисных функций, имеющих конечный носитель. При этом мера таких носителей предполагается малой (сетка h) и стремящейся к нулю. Базисные функции в МКЭ берутся в виде функций сравнительно простой структуры, как правило, полиномиальной. МКСЭ также основан на представлении решения задачи в виде разложения по системе базисных функций, имеющих конечный носитель. Однако в случае МКСЭ мера таких носителей (сетка Я) не предполагается стремящейся к нулю и столь велика, что она заведомо не позволяет (при использовании МКЭ) передать особенности решения. Другое отличие касается построения базисных функций. В МКСЭ базисные функции строятся для данной рассматриваемой задачи специальным образом так, чтобы в них самих содержалась значительная информация о решении задачи. Именно специальный, под задачу, выбор базисных функций и позволяет с помощью очень грубого разбиения исходной области получить хорошее численное решение.
Несмотря на свой возраст, теоретически метод исследован сравнительно слабо. Обоснование одного варианта метода было предложено в работах В.В. Репяха2
В настоящей работе предложен теоретический алгоритм, позволяющий строить и исследовать аппроксимации МКСЭ для достаточно широкого класса задач. Помимо МКСЭ рассмотрены комбинированные аппроксимации метода конечных элементов и метода конечных суперэлементов.
Укажем причины, по которым для исследования метода конечных суперэлементов не может непосредственно применяться техника исследования обычного метода
1 Федоренко Р.П. Введение в вычислительную физику. М.: Издательство МФТИ, 1994. 528 с. 3Ртях В.В. Применение одного варианта метода конечных суперэлементов к решению задач теории упругости // ЖВМиМФ 1986. Т.26. №11. с. 1643-165»,- Pra.ix В В, Anartm ошибок метода приближенных суперэлементов // ЖВМиМФ. 1989. Т.30^№7СУСЯ{К94«ОНАЛЬИля
УДОЙГ
БИБЛИОТЕКА С. flereptor V//.-"' 09
конечных элементов.
Алгоритмически эти методы сходны. Во всех случаях приближенное решение ищется как линейная оболочка некоторой системы базисных функций. Конечномерная задача (система линейных алгебраических уравнений для узловых значений конечных элементов/суперэлементов) во всех случаях формально получается из условия ортогональности невязки приближенного решения некоторому конечномерному подпространству. Отличие состоит в том, что в случае МКЭ базисные функции задаются (исходя из некоторых условий), а в случае МКСЭ (или комбинированного подхода) все (или некоторые) базисные функции рассчитываются как точные решения рассматриваемой задачи.
Несмотря на внешнее сходство, эти методы качественно различны с точки зрения их теории. Остановимся на этом подробнее.
Предположим, что ищется решение и следующей задачи:
ueW: a(u,v) = f(v) S/veW,
где W - некоторое гильбертово пространство, о - билинейная непрерывная положительно определенная форма в пространтве WxW, f - линейная непрерывная форма на W. Такой вид имеют слабые постановки большого количества задач математической физики.
Пусть приближенное решение задачи ищется как элемент конечномерного пространства Wk С W. Пространство W/, является линейной оболочкой той или иной системы базисных функций.
Приближенное решение uj, определяется как решение следующей конечномерной задачи:
щ Wh : o(uA) vh) - f(vh) V«fc Є Wh.
В методе конечных элементов пространство W), выбирается так, чтобы его элементы аппроксимировали произвольный элемент пространства W. Другими словами, требуется наличие опенки для ошибки наилучшего приближения вида
Vtu Є W : inf |1їо — tBfcH < e(w, h), и>»єи\
где є(ш, Л) - оценка ошибки интерполяции элемента w функцией ш/, Є Wh, e(w, Л) -v О при стремлении параметра дискретизации h к нулю. Последовательность пространств {Wh} при этом называется предельно плотной в W.
В соответствии с леммой Сеа8 ошибка приближенного решения оценивается сверху величиной ошибки наилучшего приближения решения элементом пространства Wh,
\\и - иЛ|| < С inf ||» - will ^ Сє(щ h), ч>»єи\
где uh - приближенное решение задачи. Таким образом, оценка ошибки приближенного решения сводится к оценке ошибки интерполяции произвольной функции из W элементами пространства Wh- Отсюда вытекает одно го основных требований к базисным функциям - они должны обладать аппроксимирующими свойствами.
Это требование не будет выполняться, если пытаться использовать для исследования МКСЭ рассмотренную процедуру непосредствено. В самом деле, в МКСЭ базисные функции не задаются, а рассчитываются как точные решения задачи. Поэтому
гСъярм Ф. Метод конечных элементов для эллиптических задач. М.:Мир, 1980. 512 с, с. 109
в общем случае они не будут обладать какими-либо аппроксимирующими свойствами и для них величина e(w, h) уже не будет стремиться к нулю для произвольного wW.
Описанное выше препятствие теоретического характера можно обойти. А именно, обратим внимание на то, что при построении аппроксимаций МКСЭ по существу задаются не сами базисные функции (они, как указано выше, рассчитываются специальным образом), а следы этих функций на границах суперэлементов, т.е. некоторые граничные базисные функции. Эти базисные функции уже обладают аппроксимирующими свойствами в пространстве функций, заданных на границе. Бели записать задачу относительно них, исключив из рассмотрения «внутренности» суперэлементов, то к ней уже можно применить описанный выше подход с использованием леммы Сеа и другого теоретического аппарата теории вариационных уравнений и проекционных методов. При этом ошибка приближенного решения будет определяться ошибкой интерполяции произвольной функции, заданной на границе суперэлементов, элементами соответствующего конечномерного пространства функций, также заданных на границе.
Реализации этого подхода для ряда задач и посвящена настоящая работа.
Отметим, что термин «суперэлемент» известен в теории метода конечных элементов и не в связи с методом конечных суперэлементов Р.П. Федоренко. Обычно он употребляется для обозначения группы конечных элементов, рассматриваемых совместно. Поясним это на примере.
Рассмотрим некоторую расчетную область, в которой решается задача. Пусть в этой области задана некоторая триангуляция, на которой определены финитные базисные функции. Задача определения узловых значений конечных элементов сводится к решению системы линейных алгебраических уравнений вида
Au = f,
(1)
где А - матрица жесткости задачи, и - вектор узловых значений конечных элементов, / - вектор правой части.
Выделим в этой области некоторую группу соседних конечных элементов, занимающих подобласть S расчетной области. Пусть и„ щ - векторы узловых значений конечных элементов, соответствующие внутренним и граничным узлам области S соответственно, Uo - вектор узловых значений, соответствующий оставшимся узлам расчетной области.
Тогда задача (1) может быть записана в следующем блочном виде:
(2)
Эта система уравнений совпадает с системой (1) с точностью до перестановки строк. Нулевые блоки в матрице этой системы появляются вследствии финитности базисных функций. Считая, что матрица А„ невырожденная, выразим из второй строки этой системы уравнений вектор и,. Получим
«» = KU> - А~1А,ьЩ-
Подставляя это выражение обратно в систему (2), получим:
Аьь - аь.а,;ал л» _ Щ _ /»-«*»„-/, ,^
Таким образом, можно понизить порядок системы уравнений для определения узловых значений конечных элементов, заранее исключив из нее неизвестные, соответствующие части узлов. Такой подход, когда группа конечных элементов рассматривается как одно целое, называется методом суперэлементов4, а указанная группа конечных элементов называется суперэлементом и может рассматриваться как самостоятельный объект при построении аппроксимаций задачи.
В рассмотренном выше примере был всего один суперэлемент S. В общем случае их может быть несколько, и они могут покрывать всю расчетную область. В этом случае система (3) будет связывать только неизвестные, соответствующие узлам на границах суперэлемевтов.
Рассмотренный выше подход также иногда называют методом разделения (декомпозиции) области для конечномерных задач Известны его модификации, когда системы базисных функций для узлов на границах суперэлементов выбираются специальным образом, отличным от способа выбора базисных функций для внутренних узлов суперэлементов5.
Для полноты отметим, что в последнее время активно развиваются и другие методы численного решения задач на основе представления решения в виде разложения по системе базисных функций, в свою очередь являющихся решениями специальных вспомогательных задач для исходного оператора (RFB - методы на элементах с нулевой невязкой)6. Очевидно их родство с МКСЭ.
Также к родственным подходам можно отнести методы декомпозиции области. Эти методы можно расматривать как итерационные методы решения некоторых уравнения для следов решения на границах некоторых подобластей7.
Отметим также метод наименьших квадратов на границе и метод Треффтца8. В этих методах решение ищется в виде линейной комбинации функций, каждая из которых является точным решением исходной задачи. Неизвестные коэффициенты в этой линейной комбинации определяются так, чтобы граничные условия на границе расчетной области выполнялись в некотором наилучшем смысле. В отличие от МКСЭ здесь не происходит разбиения области на меньшие подобласти, решение и базисные функции определены сразу во всей расчетной области. Функции, входящие в указанную линейную комбинацию, обычно задаются явно в виде некоторых алгебраических или тригонометрических многочленов.
1.2 Цели работы
Диссертация посвящена дальнейшей разработке и развитию метода конечных суперэлементов Р.П. Федоренко и его применению для решения задач математической физики в физически и геометрически неоднородных областях.
4Зарубин B.C., Селиванов В.В. Вариационные и численные методы механики сплошной среды. М.: Изд-во МГТУ им. Н.Э. Баумана, 1993. 360 с.
ЪВ.И. Агошков, Методы разделения области: некоторые результаты теории и приложения // М.: Отдел вычислительной математики АН СССР, 1990, 40 с.
eBrezzi F., Franca L.P., Rusao A. Futher consideration on residual-free bubbles for advective-difrusion equation // Comput. Methods Appl. Mech. Engrg. 1998. №166. p. 25-33., Franca L.P., Russo A. Approximation of the Stokes problem by Residual-Free Macro Bubbles // East-West J. Appl. Math. 1996. №4. p. 265-278., Brezzi F, Hughes T.J R, Marini L D, Ruseo A A priory error analysis of residual-free bubbles for advective-diffusion problems // SIAM J. Numer. Anal. 1999. V.36. №6. p. 1933-1948., Жуков В. T, Новикова Н.Д., Страховская Л.Г., Федоренко Р.П., Феодоритова О-В. Метод конечных элементов в задачах конвекции-диффузии. // Препринт ИПМ РАН. М., 2001. №8.
1 Марчук Г.И. Введение в вычислительную математику. М.:Наука, 1989. 608 с.
*Ректорис К Вариационные методы в математической физике и технике, М.:Мир, 1985, 590 с.
В работе разработан теоретический подход для построения и исследования метода конечных суперэлементов, а также исследована эффективность метода при решении ряда задач.
В работе МКСЭ применен для численного решения таких задач, как задача о скважине для уравнения Лапласа в двумерной области, задача о скоростном скин-слое в пространственно двумерном случае и задача теории упругости композиционных материалов в пространственно трехмерном случае. Для этих задач построены и исследованы аппроксимации метода конечных суперэлементов, проведено их численное исследование.
1.3 Научная новизна
Основой работы является идея замены исходной краевой задачи эквивалентной ей задачей для определения следов неизвестного решения на границах суперэлементов с помощью граничных операторов Пуанкаре - Стеклова.
Расчетная область при этом разбивается на некоторое количество непересекающихся подобластей-суперэлементов. Далее на основе граничных операторов Пуанкаре-Стеклова и формулы Грина, соответствующих оператору задачи, строится обобщенная постановка (вариационное уравнение) для определения следов решения исходной задачи на границах гуперэдементов. Использование операторов Пуанкаре-Стеклова позволяет исключить из рассмотрения «внутренности» суперэлементов. Для построения аппроксимаций указанного вариационного уравнения для следов могут применяться стандартные подходы и методы теории абстрактных вариационных уравнений и проекционных методов*. Это позволяет формально и единообразно рассматривать различные варианты МКСЭ, соответствующие тем или иным проекционно-сеточным методам, получить оценки ошибок и т.д. Таким образом МКСЭ «вкладывается» в известную и хорошо разработанную теорию.
При сведении исходной задачи к задаче для следов используются операторы Пуанкаре-Стеклова. Они были предложены в работах Лебедева и Агошкова10 как средство теоретического исследования методов декомпозиции области. При этом методы декомпозиции области рассматриваются как итерационные методы решения соответствующих уравнений для следов. Это позволяет использовать при их исследовании теорию итерационных методов решения абстрактных операторных уравнений. В соответствии с описанным выше подходом метод конечных суперэлементов может рассматриваться как проекционный метод решения этих уравнений, что позволяет использовать при их исследовании известную теорию проекционных методов.
Наряду с МКСЭ в работе рассмотрен комбинированный подход, в котором для аппроксимации задачи используется как МКСЭ, так и МКЭ. В этом случае в некоторых конечных элементах используются обычные аппроксимирующие базисные функции
90бэн, Ж.-П. Приближенное решение эллиптических краевых задач. М.: Мир, 1977 384 с, Съяр-ле Ф. Метод конечных элементов для эллиптических задач. М.:Мир, 1980. 512 с, Марчук ГИ, Агошков В.И. Введение в проекционно-сеточные методы. М.: Наука, 1981. 416 с, Красносельский М.А., Ваймикко Г.М., Забрейко П.П., Рутщкий Я.Б., Стеценко В Я Приближенное решение операторных уравнение. М.: Наука, 1969. 456 с.
10 Агошков В.И., Лебедев В.И. Операторы Пуанкаре-Стеклова и методы разделения области в вариационных задачах // Вычислительные процессы и системы. Т.2. М.: Наука, 1985, Агошков В.И Методы разделения области в задачах математической физики // Сопряженные уравнения и алгоритмы возмущений в задачах математической физики. М.: ОВМ АН СССР, 1989, Лебедев В И., Агошков В.И. Операторы Пуанкаре-Стеклова и их приложения в анализе. М: ОВМ АН СССР, 1983., Марчук ГИ Введение в вычислительную математику. М.:Наука, 1989. 608 с.
(как в обычном методе конечных элементов), а в некоторых - супералементные базисные функции. В этом случае в расчетной области также выделяется некоторое количество подобластей-суперэлементов, но они уже не покрывают всю расчетную область. В подобластях, занятых суперэлементами, осуществляется переход к рассмотрению следов решения на границах суперэлементов. В части области, не занятой суперэлементами, используются обычные аппроксимации метода конечных элементов.
При построении оценок приближенного решения в этом случае используется как обычный конечно-элементный подход (в тех подобластях, где используются обычные базисные функции) так и суперэлементный подход (в тех подобластях, где используются суперэлементные базисные функция). Переход от исходной задачи к задаче определения следов осуществляется локально, лишь там, где это необходимо.
Отметим, что сведение задачи к задаче для следов для МКСЭ или комбинированного подхода необходимо лишь для построения расчетной схемы и теоретического исследования метода. Это позволяет использовать при обосновании этих методов готовую и хорошо разработанную теорию. С точки зрения алгоритма построения конечномерной задачи эти методы сходны с обычным методом конечных элементов.
Подчеркнем отличия метода конечных элементов и метода конечных суперэлементов. При построении аппроксимаций МКСЭ задача записывается относительно следов решения на границах суперэлементов. При этом суперэлементная сетка не является «разностной» в обычном смысле этого слова. Это разбиение области на меньшие подобласти, оно происходит на этапе построения уравнения для следов, до построения конечномерных аппроксимаций. На границах суперзлементов задается некоторая разностная сетка, на которой задаются те или иные граничные базисные функции. Эти базисные функции должны обладать аппроксимирующими свойствами в подходящем пространстве следов. Шаг h граничной сетки является параметром дискретизации, приближенное решение задачи сходится к точному, когда h — 0. Диаметр суперэлементов не зависит от шага h разностной сетки и не меняется при fc->0.
1.4 Практическая ценность
Разработанный в работе подход позволяет единообразно и формально строить и исследовать аппроксимации метода конечных суперэлементов для большого класса задач математической физики. В работе рассмотрена модельная задача о скважине для уравнения Лапласа. Также рассмотрены задачи, имеющие важное прикладное значение. Одной из таких задач является задача о скоростном скин-слое, возникающая при математическом моделировании электродинамических ускорителей типа «рельсотрон». Помимо этого, рассмотрена задача теории упругости композиционных материалов. Предложенный в работе подход может быть формально распространен на большой класс эллиптических задач с оператором дивергентного вида. Существенным здесь является наличие для конкретной задачи формулы Грина (некоторого соотношения, связывающего интегрирование по объему с интегрированием по границе) и оператора Пуанкаре-Стеклова, описывающего реакцию решения задачи «в целом» на внешнее воздействие на границе расчетной области. Отметим, что существование формулы Грина для симметричных положительно определенных операторов в гильбертовом пространстве является следствием общей теории («абстрактная» формула Грина11). Теория же операторов Пуанкаре-Стеклова хорошо разработана в связи с
пОбэн, Ж.-П Приближенное решение эллиптических краевых задач. М.: Мир, 1977. 384 с.
обоснованием методов декомпозиции области.
Разработанный подход также применим и для построения расчетных схем МКСЭ для нестационарных параболических задач. В этом случае задача сначала может аппроксимироваться только по времени (метод прямых или метод Роте12), и возникающие на каждом временном слое эллиптические задачи решаются с помощью метода конечных суперэлементов.
Метод конечных суперэлементов входит в класс методов, в которых решение исходное задачи сводится к решению серии более простых задач, например, задач в областях более простой формы. Методы данного класса, например, методы разделения области, активно исследуются в настоящее время в связи с появлением эффективных алгоритмов решения краевых задач в областях простой формы и возможностью эффективной реализации алгоритмов этих методов на многопроцессорных и параллельных ЭВМ. В работе рассматриваются некоторые результаты по реализации МКСЭ на вычислительных машинах с параллельной архитектурой.
1.5 Аппробация работы
Результаты работы докладывались на семинарах ИПМ ям. М.В. Келдыша РАН и на следующих конференциях: «Студенческая научная веснаг2001», 2001, МГТУ им. Н.Э. Баумана, Москва; Второй международный конгресс студентов, молодых ученых и специалистов «Молодежь и наука-третье тысячелетие»/YSTM'02, 15-19 апреля, 2002, Москва, Россия; 7th International Conference «Mathematical Modelling and Analysis», May 31-June 2, 2002, Kaariku, Estonia; First International Conference «Computational Methods in Applied Mathematics CMAM-1», July 20-24, 2003, Minsk, Belarus; International Conference «Mathematical Modelling and Analysis MMA-2004», May 26-29,2004, Jurmala, Latvia; Международная конференция по вычислительной математике МКВМ-2004, 21-25 июня, 2004, Новосибирск, Академгородок, Россия.
1.6 Публикации
Результаты выполненной работы представлены в 12 печатных работах (см. раздел «Публикации автора по теме диссертации»).
1.7 Структура диссертации
Диссертация состоит из введения, трех глав и заключения. Работа изложена на _iZ& страницах и содержит Чт рисунков, а также билиографию из 62. ссылок.