CC BY
22
Численные методы расчета конструкций
УЧЕТ СМЕЩЕНИЯ КАК ЖЕСТКОГО ЦЕЛОГО ОСЕСИММЕТРИЧНО
НАГРУЖЕННОЙ ОБОЛОЧКИ ВРАЩЕНИЯ НА ОСНОВЕ МКЭ
А.П. КИСЕЛЁВ, канд. техн. наук, доцент Р.З. КИСЕЛЁВА, канд. техн. наук, доцент А.П. НИКОЛАЕВ, д-р техн. наук, профессор
ФГБОУ ВПО Волгоградский государственный аграрный университет, 400002, г. Волгоград, Университетский проспект, 26; vhs2010@yandex.ru
Для исследования напряженно-деформированного состояния оболочки вращения при смещении как жесткого целого предлагается инвариантная аппроксимация полей перемещений в векторной формулировке, реализованная при формировании матрицы жесткости кольцевого объемного конечного элемента, поперечным сечением которого является четырёхугольник с узлами ¡, ], к, I. В качестве узловых неизвестных конечного элемента выбирались перемещения узловых точек и их первые производные. На примере расчета показано, что предложенная инвариантная аппроксимация полей перемещений в векторной формулировке позволяет решение проблемы учета смещения конечного элемента как жесткого целого в общем виде.
КЛЮЧЕВЫЕ СЛОВА: МКЭ, инвариантная аппроксимация полей перемещений, векторная формулировка, осесимметрично нагруженные оболочки, кольцевой объемный конечный элемент, узловые неизвестные, смещение жесткого целого.
Для определения напряженно-деформированного состояния конструкций из тонкостенных элементов эффективно используется метод конечных элементов (МКЭ) [1,2,3,4,5,6]. Корректным расчет тонкостенных инженерных сооружений, созданных, например, на основе аналитических поверхностей [7], возможен только при использовании криволинейных систем координат. В работах [2,3] отмечается, что при использовании таких систем координат появляется проблема учета смещения конечного элемента как твердого тела.
В настоящей работе предлагается решение этой проблемы при получении матрицы жесткости кольцевого конечного элемента на основе инвариантной аппроксимации искомых величин в векторной формулировке.
1. Геометрия осесимметрично нагруженной оболочки вращения
Положение произвольной точки М срединного меридиана осесимметрично нагруженной оболочки вращения в декартовой системе координат xoz определяется радиус-вектором (рис.1)
Я = хX + г (х) к, (1)
где х, j, к - орты декартовой системы координат; г - радиус вращения точки. Векторы локального базиса точки М определяются выражениями
а1 = Я, ж = х
+ r,
к; а = а1 х j =
= (х, 6.Х + г, ) х j = -г, аХ + х, ак, (2) где s - координата, отсчитываемая вдоль дуги отсчётного меридиана; г, ж = г, х х, ж - производная радиуса
' X
Рис. 1. Базисные векторы точек оболочки
вращения по дуге s.
Радиус-вектор произвольной точки М оболочки, расположенной на расстоянии С от ее срединного меридиана, представляется выражением (рис.1):
яс = Я + &. (3)
Базисные векторы точки определяется дифференцированием (3) по координатам 5 и
§1 = = ах,^ +&„; §з = = а . (4)
С учетом (2) можно составить матричные выражения:
{§}=И{"I ЬЫ"1 (5)
2x1 2x2 2x1 2x1 2x2 2x1
где {§}Г = § §з1 ^ Г
Дифференцированием (5) можно записать матричные соотношения:
{§, Ы }=Ы Ы ]{Т }=[ы, Ы ][Ы]1{§}=И{§}; {§, , }=к , Ш=к ][Ы]1 {§}=[-]{§}. (6)
Произвольная точка М оболочки под действие м заданной нагрузки получит перемещение, которое определяется вектором V с компонентами в базисе соответствующей точки:
V = V1 §1 + у§ з. (7)
Производные вектора перемещения (7) по координатам 5 и С с учетом (6) запишутся выражениями:
V ы = §1 + V? ) Ы = ^§1 + vl ы +V, ы § 3 + V? 3, Ы =
= V ,1Ы §1 + V1 (Уп §1 + У12 §3 ) + V ы § 3 + 4^21 §1 + У 22 §3 )= (8)
= §1 ^ +^Уп + ^21 )+ §3 ^Уп + V Ы +У )= У/ § 1 + /13 §3;
V, С = (у1 §1 + у§3 \ с = V \ §1 + ^ §l, с с §3 + 3, с = /,3 §1 + V с §3 .
Векторы локального базиса точки
в деформированном состоянии Мс
определяются выражениями
§Г = ^ Ы *=Д Ы , Ы = §1 + ?, Ы ; §3 = Д с *=Д с , с = §3 + V, с. (9)
щи в точке М
ной среды [8]:
г *
Деформации в точке Mh определяются соотношениями механики сплош-
Saß ={g"aß - gaß V 2 . (10)
Принимая во внимание (8), (5) и (6), можно представить деформации через компоненты вектора перемещения выражениями:
=[gi (1 + $721)+ g3 С 722 ][gid1 + dg з J= d1 (1 + ^) + dt; 722 = v,1 (l + СГ21) +
+ vl [711 (1 + ;721 )+ 722 ; 721 ]+ v[721 (1 + ;721 )+ 722 С 722 ] + V,, С 722; езз = g3 •V, ; = g3 -fe g1 + v, ; gз) = v, ; ;
2е13 =[g1 (1 + ;721)+ gз ; 722g1 + v ; gз)+ gз(d 1 g1 + dgз)= (11)
= vk (1 + ; 721 ) + ; 722v, ; +v712 + V,. +v722 , или в матричном виде {s}=[z]{v}, (12)
4x1 4x2 2x1
где {s}r = {s11 s22 s33 2s13}-вектор-строка деформаций в произвольной точке
оболочки; {vf = {v1 v}- строка перемещений произвольной точки оболочки.
Связь между напряжениями и деформациями определяется соотношениями механики сплошной среды [8]:
CTaß = ЯВДgaß + 2/ga7gßpsyp , (13)
где Л, / - параметры Ламе; I1 (е) = £пg11 + g22£22 + е33g33 + 2е12g13 - первый инвариант тензора деформаций, gaß - контравариантные компоненты метрического тензора; gaß - ковариантные компоненты метрического тензора.
Зависимости (13) можно представить в матричном виде
М=С ]И. (14)
4x1 4x4 4x1
2. Матрица жесткости конечного элемента
Конечный элемент принимается в виде кольцевого объёмного тела вращения с поперечным сечением в форме четырёхугольника с узлами г, 1, к, 1 [9]. Для выполнения численного интегрирования произвольный четырёхугольник отображается на квадрат с локальными координатами £ П, изменяющимися в пределах -1 < £,л < 1.
Координаты внутренней точки конечного элемента выражаются через координаты узлов четырёхугольника с использованием билинейных функций локальных координат £,л соотношениями:
8 = {/ (ЫГ к }; С = {/ (ЫГ Су Ц (15)
где }т = {у',18к 81} С у }Т = СС1 СС }-матрицы-строки узловых координат четырёхугольника. Из соотношений (15) определяются производные глобальных координат 5, £ в локальной системе координат 8,^ , 8,Л ,С,% С,Л и локальных координат С, П в глобальной системе координат 8, г],с , £,с , Л,8 •
Векторы-столбцы узловых неизвестных, включающие векторы перемещений узловых точек КЭ и их первые производные в локальной и глобальной системах координат, можно представить в видеА
V/ } = V', V1, Vк,V1 ,У> ,У* ,У1 ,У* V V V } (16)
1x12
V/ } = V7' V1Vк V1 ,у' V,8 VI ^ VI}. (17)
Используя дифференциальные выражения
V, ^ = V,, 5, 5 +У ,С С,! , V Л = V,, 8,Л +У ,С С Л, (18)
между векторами (16) и (17) можно сф ормировать матричное соотношение
V, У }=М V Г }. (19)
Вектор перемещения внутренней точки объемного конечного элемента аппроксимируется через вектор -столбец узловых перемещ ений выражением
V=И^л }= ¥ и V7; }= ыт№}, (20)
1x12 12x1 1x12 1М2 12x1 1x12 12x1
где {у}т - вектор-строка аппроксимирующих функций, элементами которой являются полиномы Эрмита третьей степени:
Ф = ¥ Ф = ¥2; Ф = ¥з; Ф = ¥4; Ф = ¥5 8,\ +¥9 8,Л; Ф = ¥б ^ 1 +¥ 8,Л; ф = ¥7$ \ +¥п8,Л; ф = ¥8 $ 1 +¥12 8,Л; Ф9 = ¥5С, ^ +¥9С,Л; Фо =¥бС,\ +¥\оС,Л ; Ф11 =¥7С,1 +¥11С,Л ; Ф12 =¥8С4 +¥12С,Л•
Производные вектора (20) опр еделяются выражениями:
V';}- Ф }Т V, У} М=
Вектор- столбец узловых векторов перемещений (17) можно представить матричным выражением
& }= ШГ } (21)
-Г I у \
12x1 12x24 24x1
где [а] - матрица, ненулевыми элементами которой являются базисные векторы узловых точек конечного элемента;
У } -{АУ /1 ...V, ^ }. (22)
1x24
С использованием выражения (5) базисные векторы узловой точки т выражаются через базисные векторы внутренней точки конечного элемента:
[§т ]=[$т ][5]-1 {§}=^т ]{§}, (т = ?, k,l). (23)
2x1 2x2 2x2 9x1 2x2 2x1
При использовании выражений (21) и (23) вектор (20) запишется как
V = {§Т ШГ }, (24)
2x24 24x1
где [к] = [р ^ Г ... рА ]р ^ Г ... СР8 ] ... СР9 ^ Г ... СР12 ]
Г' ] Р5 ] ... Р8 [Zl ] ...Р9 ^ ] ... Р12 К
2x24 2x2 2x2 2x2 2x2 2x2 2x2
Производные вектора V определяется выражениями
V,5 = {§} [к,5 ]/Г } (25) V, с = {§}Т[к, с ]/Г }, (26)
1x2 2x24 24x1 1x2 2x24 24x1
где
[к, Ы ]=[Р1, Ы к [z' ]Т...Р12Ы [z'Г]; [к, с ]=Р, с ^ ([z' ]...рп4 ]"].
Сравнивая (7) и (24), можно получить аппроксимирующие выражения для компонент вектора перемещения
V1 ={'1}Т/Г } V ={'2 }Т/Г } (27)
где строки {'1 }Т и{'2 }Т являются первой и второй строками матрицы [к].
Сравнивая правые части (25), (26) и (8), можно получить аппроксимирующие выражения для производных компонент вектора перемещения внутренней точки конечного элемента:
V,! = {'3 }Т {/Г } V,Ы = {/4 }Т {/Г } V, с = {'5 }Т {/Г } V, с = {'6 }Т {/Г } (28)
Столбец неизвестных {/уГ } можно представить выражением
{/Г }=[*]& } (29)
24x1 24x24 24x1
где {¡Г } = , V11,V1, Vй, V,! ,...,V 1?,...,V, V1, V1', V1, V,! ,...v, ^ } - вектор-строка
1x24
компонент вектора узловых перемещений; [^ ] - матрица перехода.
При традиционной - скалярной аппроксимации искомых величин компоненты вектора перемещения внутренней точки конечного элемента и их производные определяются через узловые значения соотношениями
V1 =МЧУ } V =№{>, }; V,1! =Ь ! ^ } V, ! ={р, ! }Т{^ }
1x12 12x1 1x12 12x1 1x12 12x1 1x12 12x1 (30)
4 = {F, с }Т И } v, с = к с }ТК}.
1x12 12x1 1x12 12x1
В аппроксимирующих выражениях (30) отсутствуют параметры, характеризующие конкретную криволинейную систему координат, и каждая компонента вектора перемещения внутренней точки конечного элемента зависит от узловых значений только этой же компоненты.
С использованием соотношений (27), (28) и (29) формируется матрица жесткости кольцевого объемного конечного элемента [9]
[К ][< ]={Р }, (31)
У
1x24 24x1
]- матрица жесткости элемента в глобальной системе координат; {^} — вектор узловых нагрузок элемента в глобальной системе координат.
Пример 1. В качестве примера рассмотрим оболочку, срединная поверхность которой представляет собой форму усеченного эллипсоида вращения. Оболочка нагружалась равномерно распределенным внутренним давлением интенсивности q (рис. 2). В узлах сечения х=0 устанавливались опоры переменной жесткости, позволяющие смещение конструкции как твердого тела. Были приняты следующие исходные данные: А =1м; В=0.5м; D= 0,9 м; Е =2-105 МПа; V =0.3; h = 0.01 м; q = 4 МПа, (А и В являются главными полуосями эллипса). Радиус вращения срединной поверхности оболочки вычислялся по формуле:
г = Вт]А2 — X2 / A, где осевая координата х изменялась в пределах 0 < х < 0.9 м. Расчеты выполнялись в двух вариантах: в первом реализовывалась традиционная скалярная (30) аппроксимация полей перемещений; во втором - предложенная векторная аппроксимация (27), (28). Результаты расчета представлены в таблице, где в левой колонке приведены смещения конструкции как твердого тела при изменении жесткости пружинных опор.
Во второй и третьей колонках даются кольцевые напряжения в сечении х = 0 при использовании скалярной аппроксимации перемещений (вариант I) и векторной аппроксимации (вариант II).
Таблица
Величина жесткого смещения А, м Кольцевые напряжения <к, МПа
Вариант расчета
I II
0.00 126.2 125.6
0.01 99.7 125.6
0.03 42.1 125.6
0.05 -6.02 125.6
4.80 -2056.5 124.7
5.20 -5452.16 123.26
Анализ табличного материала показывает, что значения кольцевых напряжении в сечении х = 0 в обоих вариантах расчета практически совпадают при установке стержневых опор в узлах конструкции сечения х=0 и незначительно отличаются от аналитического решения <к = 125,3 МПа . Заменив шарнирное закрепление на левом конце рассматриваемой оболочки пружинной опорой, можно позволить оболочечной конструкции перемещаться вдоль оси х как жесткому целому на некоторую величину А. Величина смещения конструкции А будет зависеть от выбранной жесткости пружины. Как видно из таблицы напряжения, полученные с использованием традиционной независимой аппроксимации (вариант I) уже при незначительном смещении оболочки как жесткого целого отличаются от первоначальных и с увеличением А становятся неприемлемыми. В то время как напряжения, вычисленные с использованием векторного способа аппроксимации (вариант II), остаются постоянными независимо от А и лишь при величине смещения 5,2 м, что почти в пять раз превышает размеры самой конструкции, получают некоторую погрешность. Анализ представленно-
Рис. 2. Эллипсоид вращения, загруженный внутренней равномерно распределенной нагрузкой
го табличного материала позволяет сделать вывод о том, что предложенный способ инвариант аппроксимации в векторной формулировке полей перемещений объемного кольцевого конечного элемента обладает существенными преимуществами в сравнении с существующей методикой аппроксимации. Использование векторной аппроксимации полей перемещений конечного элемента, в отличие от существующей скалярной аппроксимации, позволяет решить проблему учета смещения конструкции как жесткого целого и эффективно исследовать их напряженно-деформированное состояние.
Л и т е р а т у р а
1. Оден Дж. Конечные элементы в нелинейной механике сплошных сред / Дж. Оден // М.: Изд-во «Мир», 1976. - 464 с.
2. Голованов А.И. Метод конечных элементов в статике и динамике тонкостенных конструкций/ А.И.Голованов, О.Н. Тюленева, А.Ф. Шигабутдинов// М.: Физматлит, 2006. -392 с.
3. Скопинский В.Н. Напряжения в пересекающихся оболочках. / В.Н. Скопинский // М.: Физматлит, 2008. - 400 с.
4. Feng S.Z. Thermo-mechanical analysis of functionally graded cylindrical vessels using edge-based smooth head finite element method/ S.Z.Feng, X. Cui, G.Y. Li, H. Fend, F.X. Xu// Int. J. Pressure Vessels and Pipe. 2013. - 111-112 . - C.302-309.
5. Belalia S.A. Nonlinear free vibration of functionally graded shear deformable sector plates by a curved singular ¿»-element/ S.A. Belalia, A. Houmat// Eur. J. Mech. A. - 2012. -35. - С. 1-9.
6. Marcio A.A. Generalized finite-volume theory for elastic stress analysis in solid mechanics. Pt I,II. Framework . Results/ A.A. Marcio, Pindera Marek-Jerzy// Trans. ASME. J. Appl. Mech. - 2012. - 79, №5. - C. 051007/1-051007/11-13.
7. Кривошапко С.Н. Энциклопедия аналитических поверхностей/ С.Н. Кривошап-ко, В.Н. Иванов// М.: Книжный дом «ЛИБРОКОМ», 2010, 560с.
8. Киселев А.П. Расчет многослойных оболочек вращения и пластин с использованием объёмного конечного элемента / А.П. Киселев, Н.А. Гуреева, Р.З. Киселёва // Изв. Вузов, сер. «Строительство». - 2010. - № 1. - С. 106-112.
9. Седов А.И. Механика сплошной среды. - М: «Наука», 1976, т. 1, 35 с.
References
1. Oden J (1976). Finite Elements in Nonlinear Continuum Mechanics. M.: Izd-vo "MIR", 464 p.
2. Golovanov AI, Tuleneva ON, Shigabutdinov AF (2006). The Finite Element Method in Statics and Dynamics of Thin-Walled Structures. Moscow: Fizmatlit, 392 p.
3 Skopinsky VN (2008). Stresses in Intersecting Shells. Moscow: Fizmatlit, 400 p.
4. Feng SZ, Cui X, Li GY, Fend H, Xu FX (2013). Thermo-mechanical analysis of functionally graded cylindrical vessels using edge-based smooth head finite element method. Int. J. Pressure Vessels and Pipe. 111-112, p. .302-309.
5. Belalia SA, Houmat A (2012). Nonlinear free vibration of functionally graded shear deformable sector plates by a curved singularp-element. Eur. J. Mech. A, 35, p. 1-9.
6. Marcio AA, Pindera Marek-Jerzy. Generalized finite-volume theory for elastic stress analysis in solid mechanics. Pt I, II. Framework. Results. ASME, J. Appl. Mech., 79, № 5, p. 051007/11-13.
7. Krivoshapko SN, Ivanov VN (2010). Encyclopedia of Analytic Surfaces. Moscow: "LIBRO-KOM", 560 p.
8. Kiselev AP, Gureeva NA, Kiseleva RZ (2010).Calculation of multilayer shells of revolution and plates using the volumetric finite element. Izv. Vuzov, Ser. "Stroitel 'stvo", № 1, p. 106-112.
9. Sedov AI (1976). Continuum Mechanics. Moscow: "Nauka", Vol. 1, p. 35.
ACCOUNT OF THE SHIFT AS RIGID BODY OF SHELL OF REVOLUTION AXIALLY SYMMETRIC LOADED ON THE BASE OF FEM
A.P. Kiselev, R.Z. Kiseleva, A.P. Nikolaev
The displacements of nodal points and their first derivatives were selected as unknown nodal finite element parameters. An example of calculation shows that the presented invariant approximation of displacement fields in the vector formulation allows solving the well-known finite element problem as a rigid body displacement in general.
KEY WORDS: FEM, invariant approximation of displacement fields, vector formulation, axisymmetrically loaded shell, circular volume finite element, rigid body displacement.
