CC BY
11
Соболев В. И., Черниговская Т.Н. УДК 621.06
АЛГОРИТМ ФОРМИРОВАНИЯ ГАРМОНИЧЕСКОГО ЭЛЕМЕНТА В МОДЕЛИРОВАНИИ КОЛЕБАНИЙ ТОНКОЙ ПЛАСТИНЫ
Для расчета на стационарные гармонические воздействия стержневых конструкций известен метод динамических податливостей, основанный на аналитическом выражении амплитуд динамических реакций в узловых связях в зависимостях от величин амплитуд узловых перемещений по возможным направлениям. На основе определения таких зависимостей осуществляется формирование уравнений динамического равновесия в соединительных узлах элементов по направлениям возможных перемещений узлов [1].
В работе [2] на основе развития методов динамической податливости разработаны дискретно-континуальные модели, включающие бесконечномерные изгибаемые элементы (балки), материальные точки, твердые тела и пружины. Создан программный комплекс, реализующий метод расчета систем, включающих произвольные комбинации таких элементов на стационарные гармонические воздействия.
Ансамбль балочных гармонических элементов (ГаЭ) формируется путем декомпозиции исходной динамической системы на изгибаемые одномерные конструктивы, ограниченные краевыми узлами. Для таких элементов осуществляется разрешение узловых динамических реакций и колебательных форм в процессе так называемого гармонического сканирования связей при различных допустимых вариантах краевых условий, задаваемых в граничных соединительных узлах элементов. Метод построения балочного ГаЭ в соответствии с принципами динамической податливости основан на аналитическом выражении величин динамических реакций в зависимости от узловых гармонических перемещений по возможным степеням свободы узловых точек элементов конструкции, в которых осу-
ществляется «сшивка» решений ансамбля гармонических элементов.
На основе применения принципа динамической податливости в сочетании с методом гармонического сканирования связей возможно также построение бесконечномерных изгибаемых гармонических элементов в виде тонких пластин, пригодных для включения в дискретно- континуальные динамические модели. Общая методика данного подхода была приведена в работе [5] , в которой было показано, что в отличие от балочных гармонических элементов [2], динамические реакции в узлах плоских элементов не удаётся определить в аналитическом виде. Для этих целей приходиться использовать приёмы аппроксимации, схожие с конечноэлементными. Однако, в отличие от классических конечноэлемен-тных моделей, при таком подходе удается избежать процедур дискретизации инерционных параметров пластины.
Для построения ГаЭ прямоугольной пластины рассмотрим вариант параметрической сшивки решения в узлах, расположенных в вершинах углов некоторой конечной прямоугольной области тонкой пластины. Используем принцип сканирования связей [2]. Для этого наложим на возможные перемещения узлов связи, обеспечивающие кинематическую узловую неподвижность изгибаемого пластинчатого элемента (рис. 1).
В каждом узле необходимы три связи; одна линейная — по направлению оси Z и две угловых — в плоскостях ZoX и ZoY. Таким образом, степень кинематической подвижности элемента, которая определяет размерность параметрического пространства функции, интерполирующей поверхность перемещений д( х, у), равна двенадцати [3,4].
ИРКУТСКИМ государственный университет путей сообщения
Рис. 1. Расчетная схема тонкостенной пластины. 1,4,7,10 - номера линейных связей, 2,5,8,11 - номера угловых связей в плоскостях параллельных ZoX, 3,6,9,12 - номера угловых связей в плоскостях параллельных ZoY.
Изложим алгоритм построения гармонического элемента прямоугольных пластин с жестким закреплением. Для этого сформируем интерполирующие полиномы для базиса граничных условий, определенного возможными вариантами поочередных единичных перемещений узловых связей.
Для варианта возбуждения первой связи имеем:
д( х, у) «Д( х, у )• Д,
где
Д х, у) = (1, х, у, х2, ху, у2, х3,
2 2 3 3 2 Ч
х у, ху , у , х у, ху ). А1 = (а1, а2, аз, а4, а^ а^,
\Т
(1) (2)
где индекс вектора коэффициентов определяется номером перемещаемой связи.
Перемещая первую связь на единичную величину и фиксируя все остальные в нулевом положении, определяем элементы вектора А1 путем решения системы уравнений:
V- А1 = (1,0,0,0,0,0,0,0,0,0,0,0)
(3)
где
где V- матрица узловых условий поля перемещений [5].
Здесь Ь, с - размеры прямоугольного элемента по осям х и у соответственно. Вычисляя значения функции Д(х, у) и её производных в узловых точках, имеем развернутую запись матрицы V:
V =
Уц V12
V 21 V 22
V» =
V12 =
V 21 =
"1 0 0 0 0 0 "
0 -1 0 0 0 0
0 0 -1 0 0 0
1 0 с 0 0 с2
0 -1 0 0 -с 0
0 0 -1 0 0 -2с
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 0 0 0
0 0 0 с3 0 0
0 0 - -с2 0 0 -с3
0 0 0 -3с 20 0
1 Ь с Ь2 Ьс с2
0 -1 0 -2Ь -с 0
0 0 -1 0 -Ь -2с
1 Ь 0 Ь2 0 0
0 -1 0 -2Ь 0 0
0 0 -1 0 -Ь 0
(4)
V 22 =
Ь3 Ь2 с Ьс2 с3 Ь 3с Ьс3
3Ь2 -2Ьс -с2 0 -3Ь2 с -с3
0 -Ь2 -2Ьс -3с2 -Ь3 -3Ьс2
Ь3 0 0 0 0 0
3Ь2 0 0 0 0 0
0 -Ь2 0 0 -Ь3 0
Выполняя операцию обращения матрицы V, имеем матрицу:
где
V» =
V-1 =
Уп У12 V 21 V 22
1 0 0 0 0 0
0 -1 0 0 0 0
0 0 -1 0 0 0
3 2
0 0 0 0
Ъ2 Ь
1 1 1 1 1 0
Ьс с Ь Ьс с
3 0 2 3 0 1
с2 с с2" с
(5)
^^ 22 =
-2 -1
0 0 0 Ь7 Ь7 0
3 1 0 -3 -1 0
Ь2 с Ьс Ь2 с Ьс
3 0 1 -3 0 2
Ьс2 Ьс Ьс2 Ь
0 0 0 0 0 0
-2 -1 0 2 1 0
Ь Зс Ь2 с Ь Зс Ь2 с
-2 0 -1 2 0 -1
Ьс3 Ьс2 Ьс3 Ьс2
Столбцы матрицы (5) являются коэффициентами полиномов д(х, у), интерполирующих функции перемещений поверхности пластины при соответствующих вариантах задания граничных условий (вариантах перемещения узловых связей).
Например, при перемещении первой связи вектор коэффициентов полинома имеет вид:
3 3 2 -2 -2
А =
-3 _1 -3 2 10 0 _3 -3
Ь2 Ьс с2 Ь3 Ь2с Ьс2 с3 Ь3с Ьс3
V =
у 21
V,, =
Ь2 с
Ьс2
2 Ь 3с 2
Ьс3
0
0
0
0 1
Ьс 3
2
с
-— 0
2
Ьс 0
0 1
Ь2 с
0 0 0 0 _ 1 с
0 -
2 1
— -— 0 Ь3 Ь2
0 0 0
0
_1_
Ьс 0
Ь2 с
Ьс 1
Ьс2
0 0" 0 0 0 0 0 0 0 0
0 0
0 2
Ьс 0
— 0
с 0
Ьс2
с 2
Ь Зс
Ьс3
Ь2 с
Ьс 1
с 0
Ьс2
Соответствующий интерполирующий полином представлен в виде:
_3 _ 1 _з 2
д(х<У) =1 + Т2х2 +— ху + — у2 + -3X3 +
Ь2 Ьс с2 Ь3 (6)
3 2 3 2 2 2 2 3 2 3 х Ух У+"ТхУ + тт~х у+ хУ .
Ь2 с Ьс2 с Ь с Ьс
Имея равенство (6), можно для сформированных граничных условий определить аналитические выражения реакций в узловых связях. Скомпонуем для этого вектор деформаций изгиба пластины в виде:
^Х ^ ^Х х X У X хУ ) ,
(7)
где X:
д д
~, Х У
д д Х =-2 д д
2 ' Х ху ^
дх ду дх ду
величины кривизны изгиба и скручивания пластины.
В данном случае вектор деформаций изгиба примет вид:
1
1
иркутский государственный университет путей сообщения
М =
6(Ь - 2х)(с - у ) Ь Зс
б(с - 2у )(Ь - х )
Ьс3
2 (^ _ 6у (с - у )_ 6х(Ь - х) Ьс с2 Ь2
(8)
Ох =--
Eh3 ( 2(с - у ) (с - 2у )( 2ц-1)
2(1 -ц2)
Ь с
Ьс3
Оу = -
Eh3 ((Ь - 2х ) 2(Ь - х)
2(1 -ц2)
Ь с Ьс3
Для определения вектора моментов {М} = (Мх Му Му) используем равенство:
{М} = С{х},
(9)
где для изотропного материала пластины
С = В •
1 ц 0 ц с 0 1 -ц
0 0
2
, ц - коэффициент Пуассона,
В =
Eh3
12(1 -ц2)'
Выполняя преобразование (9) с использованием (8), получим
{М} =
Eh3
12(1 -ц2)
6(с - у)(Ь -2х) + 6ц(Ь -х)(с -2 у)
Ь3с Ьс3
6(Ь -х)(с-2 у) + 6ц(с - у)(Ь -2х)
Ьс3 Ь3с
(ц -1)( 6х(Ь - х) + 6у(с - у)-1 Ьс ( Ь2 с2
(10)
Для определения реакций в связях воспользуемся условиями равновесия элементарного участка пластины [5,7]:
О + О
дх ду
+д( х, у) = 0,
дм дМх
ду
дх
дМх дМ
- Ох = 0,
ху
дх
ду
- Оу = 0.
(11) (12) (13)
Условие равновесия (10) выполняется на расчетной области почти всюду за исключением узловых точек, в которых величины Ох, Оу имеют разрыв в силу того, что в этих точках возникают реакции в линейных связях. Для узловых точек справедливы условия равновесия в виде:
Ох + Оу + Г{ 1 = 0, (14)
где / - номер линейной связи, в которой формируется реакция от единичного перемещения первой связи,
Eh
12(1 -ц2)
6(Ь -2х) + 12(с - у) 12(Ь -х) + 6(с -2 у)(2ц -1)4
(15)
Ь3 с
Ьс3
Величины линейных реакций г11, г41, г71, г101 определяются из равенства (15) при подстановке соответствующих координат узловых точек.
Выражения реакций в угловых связях определяются аналогичным образом из условий равновесия моментов в узловых точках:
Мх + Мху + Гх1Л = 0,
Му + Мху + Гу,., = 0.
(16) (17)
Подстановка найденных выражений компонент вектора {М}в (12) и (13) позволяет определить величины перерезывающих сил
Ох, Оу:
То есть из (16) находим величину реакций в угловых связях 2,5,7,11 в плоскости ZoX, а из (17) —в связях 3,6,9,12 плоскости ZoY.
Полученные величины г, 1, гх1 1, гу, 1 формируют первый столбец матрицы N единичных реакций, определенных при помощи параметров вектора А1 , выраженных в свою очередь через: геометрические параметры ^ Ь, с; механические параметры т, Е, ц пластины.
Сформированные компоненты матрицы N не учитывают влияния распределенной инерционной нагрузки ю2 mhg(х, у). Для ее учета воспользуемся теоремой о взаимности работ [6], согласно которой работа по преодолению внешних сил при перемещении связи
равна работе, совершенной поверхностной нагрузкой при прогибах пластины.
Исходя из того, что при перемещении используются связи с номером г, получим:
„т Ь с 2
2
1=! ах 1
ю mhg (х, у) 2
ау. (18)
Таким образом, окончательно диагональный элемент формируется в виде суммы
Г. = Гй + Г™ .
(19)
Формирование второго столбца матрицы N осуществляется аналогично при единичном перемещении второй связи. Система уравнений (3) при этом имеет вид:
V- А 2 = (0,1,0,0,0,0,0,0,0,0,0,0)7
(20)
Последовательное решение систем уравнений вида (3) с бегущей единицей в правой части равносильно обращению матрицы узловых условий V, поэтому справедливо:
V-1 = А, где А = (А1, А2,... А12 )7.
(21)
Выполнение операций (7) — (19) с последовательной подстановкой столбцов матрицы А в выражение (1) формирует матрицу динамических жесткостей (амплитуд единичных динамических реакций ) N что позволяет аппроксимировать амплитудные состояния стационарных колебаний изгиба посредством узловых соотношений вида:
N •и = Р, (18)
где и - вектор амплитуд обобщенных узловых перемещений, а Р - вектор амплитуд узловых сил.
Матрица амплитуд единичных динамических реакций N стационарных моногармо-
нических колебаний прямоугольного изгибаемого элемента с жесткими закреплениями (матрица гармонического прямоугольного элемента пластины) получена в аналитическом виде и представлена в таблице 1.
Таким образом, предложенный алгоритм реализации гармонического элемента обеспечивает устранение процедур дискретизации инерционных параметров пластины и связанных с ними погрешностей расчета. В алгоритме существование погрешностей обусловлено лишь наличием процедур аппроксимации жесткостных свойств пластины, также неизбежных в дискретных методах. Параметры модели остаются распределенными, а исходное уравнение в частных производных не заменяется системой обыкновенных дифференциальных уравнений динамики, описывающих колебания материальных точек дискре-тизированных в узлах ансамбля элементов.
БИБЛИОГРАФИЯ
1. Колоушек В. Динамика строительных конструкций. - М.: Издательство литературы по строительству, 1965.- 632 с.
2. Соболев В. И. Дискретно-континуальные динамические системы и виброизоляция промышленных грохотов. — Иркутск: Изд. ИрГТУ, 2002.-202 с.
3. Галлагер Р. Метод конечных элементов. Основы. - М.: Мир, 1984.-428 с.
4. Бате К., Вильсон Е. Численные методы анализа и метод конечного элемента. - М: Стройиздат, 1982. — 447 с.
5. Соболев В. И., Черниговская Т.Н. Построение прямоугольного гармонического элемента для моделирования колебаний тонкой пластины // Современные технологии. Системный анализ. Моделирование. Вып. №4(16). ИрГУПС. Иркутск.2007. С.28-32.
00
иркутский государственный университет путей сообщения
<
ю (в н
н
и
(в
<
С
(в н
< я
л
<
о
и
^
О «
а с
о
и
О *
и
V
а
«
и (в и
В а
н
«
* <
О <
о а с
ю «
н
