Учет граничных условий

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

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

Некоторые способы закрепления граничных узлов для плоской модели представлены на рис.15 . Запрещенные степени узловых перемещений зачеркнуты.

Наиболее просто граничные условия могут быть заданы в виде комбинации 0 и 1, при этом 0 означает отсутствие ограничений для заданной степени свободы, а 1 – наличие ограничений. Для узла А (рис. 15) граничные условия должны быть заданы в виде 01, для узла С –11.

Рис. 15

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

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

Алгоритмы построения сеток для решения задач механики деформируемых твердых тел

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

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

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

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

Нумерация узлов в квадратичных четырехугольниках производится в соответствии с рис. 16. Безразмерные координаты узла 1 всегда равны ==-1.

Для каждого четырехугольника выполняются следующие процедуры.

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

2.Производится проверка граничных узлов на предмет отсутствия присвоения им номера. Если такие узлы есть, то за ними сохраняются номера, присвоенные ранее.

Рис. 16

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

=1, =+1 и двигаясь слева направо (при изменении  от –1 до +1) и сверху вниз (при изменении  от +1 до –1).

Узлы, пронумерованные ранее, пропускаются. Номера граничных узлов сохраняются для последующего учета при формировании соседних зон.

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

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

В качестве примера рассмотрена плоская прямоугольная область с круговым вырезом (рис. 17). Для обеспечения достаточно полного описания границы кругового выреза выделено 12 четырехугольных зон (на рис. 17 обозначены римскими цифрами). Стороны четырехугольных зон, примыкающие к вырезу, разделены на два равных отрезка: это обеспечивается заданием соответствующего числа строк и столбцов узлов в каждой зоне.

В результате такого расположения узлов вдоль периметра выреза получена равномерная сетка треугольных элементов. Исходные данные для рассматриваемой области размещены в файле «DT»:

ELEMENT GENERATION FOR EIGEN PROBLEM

12 60

0.0 0.045 0.090 0.105 0.120 0.135 0.150 0.275

0.400 0.400 0.400 0.400 0.400 0.400 0.400 0.400

0.400 0.275 0.150 0.135 0.120 0.105 0.090 0.045

0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

0.045 0.041 0.045 0.090 0.120 0.150 0.275 0.280

Рис. 17

0.275 0.150 0.120 0.090 0.090 0.084 0.081 0.084

0.090 0.105 0.120 0.135 0.150 0.157 0.161 0.157

0.150 0.135 0.120 0.105

0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

0.0 0.033 0.073 0.085 0.100 0.115 0.127 0.167

0.200 0.200 0.200 0.200 0.200 0.200 0.200 0.200

0.200 0.167 0.127 0.115 0.100 0.085 0.073 0.033

0.073 0.10 0.127 0.167 0.171 0.167 0.127 0.100

0.073 0.033 0.029 0.033 0.073 0.085 0.100 0.115

0.127 0.138 0.140 0.138 0.127 0.115 0.100 0.085

0.073 0.063 0.060 0.063

1 2 0 12 0

2 0 3 1 0

3 0 4 0 2

4 0 5 0 3

5 0 0 6 4

6 5 0 7 0

7 6 0 8 0

8 7 0 0 9

9 0 8 0 10

10 0 9 0 11

11 12 10 0 0

12 1 0 11 0

1 3 7 31 33 45 46 47 34 29 30

2 5 7 1 2 3 44 45 33 31 32

3 5 3 3 4 5 43 59 60 45 44

4 5 3 5 6 7 42 57 58 59 43

5 5 13 7 8 9 10 11 41 57 42

6 3 13 57 41 11 12 13 40 55 56

7 3 13 55 40 13 14 15 39 53 54

8 5 13 53 39 15 16 17 18 19 38

9 5 3 51 52 53 38 19 20 21 37

10 5 3 49 50 51 37 21 22 23 36

11 5 7 27 35 49 36 23 24 25 26

12 3 7 29 34 47 48 49 35 27 28

В первой строке обозначена решаемая проблема, в следующей указано число четырехугольных зон и количество граничных точек. Следующие 16 строк заполнены координатами x и y граничных точек. За координатами введены 12 строк, соответствующих числу четырехугольных зон, в которых указаны номера смежных зон. Данные каждой зоны считываются в следующем порядке: номер зоны, число строк и число столбцов, глобальные номера узлов, определяющие четырехугольник.

Вычисленные значения координат узлов и массив связей треугольных элементов сохраняются для проведения дискретного моделирования методом конечных элементов. Сетка конечных элементов для рассматриваемого примера представлена на рис. 18. Число узлов построенной сетки равно 290, число треугольных конечных элементов – 496.

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

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

Рис. 18

Тип элемента вводится перед началом вычислений: «1» – для треугольного элемента с узлами в вершинах, «2» – для треугольного элемента с 6 узлами; если ввести тип «0», то происходит переход к подпрограмме вывода полученных результантов. Затем задается число участков, на которые требуется разделить каждую сторону треугольной зоны. После этого осуществляется ввод номеров узлов и их координат, при этом номера угловых узлов надо подобрать таким образом, чтобы они соответствовали введенному числу деления сторон на элементы. Если ввести отрицательный номер узла, то последующий ввод данных прекратится и произойдет переход к подпрограмме генерации сетки. После ввода номера последнего узла и его координат необходимо ввести число «0», при этом после генерации сетки осуществляется вывод результатов вычислений. При введении новых узлов используется общая система нумерации для всех треугольных зон. Вычисления координат узлов, расположенных вдоль сторон, осуществляется при помощи системы L – координат /3/.

Для угловых узлов треугольных зон вводятся функции формы:

N₁(2L₁-1)L₁; N₂(2L₂-1)L₂; N₃(2L₃-1)L₃ .

Для узлов, расположенных на серединах сторон, функции формы задаются в виде:

N4=4L1L2 ; N5=4L2L3 ; N6=4L3L1 .

Координаты узлов в глобальной системе отсчета вычисляются по формулам:

XNiXi ; YNiYi ; ZNiZi ( i=1,6) ,

где X_i, Y_i, Z_i – координаты угловых узлов, заданные на начальном этапе вычислений. При измельчении моделируемой области, состоящей из нескольких треугольных зон, проверяется совпадение координат узлов, принадлежащих разным элементам. Узлу при последовательной нумерации присваивается номер, имеющий меньшее значение для двух смежных элементов.

В качестве примера, на рис. 19 изображена плоская область с нанесенной исходной сеткой триангуляции; сетка, измельченная с использованием описанного алгоритма, показана на рис. 20.

Третья схема построения сетки конечных элементов для тонкостенных конструкций производится на основе массивов координат узловых точек, полученных при помощи генератора узлов поверхности [6].

Рис. 19

Рис. 20

Использована локальная система координат ,, связанная с исследуемой областью (рис. 21).

Рис. 21 Локальная и глобальная системы координат

Локальные координаты  и  изменяются от 0 до 1 вдоль каждой границы. Так, для границы 1 (рис. 1) =0 и =1 для границы 2.

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

Преобразование локальных координат в глобальную систему может быть определено, например, следующим образом [6]:

где определяются как граничные функции для стороны n.

Для Y и Z могут быть составлены аналогичные соотношения.

С помощью этих преобразований могут быть вычислены глобальные координаты внутренних узлов, если известны их локальные координаты и граничные функции. Локальная координата  узлов на границе 1 вычисляется по формуле:

где

Глобальные координаты, вычисленные в соответствии с описанным алгоритмом, используются для автоматизированного формирования сетки конечных элементов.

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

Сетка конечных элементов рассматриваемой оболочки представлена на рис. 22.

Рис. 22 Сетка конечных элементов

Постановка задачи предусматривает исследование собственных колебаний линейно-упругой тонкостенной оболочки, защемленной вдоль одной из сторон. Для решения неполной проблемы собственных значений использован метод итераций в подпространстве [7].

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

Некоторые из низших частот и форм собственных колебаний лопатки представлены на рис. 23–28.

Рис. 23 Низшая частота ₁ = 526 с^-1 собственных колебаний лопатки

Рис. 24 Низшая частота ₂ = 1469 с^-1 собственных колебаний лопатки

Рис. 25 Низшая частота ₃ = 2223 с^-1 собственных колебаний лопатки

Рис. 26 Низшая частота ₄ = 2386 с^-1 собственных колебаний лопатки

Рис. 27 Низшая частота ₅ = 4017 с^-1 собственных колебаний лопатки

Рис. 28 Низшая частота ₆ = 5461 с^-1 собственных колебаний лопатки

В приложении 1 представлены тексты программ триангуляции плоских областей по схеме 1. В приложении 2 представлен текст программы генерации сетки по схеме 2. (В отладке программ принимали участие студенты РД981 Житенев С.В., Гуртовой А.А., НТ991 Калядин О.В)