Модели и методы анализа проектных решений

МЕТОД КОНЕЧНЫХ ЭЛЕМЕНТОВ

Л. А. РОЗИН

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

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

Допустим, требуется построить такое семейство функций u(x) при a # x # b. Интервал ab разбивается на конечное число частей (элементов), соединяющихся между собой и с концами интервала в узловых точках (узлах) xi (рис. 1). В пределах каждого элемента задается функция, например в виде линейного полинома. Она определяется своими значениями u(xi) в узлах на концах элемента. Если отыскиваемая функция является непрерывной, то значения ее в каждом узле для соседних элементов совпадают. В результате имеем семейство кусочно-линейных непрерывных функций, которые изображаются в виде ломаных и определяются конечным числом параметров - своими узловыми значениями. На рис. 1 показана одна из функций такого семейства. Здесь 5 элементов, 6 узлов и 6 узловых параметров u(xi) = ui . В случае нескольких переменных схема метода конечных элементов в принципе не меняется. Таким образом, метод конечных элементов заменяет задачу отыскания функции на задачу отыскания конечного числа ее приближенных значений в отдельных точках-узлах. При этом если исходная задача относительно функции состоит из функционального уравнения, например дифференциального уравнения с соответствующими граничными условиями, то задача метода конечных элементов относительно ее значений в узлах представляет собой систему алгебраических уравнений.

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

1. Метод конечных элементов позволяет  построить удобную схему формирования  системы алгебраических уравнений  относительно узловых значений  искомой функции. Приближенная аппроксимация  решения при помощи простых  полиномиальных функций и все  необходимые операции выполняются на отдельном типовом элементе. Затем производится объединение элементов, что приводит к требуемой системе алгебраических уравнений. Такой алгоритм перехода от отдельного элемента к их полному набору особенно удобен для геометрически и физически сложных систем.

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

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

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

Первая работа, где рассматривалась схема типа метода конечных элементов, принадлежит известному математику Р. Куранту [2]. Построение метода с использованием физических соображений и его название "метод конечных элементов" содержатся в статье, написанной инженерами [3]. Такое сочетание специальностей авторов характерно для работ по методу конечных элементов. В последующем было опубликовано много статей и книг, посвященных этому методу и его различным модификациям. Некоторое представление об этом можно получить из списка литературы, например в [4]. Метод конечных элементов реализован в больших универсальных компьютерных пакетах программ, которые имеют широкое применение.

ПРИМЕР ДИСКРЕТНОЙ ЗАДАЧИ

Обратимся к дискретной задаче, состояние которой точно определяется конечным числом параметров. Рассмотрим упругий стержень в виде прямого кругового цилиндра, длина которого значительно больше его диаметра. Это позволяет отождествить стержень с его осью. Пусть три таких стержня расположены на оси x и соединены между собой в точках 1 и 2 (рис. 2, а). Точки a и b закреплены, что условно изображено на рисунке. К осям стержней вдоль x приложим внешнюю нагрузку. Очевидно, точки на осях стержней перемещаются вдоль x. Силы и перемещения считаются положительными, если они направлены в положительном направлении x. Задача состоит в определении перемещений точек, принадлежащих осям стержней и продольных внутренних сил в поперечных сечениях стержней.

Согласно методу конечных элементов, представим стержневую систему в виде элементов, соединенных в узлах. В качестве элементов примем отдельные стержни, а узлов - точки 1 и 2. На рис. 2, а в скобках указаны номера элементов. Обратимся к типовому для данной системы элементу r. На элемент r с узлами i, j (рис. 2, б ) может действовать распределенная нагрузка интенсивности q(r)(x) и перемещения его узлов Примем x = 0 в узле i и обозначим длину элемента l (r). Получим задачу для функции u(r)(x) перемещений точек оси r. Бесконечно малая часть элемента dx находится в равновесии под действием нагрузки q(r)(x)dx и продольных внутренних сил N (r)(x), действующих так, как показано на рис. 2, в. Из условия равновесия dx имеем

Согласно закону Гука, для упругого стержня

где c(r) > 0 носит название продольной жесткости стержня и определяется из опыта. Пусть c(r) = const для элемента r. Подставляя (2) в (1) получим задачу относительно u(r)(x) в виде дифференциального уравнения и граничных условий

В нашем примере для случая дискретной задачи положим q(r) = 0 и будем считать, что на стержневую систему действуют только сосредоточенные силы P1 и P2 соответственно в узлах 1 и 2. Тогда решение (3) примет вид

На основании (4) можно заключить, что состояние типового элемента r, то есть u(r)(x), N(r), точно определяется двумя параметрами - перемещениями его узлов что делает задачу дискретной.

Рассмотрим внутренние силы действующие в узлах i, j на элемент r. Поскольку имеет место линейная задача, то они линейно зависят от

Здесь (l = i, j; t = i, j) есть внутренняя сила действующая на элемент r в узле l и возникающая от единичного перемещения узла t. При этом перемещение другого узла равно нулю. На рис. 3, а показана сила для сжатого элемента при

Соотношения (5) можно представить в матричной форме. Введем столбцы f (r), u(r) и матрицу K(r)

Тогда (5) можно записать в виде

f (r) = K(r)u(r).

Для упругой пружины коэффициент пропорциональности между силой и перемещением называется коэффициентом жесткости пружины. Аналогично K(r) носит название матрицы жесткости элемента r.

От типового элемента перейдем к отдельным элементам данной системы. Для элемента 2 с двумя узлами 1, 2 справедливы все зависимости (4)-(7), где следует положить r = 2, i = 1, j = 2. Поскольку точки a, b неподвижны, то состояние элемента 1 определяется перемещением узла 1, а элемента 3 - перемещением узла 2. На основании (4) будем иметь для элементов 1 и 3 соотношения

В (8) координата x для элемента 1 равна нулю в точке a, а для элемента 3 - в узле 2. Зависимости (5) для элементов 2 и 3 примут соответственно вид

Сравнивая (5) или (7) для элемента 2 с (9) для элементов 1 и 3, можно заключить, что вместо матрицы жесткости для двуузлового элемента 2 фигурируют коэффициенты жесткости для одноузловых элементов 1 и 3. Теперь все известно о каждом отдельном элементе системы. Следующим шагом является соединение элементов в узлах на основе условий

где u1 , u2 - перемещения узлов системы 1, 2. Отсюда следует, что состояние соединенных элементов или системы в целом определяется двумя узловыми перемещениями u1 , u2 и рассматриваемая задача является дискретной.

Для всей системы можно записать соотношения типа (5) относительно суммарных для смежных элементов внутренних сил в узлах 1, 2. Обозначим их f1 , f2 . Очевидно, как и в случае типового элемента, они должны линейно зависеть от u1 , u2 :

f1 = f11u1 + f12u2 ,

f2 = f21u1 + f22u2 .

Введем столбцы f, u и матрицу жесткости всей системы K по формулам

Тогда матричное соотношение типа (7) для всей системы будет

f = Ku.

Здесь flt (l = 1, 2; t = 1, 2) есть суммарная для смежных элементов в узле l внутренняя сила fl , возникающая от единичного перемещения узла t. При этом перемещение другого узла равно нулю. Эти суммарные силы определяются через узловые силы в смежных элементах

На рис. 3, б показаны силы при u1 = 1, u2 = 0. При этом элемент 1 растянут, а элемент 2 сжат. В (13) поскольку узел 2 не принадлежит элементу 1, а узел 1 не принадлежит элементу 3. Из (13) следует, что матрица жесткости системы строится на основе коэффициентов жесткости для отдельных элементов. Алгоритмически выполнить это можно по-разному [5]. Например, можно для всех элементов строить матрицы жесткости одинаковой размерности равной размерности матрицы K, основываясь на столбце u перемещений всех узлов системы. Это возможно, поскольку если по крайней мере один из узлов l или t не принадлежит элементу r. В данном примере будем иметь

и, согласно (13), K = K(1) + K(2) + K(3). Из условия равновесия элемента r следует при где N (r) > 0 при растяжении и N (r) < 0 при сжатии. В результате на основании (4) и (8) будем иметь

а матрица K примет вид

Из условия равновесия узлов 1, 2 следует f1 = P1 , f2 = P2 или для столбцов f = P, где P - столбец из P1 , P2 . Подставляя сюда вместо f его выражение согласно (12), окончательно получим систему алгебраических уравнений относительно u1 , u2

После определения u1 , u2 в результате решения (15) находятся u(r)(x), N (r) во всех элементах системы при помощи (4) и (8).

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

СИСТЕМА ЛИНЕЙНЫХ

АЛГЕБРАИЧЕСКИХ УРАВНЕНИЙ

Метод конечных элементов сводит решение линейной задачи к решению системы линейных алгебраических уравнений

f11u1 + f12u2 + _ + f1nun = P1 ,

.............................................

fn1u1 + fn2u2 + _ + fnnun = Pn .

Здесь ui (i = 1, 2, _, n) - неизвестные, Pi (i = 1, 2, _, n) - заданные свободные члены, fi j (i, j = 1, 2, _, n) - коэффициенты при неизвестных. По аналогии с (11) коэффициенты fi j образуют квадратную матрицу, состоящую из n строк и n столбцов:

Если обозначить столбец неизвестных u, а столбец свободных членов P, то (16) принимает матричную форму (15). Система алгебраических уравнений должна быть невырожденной, то есть иметь единственное решение. Казалось бы, дальнейшее ясно. Можно воспользоваться для решения (16), например, методом исключения Гаусса. Однако при применении приближенных методов обычно приходится иметь дело с системами большого порядка n, и матрица, вообще говоря, может иметь такую структуру, которая затрудняет получение решения. При этом на точности результата в той или иной степени сказываются неизбежные в процессе вычислений ошибки округления. Одним из важных достоинств метода конечных элементов является то, что он обычно приводит к таким системам алгебраических уравнений, матрицы K которых позволяют эффективно строить решение.