Точечно-узловая динамическая матрица жесткости для прямоугольной пластины Текст научной статьи по специальности «Физика»

Точечно-узловая динамическая матрица жесткости для прямоугольной

пластины

С.О. Папков, Ю.И. Папкова, В.А. Пасечник Севастопольский государственный университет, Севастополь

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

Ключевые слова: пластина, колебания, собственные частоты, метод динамической жесткости, бесконечная система линейных уравнений.

Введение

Как известно, пластина является важнейшим элементом в структурной механике и инженерных приложениях, на сегодняшний день имеются сотни работ, посвященных данной тематике, например [1, 2]. Метод спектральной динамической жесткости (Dynamic Stiffness Method или сокращенно DSM) изначально был предложен Wittrick и Williams [3] для анализа колебаний и устойчивости ансамблей прямоугольных пластин в случае, когда две противоположные стороны каждой из пластин являются шарнирно опертыми. Заметим, что только для этого случая граничных условий имеется точное решение задачи о поперечных колебаниях прямоугольной пластины в форме тригонометрического ряда. Данное решение позволяет построить точные аналитические соотношения между векторами нагрузок и векторами смещений на границах пластины, которые определяют динамическую

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

В последующие годы круг применения DSM был существенно расширен. В статьях Banerjee и Boskolo [4 - 5] данный метод развит для анализа продольных колебаний и уточненной теории поперечных колебаний прямоугольных пластин. В [6] была предложена эффективная модификация данного метода для анализа колебаний тонкостенных конструкций и многослойных пластин. В работах [7, 8] представлено развитие DSM для элементов в форме прямоугольных пластин. Хотелось бы заметить, что преодолеть существенное ограничение метода для прямоугольных пластин -шарнирное опирание противоположных сторон, можно лишь на основе новых аналитических решений теории упругости. В [7, 8] для этих целей используются решения для пластин на основе метода суперпозиции, развитого Gorman [9]. При использовании данного подхода ряды, представляющие решения, усекаются, следовательно, с возрастанием частоты колебаний приходится удерживать большее число членов ряда. В статье [10] используется декомпозиция задачи на четыре блока в зависимости от типа симметрии от координатных осей, что позволяет построить спектральную матрицу жесткости пластины для каждого случая симметрии по отдельности. Более того, в упомянутой статье для явной записи спектральной матрицы жесткости используется точное обращение бесконечных матричных блоков на основе асимптотической теории бесконечных систем линейных уравнений, то есть учитываются бесконечные остатки рядов, представляющих решение.

Заметим, что развитие DSM для точечно заданных граничных характеристик открывает возможности для построения гибридных

вычислительных схем и более широкого использования DSM в моделировании механических структур.

Математическая модель

Рассмотрим колебания прямоугольной ортотропной пластины {(х,y)e[-a;a]х[-b;b]} толщины h. Как известно, в рамках теории тонких пластин, уравнение колебаний пластины относительно прогиба w(х, y, t) = W(х, y)eimt имеет вид

d4W d4W d4W D ++D2 - D1Q4W = 0 (1)

1 dx4 3 dx2dy2 2 ay41 w

где Q = a 4ю2р h / D1 -безразмерный частотный параметр, р - плотность материала, ю - круговая частота, D1, D2, D3 -упругие постоянные.

Пусть входящие граничные смещения и нормальные углы поворота на каждой из сторон прямоугольника Wakj(y), W¡k(x), y), $>(x) также как силовые характеристики Mka(y), Mbk(x),Vk(y),Vk(x) заданы точечно. В частности, предположим, что в зависимости от типа симметрии (индекс к =0 соответствует четной функции координаты*, а к = 1 - нечетной) функции координаты х заданы на [0; a] в (N - 1+ к) равноотстоящих узлах с шагом

. , _ a 2a , ,

a/N. А именно, при к =0 в точкахx0 = 0,xJ = —,x2 =—,...,xN = a и при к = 1 в

a 2a

точкахxJ = —,x2 = —,...,xN = a. Аналогично, пусть функции координаты y также заданы точечно на [0;b] в (N - 1+ j) равноотстоящих узлах

y = (l = J,J+ J,. .,N).

Построим по заданным точечным значениям интерполяционные тригонометрические полиномы координаты х по системе функций

{Tk (<^nkx)Yn^l k и для функций координаты y по системе функций \rj($njy)Ynll~},

где обозначено:

II

/ • л

Tj.(z) = cos П-z " ;, ank = - n-1 + 1 , ßnj = n

2

cos z, j = 0 n

r ir\

V 2

к

2

n

n -1 + j

V 2 j

• . л •> nk

[sin z, j = 1 a

Так, приближенное представление прогиба Wk( y), будет иметь вид:

N+1- j

Wak(y) - WN(y) = YWWOLnTTí(ß*y) (2)

m=1

где коэффициенты Wam, N определяются по методу наименьших квадратов и представляют наилучшее приближение точечно-заданной функции {Wk( y)} в пространстве интегрируемых с квадратом функций L2 по заданной системе. Фактически можно говорить об отображении последовательности действительных чисел {Wak(y)} на действительную последовательность

Wam,N}. Подобное отображение в математической литературе принято называть дискретным преобразованием Фурье. Для введенной системы тригонометрических функций T. (ßn]y)fnll~1 дискретное преобразование

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

Для отыскания по методу наименьших квадратов коэффициентов WjN в формуле (2) решаем задачи минимизации:

£ W (Уп) - WON(Уп $ ^ min (3)

п=j

Далее, следуя стандартной технике решения подобных задач, получаем из (3) систему уравнений для отыскания коэффициентов Фурье:

N f N+1-j \

£ £Wakm,NTj j ) - W1 (Уп) T j ) = 0 ( l = 1,2,..., N + j -1) (4)

n=j V m=1 J

N+1- j N

или £MimWlmN = £ T} (ßy W (Уп ) (5)

m=1 n=j

и

N

где М{т = Т Т] ] )Т] ф1ту„)

(6)

п= ]

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

М = Мгт Хт=\.] допускает точное суммирование в представлении своих

элементов и позволяет построить матрицу обратную к ней с = М-1. Зная элементы <з]т матрицы о, можно явно записать, что:

N Г N+1" ]

^ = Т Ъ<Т]](Уп)

п=] V I=1 )

(7)

Так как интерполяционный тригонометрический полином (у) в

узлах принимает заданные значения (уп) = ^ (уп), то обозначая, как:

N+1-]

?тп = ($тУп)

I=1

(8)

получаем формулы прямого и обратного перехода для данного типа преобразования:

N+1-]

Жк (Уп ) = ТТ Фт]Уп

т=1

N

Т ^тг^а (Уп ) п=]

(9)

(10)

Аналогичным образом могут быть построены интерполяционные тригонометрические полиномы и для остальных граничных кинематических и силовых характеристик. Обозначим интерполяционные многочлены, как:

N+1-] N+1^

(У) = ТЖ^Т] (КУ)> Ж] (х) = (а ткХ)

т=1

т=1 N+1-к

N+1"]

Ф] (У) = Т]Т] ($т]У\ Ф] (х) = Т^Тк (а ткХ)

т=1

N +\-]

т =1

N +\-к

(11)

м] (У) = Т Мкт,Т ®]У\ М] (X) = Т М]^Тк (а ткХ)

т=1

N +1- ]

V] (У) = Туа]л- ®т,У)> < (х) = ТУь^Тк (а ткХ)

т=1 N +1-к

т=1

т=1

(12)

>

>

и

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

N

N

К

ат^

V (у ),

/ 1 тп а п /'

К

к

Ьт^

V рткк (Хп)

п=}

N

п=Г

N

Фa'm,N = V (уп X фЬт^ = V ^тп^Ь (Хп )

п=Г

п=Г

(13)

и для силовых характеристик:

М*ат* =£ КМ (Уп ),

п=]

N

ук\г = V —ГУк (у ),

ат,м ^^ тп а п Л

N

МЬтN =1 КМ (Хп )

п=Г

N

УЬтN V —тг^Ь (Хп )

п=Г

п=Г

(14)

Из теории рядов Фурье известно, что коэффициенты интерполяционных тригонометрических полиномов (11), (12) будут стремиться к точным значениям коэффициентов Фурье этих функций при увеличении N. Тогда, подставляя (13), (14) в соотношения спектральной динамической жесткости для прямоугольника [10], получаем 2(2^ 2 - к-/) уравнений, связывающих значения точечно заданных граничных характеристик пластины. Данные уравнения имеют вид (к, / = 0, 1):

N 2 — Й Й N+1—к (л В2 + Л а2 )—к N

2 Й-"Й1т у (—\)1+т л 2Втг + П12^1к/Г1п п=к ълх

В2

V2 Й0Й1т V (-1)1+^П22рт] 2' Л12Гй-2\ ФЬ] (х„) +хлт,^тпФаг ^)+

(ак + Ят] )(ак + Чтг) пг]

2 \Т7к 1п

т: фЬ (хп,

1=1

N 2 — Й Й N+1-к ^ 2 Й г 0Й 1т ^

V т V

(—1)1+т-1пк

1 1=1 « + ч1г )(«2к + )

N N

^ V* (Х„ ) 'V —тп ^а (Уп ) = V —Г К* (Уя )

(15)

2 — 8,.п8

Г 0 1т

N+1—к

п=к ьлх

V (—1)

((л2 — 4П626 —лл Вг а* — ДЛО4—

1+т " 3

п=Г N

(« к + 4 )(«/1 + ^)

^ ф* (Хп ) + 1А3т—тГпФа (Уп ) —

КЬ (Хп ) + V т тп а п

« 2 — 8 г о8^^ (л^ + )—1п V ^ V (—1)

/= (^/к 1 Ятг к 1 ЯтГ

(16)

п=к ЬП;

/2 2 \/ 2 — 2 \ Ь \ п) (< + ^Г )(«£ + Ят, )

У? (Х ) + V А4 —г Ук (у ) = — V —г Мк (у )

Ь У п/ / * т т^ а Куп; / > т^ а \ у п /

п=Г

п=Г

(т = 1, 2, ..., N + 1 -/)

кг

>

>

к

п

п=г

1=

N 1 Я Я N+1_j

X £ (_1)

aD-

N

l+m (AaL + №2Мп ф^ (y ) +ХД рк ^ ) +

—2 \ т a У-' n / / i m mn^b \ n'

n=k

2

l =1

(ßj + p2k )(ßj+p2k)

N

N

N 9_Я Я N+1_] (_1l+mF]

+XX ,02 ,(_Л,П22 чVk(yn) +X*mFkmv?(xn) = XFkmnwk(x„)

n=k

(17)

n=]

aD

2

l =1 (ßl] + P2mk Ш] + P2mk )

n=k

N 2 _ö,,nö

n=] aD2

N+1_ j

X (_1)

l=1

N

k0Ö1m ^ i_nl+m ((D3 4D66 D1D2 jfyj0,mk D} DO ф^ y ) +у д? Fk фЪ (x ) _

фа (yn) X i m mnTb ( n)

(ß]+pik )(ß]+Pik)

n=k

N П Я Я N+1_j

j^2 _Ök0Ö1m

n=j aD2

(d a2 + D ß2 n N

X (_1)l+m L 2 ih ^ (Уп ) +x (xn ) = _XFL№bb' (x )

(18)

l=1

(ß2- + pik )(ßt + PL)

(m = 1, 2, ..., N + 1 -k)

Заметим, что учитывая соотношения симметрии:

4Wk (Уп ) = W (a, Уи ) + (_1)kW (_a, Уп ) + (_1) JW (a,_ Уп ) + (_1) k+JW (_a,_yn) W (xn, b) + (_1) ]W (xn ,_b) + (_\)kW (_xn, b) + (_1) k+]W (_xn ,_b)

W (Xn)

4ф1 (Уи) = фy (a, Уи) _ (_1)kф У (_a, Уи) + (_1)] ФУ (a,_Уи) _ (_1)k+] ФУ (_a,_Уи)

4ф* (Хи) = ф X (Хи, b) _ (_1)] ф X (Хи ,_b) + (_1)k ф X (_ Хи, b) _ (_1)k+] ф X (_ Хи ,_b)

4Mk ( Уп ) = Mx (a, Уп ) + (_1) kMx (_a, Уи) + (_1) ]MX (a_yn) + (_1)k+JMx (_a_yn)

4Mb ( хп ) = My (Xn, b) + (_1) JMy (Xn ,_b) + (_1)kMy (_ Xn, b) + (_1)k +My (_ Xn ,_b)

4Vk (Уи) = Vx (a, Уи) _ (_1) kVx (_a, Уи) + (_1) JVx (a,_y„) _ (_1)k+JVx (_a_yn) 4Vbk(Xn) = Vy(Xn,b) _ (_1)Wy (xn,_b) + (_1)kVy (_Xn,b) _ (_1)k +Vy(_Xn,_b)

(19)

уравнения (15) - (18) фактически определяют связь между заданными в узловых точках на границе пластины смещениями и усилиями. Действительно, подставляя (19) в уравнения (15) - (18), и рассматривая четыре случая симметрии (k, j), получаем связь между точечно заданными компонентами векторов смещений и усилий в общем случае. Таким образом, представленная система соотношений (15) - (18) задает динамическую матрицу жесткости (DSM) ортотропной прямоугольной пластины

Численные результаты

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

k

n

и

примера. Первый пример относится к расчету квадратной пластины, защемленной к углам (рис.1) при у=0.3 и с = 0.5, й =0.75. Собственные частоты и формы колебаний для данной пластины не представлены в литературе и здесь приводятся впервые (таблица 1).

Рисунок 1. Пластина с защемленными углами

Таблица №1

Собственные частоты для квадратной пластины, защемленной по углам

при у=0.3, с = 0.5 и й =0.75

1 2 3 4 5 6 7

п 1.5708 2.1878 2.2214 3.1415 3.5124 4.2082 4.3136

Рисунок 2. Составная ортотропная пластина с защемленными краями

В качестве второго примера рассмотрим составную ортотропную пластину (рис.2) с защемленными краями. При этом левая половина пластины имеет параметры Е1 = 60.7 ОРа, О12 = 12 ОРа, у12 = 0.23, у21 = 0.094 (кварцевое стекло), а правая половина сделана из того же ортотропного материала, но повернута на 90°, то есть имеет параметры Е2 = 60.7 ОРа, Оп = 12ОРа, у12=0.094, у21=0.23. В табл. 2 даны первые семь собственных

частот такой пластины в сравнении с собственными частотами однородной пластины с защемленными краями из того же материала.

Таблица №2

Собственные частоты для квадратной составной ортотропной пластины, защемленной по краям

N 1 2 3 4 5 6 7

Составная пластина (рис.2) 2.6823 3.8172 3.8367 4.6408 5.0246 5.0389 5.6394

Однородная защемленная пластина 2.6975 3.5649 4.1015 4.6283 4.6678 5.4446 5.6157

Из данных табл.2 можно увидеть, что некоторые частоты сплошной и составной пластин достаточно близки (например, первая, четвертая, седьмая и др.), но есть и существенные отличия. Так пятая и шестая собственные частоты составной пластины оказываются отличными лишь на 0.3 %, и существенно отличаются от соответствующих собственных частот сплошной пластины. На рис. 3 представлены собственные формы, отвечающие первым собственным частотам составной пластины. Можно увидеть, что если первая и вторая собственные формы мало отличаются от первой собственной формы квадратной сплошной пластины с защемленными краями, то пятая и шестая формы имеют иной качественный характер, наблюдается выраженная асимметричность.

(а) (Ь)

(с) (<!)

Рисунок 2. Первая (а), вторая (Ь), пятая (о) и шестая (ё) собственные формы защемленной составной пластины с параметрами ^ = 60.7 GPa,

^ = 12 GPa, у12 = 0.23, у21 = 0.094 для левой части и Е2 = 60.7 ОРа, ^ = 12 GPa, у12=0.094, у21=0.23 для правой части

Выводы

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

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

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

Литература

3. Дородов П.В. Исследование напряжений на линии сопряжения ступенчатой пластины // Инженерный вестник Дона, 2013, №2. URL: ivdon. ru/ru/magazine/archive/n2y2013/1636

2. Овчинников М.А., Сокол В.А., Соловьева О.Ю., Тарасова Т.А., Грецова Н.В., Клячина Н.В., Лагунов Е.Н. Математическое моделирование изделий из композитных материалов с заданными свойствами // Инженерный вестник Дона, 2023, №1. URL: ivdon.ru/ru/magazine/archive/n1y2023/8131

3. Wittrick W.H., Williams F.W. A general algorithm for computing natural frequencies of elastic structures // Quart J Mech Appl. Math. 1970. 24 (3). pp. 263-284.

4. Boscolo M., Banerjee J.R. Dynamic stiffness formulation for composite Mindlin plates for exact modal analysis of structures. Part I: theory // Comput. Struct.2012. 96 (97). pp. 61-73.

5. Boscolo M., Banerjee J.R. Dynamic stiffness elements and their applications for plates using first order shear deformation theory // Comput. Struct. 2011. 89. pp. 395-410.

6. Pagani A., Boscolo M., Banerjee J.R., Carrera E. Exact dynamic stiffness elements based on one-dimensional higher-order theories for free vibration analysis of solid and thin-walled structures // Journal of Sound and Vibration. 2013. 332. pp.6104-6127.

7. Nefovska-Danilovich M., Petronijevic M. In-plane free vibration and response analysis of isotropic rectangular plates using the dynamic stiffness method // Comput. Struct. 2015. 152. pp. 82-95.

8. Nefovska - Danilovich M., Kolarevich N., Marjanovich M., Petronijevich M. Shear deformable dynamic stiffness elements for a free vibration analysis of composite plate assemblies - Part I: Theory //Comput. Struct. 2017. 159. pp. 728 -744.

9. Gorman D.J. Free vibration analysis of rectangular plates. - New-York: Elsevier, North Holland; 1982, 324 p.

10. Banerjee J.R., Papkov S.O., Liu X., Kennedy D. Dynamic stiffness matrix of a rectangular plate for the general case // Journal of Sound and Vibration. 2015. 342. pp. 177 - 199.

References

1. Dorodov P.V. Inzenernyj vestnik Dona, 2013, №2. URL: ivdon. ru/ru/magazine/archive/n2y2013/1636

2. Ovchinnikov M.A., Sokol V.A., Solovyova O.Yu., Tarasova T.A., Gretsova N.V., Klyachina N.V., Lagunov E.N. Inzenernyj vestnik Dona, 2023, №1. URL: ivdon.ru/ru/magazine/archive/n 1y2023/8131

3. Wittrick W.H., Williams F.W. Quart J Mech Appl. Math. 1970. 24 (3). pp. 263-284.

4. Boscolo M., Banerjee J.R. Comput. Struct. 2012. 96 (97). pp. 61-73.

5. Boscolo M., Banerjee J.R. Comput. Struct. 2011. 89. pp. 395-410.

6. Pagani A., Boscolo M., Banerjee J.R., Carrera E. Journal of Sound and Vibration. 2013. 332. pp.6104 - 6127.

7. Nefovska-Danilovich M., Petronijevic M. Comput. Struct. 2015. 152. pp. 82-95.

8. Nefovska - Danilovich M., Kolarevich N., Marjanovich M., Petronijevich M. Comput. Struct. 2017. 159. pp. 728 - 744.

9. Gorman D.J. Free vibration analysis of rectangular plates. - New-York: Elsevier, North Holland; 1982, 324 p.

10. Banerjee J.R., Papkov S.O., Liu X., Kennedy D. Journal of Sound and Vibration. 2015. 342. pp. 177 - 199.

Дата поступления: 29.12.2023 Дата публикации: 7.02.2024