Метод конечных элементов

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

Метод конечных элементов - это численный метод решения дифференциальных уравнений. В данной статье будет показан метод подсчёта и приведён код соответствующей программы на C++ (g++ 5.1.0; linux 4.0.x-ARCH). Для работы также потребуется MATLAB 4.0+ с возможностью вызова из терминала.

Суть[править]

Суть метода следует из его названия. Область, в которой ищется решение дифференциальных уравнений, разбивается на конечное количество подобластей (элементов). В каждом из элементов произвольно выбирается вид аппроксимирующей функции. В простейшем случае это полином первой степени. Вне своего элемента аппроксимирующая функция равна нулю. Значения функций на границах элементов (в узлах) являются решением задачи и заранее неизвестны. Коэффициенты аппроксимирующих функций обычно ищутся из условия равенства значения соседних функций на границах между элементами (в узлах). Затем эти коэффициенты выражаются через значения функций в узлах элементов. Составляется система линейных алгебраических уравнений. Количество уравнений равно количеству неизвестных значений в узлах, на которых ищется решение исходной системы, прямо пропорционально количеству элементов и ограничивается только возможностями ЭВМ. Так как каждый из элементов связан с ограниченным количеством соседних, система линейных алгебраических уравнений имеет разрежённый вид, что существенно упрощает её решение.

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

Матрица жёсткости[править]

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

Обычно матрица жёсткости получается разреженной, то есть содержащая большое количество нулей. Для работы с подобным типом матриц созданы специальные библиотеки (mtl4, SparseLib++, SPARSPAK и другие)

Определение[править]

Элементы матрицы жёсткости в общем случае равны

Например, если дано уравнение Пуассона

в пространстве и граничные условния — это

Представим функцию как ряд:

— это известные значения функции в узлах, а — некие базисные функции

то

Создание матрицы[править]

Для одного треугольника[править]

Пусть дан один конечный элемент, для простоты — треугольный. Матрица жёсткости, по сути, задаёт связи между узлами. Так как у элемента три узла (в локальной нумерации — 0, 1 и 2), то матрица будет иметь вид

В дальнейшем матрицу для одного треугольника будем называть локальной, для всей сетки сразу - глобальной.

В общем случае, элементы определяются через линейные функции

где
— площадь треугольного элемента.

и получаются из цикличной перестановкой индексов. Удобно искать как определитель матрицы

Сами

В описываемом случае для каждого треугольника составляется такая матрица:

Первый вид обобщения на несколько треугольников - сшивание[править]

Stiffness matrix - adding triangles.png

Для того, чтобы сделать из многих раздельных матриц, полученных выше, одну большую матрицу, описывающую отношения между узлами всей области расчёта, необходимо произвести процедуру объединения матриц. Пусть символ обозначает разделённые элементы (рис. а), а символ — объединённые элементы (рис. б).

Обозначим — вектор-строку значений функции в вершинах двух треугольников (см.рис). Символ обозначает транспонирование матрицы То есть, это вектор значений функции в шести узлах треугольников. Очевидно, что при объединении оных получится вектор , содержащий только четыре компоненты.

Преобразование происходит по схеме

Нумерация, конечно же, произвольная: необходимо равенство функции в соответствующих вершинах. Матрицу называют матрицей преобразования, а само уравнение называют связанной системой.

Запишем теперь матрицу жёсткости для двух треугольников:

Результирующая матрица

То есть, на каждом следующем шаге необходимо добавлять новые элементы к уже существующим.

Второй вид обобщения на несколько треугольников - дозаписывание[править]

Mesh FEM.png

Пусть есть область, представленная и разбитая на треугольники так, как преставлено на рисунке. Пусть данная сетка содержит узлов. Создадим глобальную матрицу (размера, очевидно, ) и заполним её нулями. Начнём строить локальные матрицы для треугольников, например, для

Введём локальную нумерацию для данного треугольника: пусть его верхняя вершина имеет локальный номер , далее по часовой стрелке и . Иначе говоря, пусть глобальным номерам соответствуют локальные номера соответственно.

Составим матрицу для этого треугольника так, как описано выше, получив что-то типа

Теперь заменим локальную нумерацию на глобальную. То есть запишем локальное число как глобальное число , - как , - как и так далее.

Получим

С остальными треугольниками поступаем аналогично. Необходимо помнить, что надо "дописать" число в глобальную ячейку, то есть прибавить к уже существующему.

Матрица масс[править]

Матрица масс собирается по тем же правилам, но чуть по другим формулам. Создаётся матрица размеров три на три, затем говорится, что

где
  • - площадь данного треугольника, которая считается, как в предыдущей главе.
  • и получаются из циклической перестановкой, равно как второй и третий блок элементов матрицы из первого.

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

Учёт граничных условий[править]

Условия Дирихле[править]

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

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

Условия Неймана[править]

В случае граничных условий второго рода глобальная матрица не меняется.

Код[править]

Из-за большого объёма кода он был помещён на подстраницу.

Литература[править]

  • Zienkiewicz, K. Morgan — Finite elements and approximation, 1983.
  • P.P. Silvester, R.L. Ferrari — Finite elements for electrical engineers, 1986.
  • E. Suli — Finite Element Methods for Partial Differential Equations, 2011
  • K.J. Bathe, E.L. Wilson — Numerical methods in finite element analysis, 1982