Формулировка метода конечных элементов

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

Рассмотрим трехмерный объект произвольной формы, находящийся в равновесном состоянии под воздействием некоторой нагрузки (рис. 8.4). Силы трения, действующие на поверхность, обозначим f ^s, массовые силы — f ^B, а сосредоточенные внешние силы — F ⁱ. В общем случае эти силы раскладываются на компоненты, параллельные осям координат:

Обозначим смещения произвольной точки объекта (X, Y, Z) по сравнению с конфигурацией в отсутствие нагрузки символом U. Тогда

где индекс ^Т означает транспонирование. Смещения U приведут к возникновению деформаций

и соответствующих напряжений

Наша задача состоит в том, чтобы рассчитать U, ε, τ в точке (X, Y, Z) по заданным внешним силам. Возможно, вы знакомы со следующим подходом к этой задаче: основные дифференциальные уравнения равновесия записываются путем наложения условия равновесия на элементы объекта, после чего эти уравнения решаются с учетом граничных условий и условий совместности.

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

Левая часть уравнения (8.5) описывает виртуальную внутреннюю работу, выполняемую реальными напряжениями на виртуальных деформациях, вызванных виртуальными смещениями U. В этом выражении

Слагаемые в правой части выражения (8.5) описывают внешнюю работу, выполняемую реальными силами f ^В, f ^s и F ⁱ на виртуальных перемещениях U, где

Верхний индекс S у вектора U означает виртуальное смещение на поверхности, а верхний индекс i — смещение в точке приложения сосредоточенных сил F ^i. Уравнение (8.5) включает также требования на совместимость и конститутив- ность непрерывных функции смещений, которые удовлетворяют граничным условиям. Напряжения вычисляются через деформации по соответствующим материальным уравнениям. Поэтому принцип виртуальных перемещений включает все требования, которым должно удовлетворять решение задачи строительной механики.

Посмотрим теперь, как из уравнения (8.5) получаются уравнения метода конечных элементов. Начнем с аппроксимации объекта, изображенного на рис. 8.4, сеткой конечных элементов. Элементы соединяются друг с другом в узловых точках, которые находятся на их границах. Смещение в любой точке с координатами (х, у, z) в локальной системе координат элемента считается функцией смещений в узловых точках¹. То есть для элемента т высказывается предположение, что

где Н^(т) — интерполяционная матрица смещений, a U — вектор смещений на всех узлах. Если общее количество узлов равно N, вектор U запишется следующим образом:

Это выражение можно переписать так:

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

Хотя в уравнении (8.10) перечисляются смещения всех узлов, а следовательно, эти смещения входят и в выражение (8.8), для каждого конкретного элемента смещения внутри него определяются только смещениями в его собственных узлах. В уравнение же (8.8) все узлы вошли потому, что это облегчает процесс объединения матриц отдельных элементов в матрицу структуры в целом, как будет показано ниже.

Уравнение (8.8) позволяет вычислить деформации:

Строки матрицы деформаций-смещений В^(т) из уравнения (8.11) получаются дифференцированием и объединением строк матрицы Н^(т). Производные матриц Н^(т) и В^(т) рассматриваются в примере 8.1.

Теперь мы можем записать и выражения для напряжений внутри каждого элемента:

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

Перед тем как подставить выражения (8.8), (8.11) и (8.12) в формулу принципа виртуальных смещений, мы перепишем уравнение (8.5) в виде суммы интегралов по объемам и поверхностям отдельных элементов:

где т изменяется от 1 до полного количества элементов в системе.

Подставляя выражения (8.8), (8.11) и (8.12) в (8.13), будем предполагать, что виртуальные смещения в элементе ū^(т) связаны с виртуальными узловыми смещениями Ū той же матрицей Н^(т) из (8.8). Эта подстановка даст следующее выражение:

где поверхностные интерполяционные матрицы смещений Н^s(т) получаются из объемных интерполяционных матриц смещений Н^(т) подстановкой координат поверхности элемента. F — вектор внешних сосредоточенных сил, действующих на узловые точки. Компоненты вектора F соответствуют компонентам смещений вектора Ū. Обратите внимание, что в уравнении (8.14) вектор Ū вынесен за знак суммирования, потому что он не зависит от рассматриваемого элемента.

Чтобы выражение (8.14) было верным для произвольного виртуального смещения (а это и есть условие равновесия), должно выполняться следующее равенство:

Будем отныне обозначать смещения в узлах просто буквой U. Перепишем уравнение (8.15) в виде

где

Обратите внимание, что суммирование интегралов по объемам отдельных элементов в формуле (8.17) выражает тот факт, что матрица жесткости набора элементов как целого получается сложением матриц жесткости элементов К^(т). Аналогичным образом, вектор R_B объемной силы, действующей на все тело, получается суммированием векторов объемных сил, действующих на отдельные элементы. Тем же путем вычисляются и векторы прочих сил.

Выражение (8.16) описывает статическое равновесие. Если приложенные силы изменяются во времени, это выражение применимо к любому конкретному моменту. Однако при быстром приложении нагрузки необходимо учитывать силы инерции. По принципу Даламбера силы инерции отдельных элементов могут быть добавлены к массовым силам. Если предположить, что ускорение в любой точке элемента связано с ускорениями в узловых точках матрицей Н^(т) подобно смещениям, вклад массовых сил в вектор нагрузки R будет выражаться так:

где Ū — ускорения узловых точек, а р^(т) — массовая плотность элемента т.Слагаемое f ^В(т) в выражении (8.23) больше не включает никаких сил инерции.

Подстановка (8.23) вместо (8.19) в (8.15) дает новое уравнение равновесия:

где М — матрица масс, определяемая следующим образом

Обратите внимание, что U и R в уравнении (8.24) являются функциями времени.

Демпфирующие силы могут быть учтены как дополнительный вклад в массовые силы, что позволяет описать эффект демпфирования (затухания). Уравнение (8.23) при этом принимает новый вид

где Ū — вектор скоростей узловых точек, а к^,(т) — демпфирующий коэффициент для элемента т.

Уравнение равновесия приобретает вид

где С — матрица демпфирования, структура которой описывается выражением

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

Приведенный ниже пример иллюстрирует изложенную процедуру вывода системных уравнений.

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

Анализ методом конечных элементов является мощнейшей технологией, позволяющей моделировать распределение напряжений, температур, потоки жидкостей и распространение электромагнитных полей, однако до сих пор нерешенной остается проблема подготовки данных для проведения анализа: выбор геометрии, построение сетки конечных элементов, добавление граничных условий и нагрузок, задание свойств материалов и выбор типа анализа (статический или динамический, линейный или нелинейный, анализ деформаций, напряжений и т. д.). Действия, относящиеся к подготовке данных, обобщенно называют моделированием конечных элементов (finite-element modeling). Выполняются эти действия чаще всего препроцессором, рассчитанным на работу с какой-либо конкретной программой анализа методом конечных элементов (finite-element analysis — FEA).

Работа с препроцессором начинается с выбора геометрии объекта или области задачи. Традиционные системы FEA обладают лишь зачатками функций моделирования конечных элементов, тогда как большинство современных систем либо снабжаются расширенными средствами моделирования, либо позволяют обмениваться данными с системами геометрического моделирования (а иногда предлагают пользователю и то, и другое вместе). Системы, рассчитанные на подготовку геометрической модели в системах автоматизированной подготовки чертежей, либо работают непосредственно с данными CAD. либо преобразуют и импортируют их. Вариант «поверх CAD» (direct on CAD) становится в последнее время все более популярным, поскольку он устраняет преобразования (которые могут повлечь потерю данных) и сокращает длительность цикла «проектирование — анализ — изменение». Более того, использование CAD упрощает моделирование и дает возможность работать с более сложными функциями создания и изменения геометрических форм. Современные гибридные системы моделирования (интегрирующие объемное, поверхностное и каркасное моделирование с параметрическим и объектно-ориентированным подходами) позволяют создать практически любую нужную для анализа геометрию. Большинство систем FEA могут также импортировать геометрические данные либо через промежуточные файлы стандартных форматов (типа IGES²), либо непосредственно из конкретных CAD. Однако использование геометрических моделей, подготовленных в CAD, не всегда оказывается простым делом. Модель, которую конструктор сочтет идеальной, может на самом деле содержать недопустимые в FEA элементы. Особенно это касается построения сеток. Некоторые системы уже предлагают функции проверки импортированных моделей. Более того, даже если построенная в CAD модель свободна от недостатков, она может быть чересчур подробной. Например, такие характерные детали, как фаски, в некоторых случаях вполне могут быть исключены из модели для анализа методом конечных элементов. Подобные решения принимаются конструктором исходя из ожидаемого размера ячеек сетки, а также из интуитивных предположений о важности отдельных участков объекта. Некоторые программы обладают функциями удаления элементов (defeaturing), то есть временного скрытия деталей, нс влияющих на точность анализа. Абстрагирование является основной причиной различий между моделями одного и того же объекта, используемыми проектировщиками и аналитиками. Изменения, предлагаемые одними из них, не могут непосредственно воплощаться в модели других. В настоящее время ведутся исследования возможности автоматического абстрагирования объемных моделей.

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

От сложности сетки зависит размер глобальной матрицы жесткости, численная сложность задачи и объем требуемых вычислительных ресурсов. Точность решения можно повысить увеличением количества ячеек или использованием функций формы более высоких порядков. Конечные элементы должны удовлетворять определенным требованиям. Во-первых, размерность элементов должна совпадать с размерностью области задачи. Для одномерных задач используются одномерные элементы, для двумерных — двумерные, и т. д. Во-вторых, конечные элементы должны поддерживаться выбранной программой FEA. Другими словами, программа должна уметь рассчитывать вклад конкретного элемента в матрицу жесткости. Все элементы, поддерживаемые пакетом анализа, составляют его библиотеку (elementlibrary). Чем больше элементов в библиотеке, тем большее число задач может решать программа. Наиболее типичные конечные элементы, поддерживаемые большинством программ анализа, демонстрирует рис. 8.7. Обратите внимание, что одна и та же ячейка может становиться элементами разных типов в зависимости от количества узлов на ее границах. Наконец, в зонах, где ожидаются резкие изменения неизвестных (напряжения, например, сосредоточиваются в окрестностях отверстий), плотность узлов и ячеек должна быть выше, чем в областях с плавным изменением параметров.

Другой подход к решению проблем формирования сетки предлагает р-версия конечноэлементного анализа. P-версия использует простые сетки, формируемые автоматическими методами, но зато в этой версии может изменяться степень функции формы (также автоматически). Существует достаточно много программ FEA, поддерживающих р-версию анализа, но только две программы были разработаны специально для этой версии: Pro/MECHANICA фирмы РТС и PolyFEM фирмы CADSI. Преимущества этого подхода не ограничиваются простотой сеток. P-версия позволяет задавать конкретные ограничения на точность, а также лучше аппроксимировать геометрические модели из программ CAD. Низкий уровень точности позволяет конструктору быстро получить результаты анализа на предварительном этапе разработки.

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

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

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