CC BY
17
Численные методы расчета конструкций
МАТЕМАТИЧЕСКАЯ МОДЕЛЬ ДИСКРЕТИЗАЦИИ ЦИЛИНДРА С ЭЛЛИПТИЧЕСКИМ ДНИЩЕМ ИЗ РАЗНОРОДНЫХ МАТЕРИАЛОВ1
Ю.В. КЛОЧКОВ, доктор техн. наук, профессор Т.А. КИСЕЛЕВА, кандидат техн. наук, старший преподаватель О.В. ВАХНИНА, кандидат техн. наук, доцент Волгоградский государственный аграрный университет 400002, Волгоград, пр. Университетский, 26, Klotchkov@bk. ги
Изложен алгоритм расчета произвольных сочлененных оболочек с различными физико-механическими характеристиками на основе МКЭ с использованием скалярной и векторной интерполяции полей перемещений. В качестве элемента дискретизации используется четырехугольный криволинейный конечный элемент с восемнадцатью степенями свободы в узле. Выполнен анализ НДС тонкостенной оболочки из разнородных материалов в форме цилиндра с эллиптическим днищем.
КЛЮЧЕВЫЕ СЛОВА: сочлененные оболочки, скалярная интерполяция, векторная интерполяция, четырехугольный конечный элемент, цилиндр с эллиптическим днищем.
Тонкостенные конструкции из сочлененных оболочек находят широкое применение в различных отраслях народного хозяйства. В частности, цилиндры с эллиптическими днищами используются в качестве сосудов, ёмкостей, вакуумных аппаратов, воздухосборников, котлов и фрагментов магистральных трубопроводов. Геометрия современных тонкостенных конструкций отличается широким многообразием форм [1]. Однако расчет такого рода конструкций возможен лишь при использовании численных методов анализа [2, 3], в частности метода конечных элементов (МКЭ) [4, 5, 6, 7]. Поэтому возникает необходимость совершенствования методов определения напряженно- деформированного состояния (НДС) не только оболочек вращения, но и произвольных оболочек, а также их сочленений. В связи с этим в статье для сочлененных произвольных оболочек из различных материалов предложены кинематические и статические условия сопряжения оболочек, необходимые для исследования напряженно-деформированного состояния такого рода конструкций с помощью высокоточных конечных элементов четырехугольной формы [8], матрицы жесткости (72x72) которых формировались на основе скалярного и векторного способа аппроксимации полей перемещений [9].
Рассмотрим реализацию алгоритма на примере оболочки в форме эллиптического цилиндра, сочлененного с половиной трехосного эллипсоида.
Срединные поверхности цилиндра и эллипсоида (рис. 1) в исходном состоянии определяются соответственно радиус-векторами
Я0 = х' + Кв^тв; + г{в)со$,вк; (1)
Я 0' = х' Г'+г (х\в)$,тбЦ'+г (х\в)со$,вк\ (2)
где 1, ;, k - орты декартовой системы координат основной оболочки (цилиндра); 1', ;', k' - орты декартовой системы координат примыкающей оболочки
Работа выполнена при финансовой поддержке РФФИ и Администрации Волгоградской области по проекту «Разработка математических моделей для расчета на прочность тонкостенных конструкций, входящих в структуру строительных, машиностроительных, нефтехимических и водохозяйственных объектов Волгоградской области» № 15-41-02125 р_Поволжье_а.
(эллипсоида); г(х,в) = ^ Х / ° Ьс , в - угол, отсчитываемый от оси
л1с2 2 в + Ьв Ог против хода часовой стрелки в плоскости, перпендикулярной оси Ох, тогда
Рис. 1
Г,х (Х,в) =-
bc
a 2Л--
Í
c2 sin2 в + b2 cos2 в
x
rn( Х,в) = . 1 - —
bcip2 - c2 )sin 2в
2 3 '
a 2¡c2 sin2 в+ b2 cos2 в)2
(3)
На кривой пересечения срединных поверхностей рассматриваемых произвольных оболочек справедливы следующие геометрические соотношения
г = z,
У = У, Х = х - D,
(4)
где D - длина цилиндра.
Так как рассматривается задача пересечения половины эллипсоида с эллиптическим цилиндром, то х' = 0 на кривой пересечения. В силу этого справедливы равенства
Г(в) = , _ _ЬС _ , ; (5)
Í
c2 sin2 в + b2 cos2 в
2___2,
Г, Х (в) = 0;
г,в(в) =
bc(b2 - c2 )sin 2в
2\c2 sin2 в+ b2 cos2 в)2
3 '
(6)
(7)
Касательные векторы локального базиса определяются выражениями
a10 = R0x = 7 ;
aj = Re = (r, в (в) sin в + r(e)cos в)] + (г,в (в) cos в - r(e)sine)fc. Отсюда векторное произведение
x
2
x
2
a
- 0 - 0 oj х a2 =
i j k 1 0 0 0 re(6)sin6 + r(6)cos6 r6(6)cos6- r(6)sin6
= (r(6) sin 6 - r6 (6)cos 6)j + (r(6) cos 6 + r6 (6)sin 6)k. (9)
Ковариантные компоненты метрического тензора определяются скалярными произведениями базисных векторов
ail = 1; ai2 = 0; 022 = r,6(6) + r 2(6). (10)
На основе этого орт нормали к срединной поверхности основной оболочки в исходном состоянии равен
-0 o0 хo2 r(6)sin6-r0(6)cos6 - r(6)cos6 + r0(6)sin6 -
a = i—2 =-, ' =-j +-, ' =-k, (11)
Va 0 vr,6(6) + r 2(6) Vr,6(6) + r 2(6)
где
0 a11022 (a12 ^
- детерминант метрического тензора, отнесенного к срединнои поверхности основной оболочки.
В точках, расположенных на кривой пересечения срединных поверхностей основной и примыкающей оболочек, введем следующие столбцы узловых неизвестных
и Г = ÍV'V2vv,11v,2v,1v,12v,22v,2v»,11v>,22V,22V»,12 } (12)
1x18
} — íu и ии 2 и 2 и 2 и 22U 22U 22U12^ 12U12 j; (13)
1x18
fry Í = ÍV' 'V 2'V' V11 'V,2 ' V,1' V2 ' V,22 ' V,2 ' ^11 ' V,U ' V,11 ' V22 ' V,222 ' V,22 ' ^12 ' V,22 ' V,12 '} (14) 1x18
"^Uyv } — íu и и и 1 и 1 и 1 и 2 и 2 и 2 и 11 u 11 u 11 и 22 и 22 u 22 u 12 и 12 u 12 }, (15)
1x18
где
vР, V р\ vр, ^, ^ ^р \ уар, уар'
- компоненты вектора перемещения и их первые и вторые производные в базис-
- 0 - 0 I
ных векторах ар и ар ;
ир, ир\ ир, ирр\ иа, иа \ ирр, иррр \ uaр, ир'
- компоненты вектора перемещения и их первые и вторые производные, отнесенные к базисным векторам ер и ер'.
При формировании матрицы жесткости конечного элемента в узлах, расположенных на кривой пересечения, в качестве неизвестных принимаются компоненты и многочлены столбца (13). После составления матрицы жесткости дискретного элемента основной оболочки, узлы ], k которого расположены на кривой пересечения выполняется, преобразование указанной матрицы жесткости и столбца усилий рассматриваемого элемента, обусловленное переходом от столбца узловых неизвестных (12) к столбцу (13). При формировании матрицы жесткости конечного элемента примыкающей оболочки, узлы которого V, I' расположены на кривой пересечения необходимо выполнить преобразования,
обусловленные переходом от столбца узловых неизвестных (14) к столбцу (15), а затем к столбцу (13).
Связь между столбцами (12) и (13) в узловой точке кривой пересечения произвольных оболочек может быть представлена в матричном виде
Ьг)=и]к), (16)
для основной оболочки и в матричном виде
ЬГ 'МР!']^'} (17)
для примыкающей оболочки.
В узловой точке на кривой пересечения получены 14 соотношений для выражения 14 компонент столбца узловых неизвестных примыкающей оболочки через компоненты столбца основной оболочки. Столбец неизвестных узла кривой пересечения может быть представлен в виде
|| ^ * II 1111 1 1 22222 2 11 2 1 1
уу } = ? Ч {2{11{22{12и Н ^2 ^11^22^12^ t2 {11{22{12{1 Н 117"2'^2', (18) _1 _1 _2 _1
где многочлены tl', tll', tl ' и tl2' остаются свободными и будут варьироваться как неизвестные примыкающей оболочки.
Можно составить матричное соотношение между столбцами (15) и (18)
к 'М р2 *} (19)
18x1 18x22 22x1
Многочлены , , могут быть выражены через компоненты
соответствующих векторов перемещений основной и примыкающей оболочки и их первые и вторые производные. Тогда связь между столбцами (13) и (15) может быть представлена в матричном виде
к'}=Мк}. (20)
Для получения зависимости (20) используются следующие условия сочленения.
1. Инвариантность векторов перемещения точек срединных поверхностей V' = V в узлах на кривой пересечения позволяет получить зависимости
V1' =-у!; V2' =-у2; V'=V. (21)
2. Равенство первой и второй производных по координате в векторов перемещений
Кв' = ^ ^вв' = Квв
в узлах на кривой сочленения определяет соотношения
1.12,2 , 1, 1 2, 2 ,
V,в = ^в; V,в = V,в; ^в = -^в; V,вв = -V,вв; V,вв = -^вв; ^вв = в (22)
3. Равенство углов поворотов нормалей к срединным поверхностям основной и примыкающей оболочек в точках кривой пересечения позволяет получить
V х' = х. (23)
4. Равенство нормальных усилий на кривой сочленения
" (£11 '+у£22 ')= (£11 + у£22 )
1 -у'^11 ' 1 -V2
дает возможность получить зависимость
■3 '
х , ЕН3 (1 -V'2 \ / 1 2 \ , 2 ,
V1« /Г 2)^ + )-У v20■ (24)
Е' К3 (1 -V
5. Равенство касательных усилий в узловых точках кривой пересечения
Е'# ,_ Ек_ 1 + у'_ 1 + у^12
приводит к выражению
V2,' _ ^^. (25)
'* Еh'(1 + у) '*
6. Равенство моментов Мхх _ Мхх' дает возможность выразить вторую производную V,,' функцией узловых неизвестных основной оболочки
/ (к})
.у„ (26)
Матрица жесткости и столбцы узловых усилий конечных элементов примыкающей оболочки, расположенные на границе пересечения произвольных оболочек должны быть преобразованы в результате матричных произведений
^гранич
_ [Р2 IКIР2 I [/Г"41_ [Р2 I {/у } (27)
Г2Ь У у
Пример расчета.
В качестве примера решим задачу об определении НДС цилиндра с эллиптическим днищем, загруженного внутренним давлением интенсивности q .
Оболочка представляет собой цилиндр длиной D _ 4,0 м и радиусом R _ 0,9 м (рис. 1), сочлененный с эллипсоидом вращения с параметрами а _ 1,3 м, Ь _ с _ 0,9 м ; толщина оболочки t _ 0,008 м ; коэффициент Пуассона у =0,3; д _ 1 МПа.
Расчеты выполнялись в двух вариантах: в первом варианте при формировании матрицы жесткости четырехугольного элемента дискретизации использовалась общепринятая интерполяция компонент вектора перемещения как скалярных величин; во втором варианте применялась векторная интерполяция полей перемещений четырехугольного конечного элемента. В обоих вариантах были реализованы разработанные условия сочленения оболочек.
Результаты повариантных расчетов представлены в таблице 1, в которой приведены численные значения меридиональных и кольцевых напряжений в опорном, сочлененном и концевом сечениях оболочки в зависимости от отношения модулей упругости сочленяемых оболочек, при этом модуль упругости цилиндра Е1 последовательно уменьшался, а модуль упругости эллипсоида
Е2 _ 2 • 10 МПа оставался неизменным.
Таблица 1
Интерполяция Скалярная Векторная
Е1 1 1 1 1 1 1 1 1
сечения 1 10 100 1000 1 10 100 1000
< 47,92 47,9 47,9 47,4 47,9 47,9 47,9 47,9
а о < 47,92 47,9 47,9 47,4 47,9 47,9 47,9 47,9
о
X < 112,5 112,5 112,5 112,5 112,5 112,5 112,5 112,5
ак 112,5 112,5 112,5 112,5 112,5 112,5 112,5 112,5
< 461 128,9 199,7 210,8 46,2 129,3 200,2 211,6
а о < 49,9 -31,3 -101,1 -113,1 49,8 -31,6 -101,6 -112,8
X < 100,5 63,7 63,5 63,7 100,5 63,8 63,6 63,8
< 101,6 15,7 -26,7 -33,5 101,6 15,5 -26,9 -33,5
< -0,08 -11,64 -128,6 -1256,1 0,03 0,03 0,03 0,03
а сч < 0,12 11,69 128,8 1257,4 0,05 0,05 0,05 0,05
ич" II X < 84,2 72,5 -45,0 -1177,9 84,5 84,5 84,5 84,5
< 82,3 76,3 5,0 -682,2 83,6 83,6 83,6 83,6
Анализ табличного материала показывает, что результаты повариантного расчета практически совпадают при равенстве модулей упругости основной (цилиндра) и примыкающей (эллипсоида) оболочек. Однако с уменьшением значения модуля упругости цилиндра, величины контролируемых параметров НДС в концевом сечении оболочки существенно отличаются друг от друга в зависимости от варианта расчета. Этот факт более наглядно можно проследить на диаграммах, представленных на рис. 2.
х=5,2 м
св
С
и к к
о
я
и Л
е
св
К
200 0
-200 -400 -600 -800 -1000 -1200 -1400
ей
С
и к к
о
я
и Л
е
ей
К
Величина отношения Е1/Е2 ♦ скалярная интерполяция векторная интерполяция
х=5,2 м
200 0 -200 -400 -600 -800 -1000 -1200 -1400
Величина отношения Е1/Е2 ♦ скалярная интерполяция векторная интерполяция
Рис. 2
Поскольку цилиндрическая часть оболочки становится все более податливой, то эллипсоидная часть оболочки получает дополнительную возможность смещаться как абсолютно твердое тело. Такое смещение не должно оказывать влияние на НДС конструкции, что наблюдается во втором варианте расчета, в котором реализована векторная интерполяция полей перемещений. В первом же варианте погрешность вычислений в концевом сечении стремительно нарастает. Меридиональные напряжения, которые по физическим соображениям должны быть равны нулю, так как концевое сечение незагружено, достигают неприемлемо больших значений. Кольцевые напряжения изменяют свой знак, что также является недопустимым. Таким образом, можно сделать вывод, что при расчете сочлененных оболочек с различными физико-механическими характеристиками материала необходимо использовать векторную интерполяцию перемещений в сочетании с разработанными корректными условиями сочленения.
H h t e p a t y p a
1. Krivoshapko S.N., Ivanov V.N. (2015) Encyclopedia of Analytical Surfaces. - Springer International Publishing Switzerland, 2015. - 752 p.
2. Gout C., Romani L. (2013) Multivariate Approximation: Theory and applications// Journal of Computational and Applied Mathematics, Vol. 240, pp. 1-4.
3. Guyader C. Le, Apprato D., Gout C. (2012) On the construction of topology-preserving deformation fields// IEEE Transactions on Image Processing, 21 (4), art. no. 6093964 , pp. 1587-1599.
4. Munjiza A., KnightE.E., Rougier E. (2015) Large strain finite element method: A practical course. - Wiley, 2015. - 486 p.
5. Kim D.-J., Duarte S.A., Sobh N.A. (2011) Parallel simulations of three-dimensional cracks using the generalized finite element method// Comput. Mech. Vol. 47. No. 3, pp. 265282.
6. Chen X., Liu Y. (2015) Finite Element Modeling and Simulation with ANSYS Workbench. - CRC Press, 2015. - 411 p.
7. Dhatt G., Lefrangois E., Touzot G. (2012) Finite Element Method. - John Wiley & Sons, 2012. - 612 p.
8. Klochkov Yu.V., NikolaevA.P., Kiseleva T.A. (2015) A comparative evaluation of the scalar and vector approximations of sought quantities in the finite-element method of arbitrary shells// Journal of Machinery Manufacture and Reliability, Vol. 44, No. 2, pp. 166-172.
9. Kiseleva T.A., Klochkov Yu.V., NikolaevA.P. (2015) Comparison of scalar and vector FEM forms in the case of an elliptic cylinder// Computational Mathematics and Mathematical Physics, Vol. 55, No. 3, pp. 422-431.
R e f e r e n c e s
1. Krivoshapko, S.N., Ivanov, V.N. (2015). Encyclopedia of Analytical Surfaces, Springer International Publishing Switzeland, 752 p.
2. Gout, C., Romani, L. (2013). Multivariate Approximation: Theory and applications, Journal of Computational and Applied Mathematics, Vol. 240, pp. 1-4.
3. Guyader, C. Le, Apprato, D., Gout, C. (2012). On the construction of topology-preserving deformation fields, IEEE Transactions on Image Processing, 21 (4), art. no. 6093964, pp. 1587-1599.
4. Munjiza, A., Knight, E.E., Rougier, E. (2015). Large Strain Finite Element Method: A practical course. - Wiley, 486 p.
5. Kim, D.-J., Duarte, S.A., Sobh, N.A. (2011). Parallel simulations of three-dimensional cracks using the generalized finite element method, Comput. Mech., Vol. 47. No. 3, pp. 265-282.
6. Chen, X., Liu, Y. (2015). Finite Element Modeling and Simulation with ANSYS Workbench, CRC Press, 411 p.
7. Dhatt, G., Lefrangois, E., Touzot, G. (2012). Finite Element Method, John Wiley & Sons, 612 p.
8. Klochkov, Yu.V., Nikolaev, A.P., Kiseleva, T.A. (2015). A comparative evaluation of the scalar and vector approximations of sought quantities in the finite-element method of arbitrary shells, Journal of Machinery Manufacture and Reliability, Vol. 44, No. 2, pp. 166-172.
9. Kiseleva, T.A., Klochkov, Yu.V., Nikolaev, A.P. (2015). Comparison of scalar and vector FEM forms in the case of an elliptic cylinder, Computational Mathematics and Mathematical Physics, Vol. 55, No. 3, pp. 422-431.
MATHEMATICAL MODEL OF DISCRETISATION OF CYLINDER WITH ELLIPTIC BOTTOM FROM HETEROGENEOUS MATERIALS
KLOCHKOV Yu.V., KISELEVA T.A., VAKHNINA O.V.
The algorithm for calculating arbitrary articulated shells with different physical - mechanical characteristics based on the finite element method with using of the interpolation of scalar and vector displacement fields is presented. As part of the sampling, rectangular curvilinear finite elements of eighteen degrees of freedom in the node were used. The analysis of VAT thin shell of dissimilar materials in the form of a cylinder with elliptical bottom was fulfilled.
KEY WORDS: articulated shell, scalar interpolation, spatial interpolation, quadrangular finite element, a cylinder with elliptical bottom.
