CC BY
15
Численные методы расчета конструкций
СРАВНИТЕЛЬНЫЙ АНАЛИЗ СКАЛЯРНОЙ И ВЕКТОРНОЙ ФОРМ АППРОКСИМАЦИЙ В МКЭ НА ПРИМЕРЕ СООТНОШЕНИЙ В.В. НОВОЖИЛОВА ДЛЯ ЭЛЛИПТИЧЕСКОГО ЦИЛИНДРА
Ю.В. КЛОЧКОВ, доктор техн. наук, профессор А.П. НИКОЛАЕВ, доктор техн. наук, профессор ТР. ИЩАНОВ, аспирант
Волгоградский государственный аграрный университет 400002, г. Волгоград, пр. Университетский, 26, e-mail: ishchanov. volgau@yandex. ru
Изложен сравнительный анализ скалярной и векторной аппроксимаций полей перемещений на примере конечно-элементного расчета эллиптического цилиндра при использовании теории тонких оболочек Новожилова В.В. В качестве элемента дискретизации применяется четырехугольный криволинейный конечный элемент в виде фрагмента срединной поверхности эллиптического цилиндра с восемнадцатью степенями свободы в узле. На численных примерах показано, что векторная аппроксимация, являясь инвариантной, обладает принципиальными преимуществами по сравнению со скалярной аппроксимацией при расчете произвольных оболочек.
КЛЮЧЕВЫЕ СЛОВА: векторная аппроксимация, скалярная аппроксимация, конечный элемент, эллиптический цилиндр.
В настоящее время широкое распространение получили численные методы расчета, среди которых наибольшую известность приобрел метод конечных элементов (МКЭ) [1 - 15]. В большинстве современных вычислительных комплексов (ANSYS, ABACUS, NASTRAN и др.), реализована скалярная интерполяционная процедура, суть которой заключается в аппроксимации отдельной компоненты вектора перемещения через узловые значения этой же компоненты. Для расчета оболочек произвольной формы со значительными значениями кривизн срединной поверхности или при наличии жестких смещений необходимо применять векторную аппроксимацию, суть которой заключается в использовании аппроксимирующего выражения непосредственно для вектора перемещения произвольной точки конечного элемента через векторы перемещений узлов и их производные [16].
1. Геометрия эллиптического цилиндра
Радиус-вектор точки срединной поверхности эллиптического цилиндра
может быть задан в декартовой системе координат в виде
^ ^ ^ ^
R0 = xi + b • cos(t) j + c • sin(t) k, (1)
где b и c - полуоси поперечного сечения эллиптического цилиндра (y2 / b2 + z2 / c2 = 1), t - параметр, отсчитываемый от оси OY против хода часовой стрелки.
Дифференцируя (1) по x, t и выполняя нормализацию векторов локального
базиса, можно получить орты, касательные к линиям главных кривизн:
^ ^ ^ ^ ^ ^ ^
ej0 = R,0x = i; e20 = R = (-b • sin(t) / d) j + (c • cos(t) / d) k;
^ ^ ^
e0 = (-c • cos(t)/d) j + (-b • sin(t)/d) k, (2)
где d = -Jc2 cos2(t) + b2 sin2(t).
Деформации и искривления в произвольной точке срединной поверхности могут быть выражены через компоненты вектора перемещения и их производные [17]:
е - —• е - — — -ку е - — — | ду2 • 11 дх 22 А д 12 А2 д дх
Хм —'
д 2у
1 д 2у ду
Х22 — '
дх2 д(1/А2) ду
А д дг
; Х12 — '
2
А2 дхд» дх
-к;
(3)
2
V А2
д у 1 ду
дг2
. . 1 к
----к----у
А дг А2 дг
где е„, е
11, е22, е12 - относительные удлинения и сдвиг срединной поверхности; Х11, Х22, Х12 - искривления и кручение срединной поверхности оболочки в процессе деформирования; х - осевая координата; А2 — ^(ср, )2 + {у/, )2 ; р — Ь ■ ео8(»),
у — с ^т(г); у1, у2, у - тангенциальные и нормальные компоненты вектора перемещения точки срединной поверхности; к - кривизна дуги поперечного сечения эллиптического цилиндра.
2. Элемент дискретизации и аппроксимация искомых величин В качестве элемента дискретизации выбирается четырехугольный фрагмент срединной поверхности эллиптического цилиндра с узлами
г, у, к, I, отображаемый для удобства численного интегрирования на квадрат с локальной системой координат -1 < £,ц < 1.
Глобальные криволинейные координаты в1 — х и в2 —г произвольной точки срединной поверхности криволинейного конечного элемента определяются через их узловые значения с использованием билинейных функций координат £ и Г)
ва —№л)У ву } (4)
где в У } —{вт вт вак ва1} - строка узловых значений координат ва .
Дифференцированием (4) определяются производные глобальных координат в локальной системе координат в", ву, в^ и локальных координат в глобальной системе , )вт , £ваве , ) ^ .
Столбец узловых варьируемых параметров в локальной £, п и в глобальной х, г системах координат выбирается в виде
и; Г—у }т у? } {у; }); Ш—у } {уГ } у}
У
1x24
У 1x24
У
1x24
У
1x24
У
1x24
У 1x24
(5)
где
) л ^ ) г У к I г У к I г I г I г I
Цу) — Ч Ч Ч Ч Ч,I Ч,у Ч,к Ч,ц Чщ ...Чщ Ч^ -Ч^ Чщ) ..Чщ) Ч,щ ..Ч-
1x24
..Ч,|Л
{чг } — „г „г „г „г „г
у
1x24
Ч Ч Ч Ч, xЧ, xЧ, xЧ, хЧ в ...Ч,в Ч, хх ..Ч, хх Ч, и .. Ч,и Ч, х» Ч
,хг .Ч,хг }
Здесь и ниже под Чт (т — г, У, к, I) понимается компонента узлового вектора перемещения У1т, V2т или Ут .
2
При скалярном варианте аппроксимации искомых величин каждая компонента вектора перемещения точки внутренней области элемента дискретизации и ее производные, аппроксимируются через узловые значения этой же компоненты и ее производные, оставаясь независимой от остальных двух компонент [1,12,13]:
q = И }, q,4 = KF te7}, = fo}. (6)
1x24 24x1 1x24 24x1 1x24 24x1
При вектором способе интерполяции искомых величин в качестве узловых варьируемых параметров принимаются непосредственно векторы перемещения узлов элемента дискретизации и их производные в локальной и глобальной системах координат [16]:
U Л 1 Ui aj Ak Аl ^^^^^ А А A A III i l i l i l i l i l I
1 Vy V = 1 V V V V V,i..V,i ..v^..v,ii v,w..^„..^ \,
i j k l
WT. i WT. l WT. i WT. / т i WT. l WT. i WT. l WT. i
(7)
1x24
При векторном варианте аппроксимации искомых величин вектор перемещения точки внутренней области конечного элемента аппроксимируется через векторы узловых точек соотношением [16]:
' 1
V |„ I („ЛМг 1 |„ I С,.Лт Г ГАП
V = ] Vу I = Ыт И 1 Vу 1 = [2] [С] [Щ {у^ },
1x24 I I 1x24 24x24 I I 1x24 24х72 72 х72 72 х72 72х1 24x1 24x1
у"Гх = V,х }т А [С] [Щ Ц/уГ }, уТх, = V,„ }т А [С] [N] Ц/' }, (8)
1x24 24х72 72 х72 72 х72 72xl 1x24 24х72 72 х72 72 х72 72xl
А
где [А] - квазидиагональная матрица, содержащая векторы локальных базисов
24 х72
узловых точек {ет} = (е^" е!Т е3т}; {у^} - столбец узловых варьируемых пара-
72x1
метров в глобальной системе координат, [С] , [Щ] - матрицы, содержащие уз-
72 х72 72 х72
ловые значения геометрических характеристик.
Из (8), используя выражения базисных векторов узловых точек через ба-
А А
зисные векторы внутренней точки конечного элемента (ет } = [X] { е }, (т= i, ], к, I) можно получить выражения
А = |А) [к] }, С = |АI кх ] }, у!, = |А) [h,„ ] и }. (9)
I I Зх72 72x1 I ) 3х72 72x1 I ) 3х72 72x1
1x3 1x3 1x3
Из (9) определяются матричные выражения для компонент вектора перемещения и их производных
V1=мт иг}, V2=к }т к}, V={«}т иг },у>={/*т и},
1x72 72x1 1x72 72x1 1x72 72x1 1x72 72x1
V, = К}т {К }, V*; = & }т {игу V = }т {игу } . (10)
,xt
1x72 72x1 1x72 72x1 1x72 72x1
Как видно из (9), при векторном способе аппроксимации искомых величин каждая компонента вектора перемещения точки внутренней области конечного
элемента зависит от полного набора узловых варьируемых параметров {у^ }, в
структуру которых входят узловые значения всех трех компонент вектора перемещения и их производные. При скалярном варианте аппроксимации искомых величин отдельная компонента вектора перемещения внутренней области конечного элемента зависит от узловых значений только этой же компоненты и не зависит от узловых значений остальных двух компонент.
Матрица жесткости и столбец внешней нагрузки формируются на основе функционала Лагранжа стандартным образом [1]
(11)
} Наг = {{иГ \F\dF,
V F
где М= [С]{^ }, ^ }= [Г]М, {и Г = {уУ2у}.
Пример расчета 1.
В качестве примера была решена задача по определению напряженно-деформированного состояния бесконечно длинного эллиптического цилиндра, нагруженного вдоль образующей линейно распределенной нагрузкой, интенсивности q и имеющего на диаметрально противоположной образующей цилиндра шарнирные опоры, препятствующие вертикальному смещению (рис. 1). Были приняты следующие исходные данные: q=5• 10-5 МПа; Ь =0,1 м; Е=2-105 МПа; v=0,3; h=0,001 м. Величина параметра с варьировалась от 0,1 до 0,025 м. Вследствие наличия плоскостей симметрии рассчитываемая оболочка моделировалась полоской конечных элементов шириной Дх=0,01 м. Расчеты проводились в двух вариантах: в первом варианте использовалась скалярная аппроксимация искомых величин; во втором варианте была реализована векторная аппроксимация полей перемещений.
Рис. 1.
Первоначально величина параметра с был принята равной 0,1 м, то есть цилиндр был круговым. Результаты повариантного расчета при с =0,1 м представлены в таблице 1, в которой приведены численные значения нормальных
„ В ~ „я
напряжений на внутренней ои и наружной ои поверхностях цилиндра в точках 1 и 2 в зависимости от густоты сетки дискретизации оболочки.
Анализ табличного материала, представленного в таблице 1, показывает, что результаты повариантного расчета кругового цилиндра, практически совпадают между собой. Кроме того, следует отметить полное совпадение значений напряжений в точках 1 и 2 в обоих вариантах расчета.
Таблица 1
Координаты точек (х, Напряжения, МПа Варианты аппроксимации
скалярная векторная
17x2 33x2 65x2 17x2 33x2 65x2
Точка 1 (1) оВ 190.98 190.97 190.99 190.98 190.99 190.99
оН -190.99 -190.97 -190.99 -190.99 -190.99 -190.99
Точка 2 («НО оВ 190.98 190.97 190.99 190.98 190.99 190.99
оН -190.99 -190.97 -190.99 -190.99 -190.99 -190.99
Пример расчета 2.
Если параметр с последовательно уменьшать, то цилиндр будет принимать все более выраженную эллиптическую форму. Результаты повариантных расчетов эллиптического цилиндра при сетке дискретизации (81^2) представлены в формах диаграмм (рис. 2 и 3), в которых приняты следующие обозначения:
п = Ь/с; Х = о1 /о^в; ц = о1 /о2, верхние индексы 1 и 2 указывают на номер точки оболочки, в которой контролируется НДС эллиптического цилиндра.
Я 3
- Скалярная I гнт ерполящ 1я
(81x2) -Векторная I шт ерполящ 1я (81^2)
Из анализа диаграмм видно, что по мере увеличения отношения параметров Ь / с значения коэффициентов X и ц при векторной аппроксимации остаются постоянными и равными единице, а при скалярной аппроксимации искомых величин значения коэффициенты X и ц возрастают с увеличением п.
Согласно общепринятому в МКЭ положению, погрешность конечно-элементных решений должна уменьшаться при последовательном сгущении сетки, поэтому значительный интерес представляет попытка получить удовлетворительные результаты при скалярном варианте аппроксимации при весьма значительном сгущении сетки. В табл. 2 представлены значения напряжений при густоте сетки большей, чем представленной на диаграмме. Из таблицы
видно, что, несмотря на значительное сгущение сетки приемлемых результатов
Таблица 2
Координаты точек (x, t) Напряжения, МПа Варианты аппроксимации
скалярная векторная
65x2 129x2 257x2 65x2 129x2 257x2
Точка 1 (0; I) < 223.29 223.29 223.29 158.51 158.51 158.51
-226.98 -226.98 -226.98 -158.51 -158.51 -158.51
Точка 2 N) < 44.61 44.61 44.61 158.51 158.51 158.51
-48.45 -48.45 -48.45 -158.51 -158.51 -158.51
Анализируя данные таблиц и диаграмм, можно сделать вывод, что при использовании соотношений Новожилова В.В. векторная аппроксимация эффективна как для оболочек вращения, так и для оболочек произвольного вида. Скалярная аппроксимация дает удовлетворительные результаты только при круговом цилиндре. При расчете эллиптического цилиндра результаты не удовлетворяют условию равновесия, сгущение сетки дискретизации не устраняет погрешности данного способа аппроксимации.
Л и т е р а т у р а
1. Постнов В.А., Хархурим И.Я. (1974). Метод конечных элементов в расчетах судовых сооружений. - Л.: Судостроение. - 344с.
2. Tang C.B., Liul S.Y., Zhou G.X., Yu J.H., Zhang G.D., Bao Y.D., Wang Q.J. (2012). Nonlinear finite element analysis of three implant- abutment interface designs// Int J Oral Sci, 4, №2, р. 101-108.
3. MangS., Vipperman J.S. (2013). Mutual validation of muffler performance evaluation methods using field measurements and finite element analysis// J Acoust Soc Am, 134, №5, р. 4217.
4. Li S.Y., Xu B.G., Tao X.M., Feng J. (2011). Numerical analysis of the mechanical behavior of a ring-spinning triangle using the Finite Element Method// Text Res J, 81, р. 959971.
5. Peng Shang, Xin Ye, Xueling Bai, Zhentao, Hou, Zhaobin Xu and Linan Zhang (2014). Hub bearing fatigue life prediction using finite element method// WIT Transactions on Information and Communication Technologies, 48, p. 1227-1233.
6. Bing B., Pang Y., Wang B., Yang A. (2014). Finite Element Analysis of Strip Steel Jitter For Cooling Tower// Mechanical Research and Application, 1, p. 30-32.
7. Zheng-wei L. U. (2014). Static and Dynamic Characteristic Finite Element Analysis of Slide Bearing Block Based on ANSYS// Mechanical Research and Application, 1, p. 33-35.
8. Jing J. (2014). Finite Element Analysis of Wheels// Mechanical Research and Application,!, p. 102-104.
9. Kyoung-Sik Chun, Samuel Kinde Kassegne, Berhanu Kebede Wondimu. Hybrid/mixed assumed stress element for anisotropic laminated elliptical and parabolic shells// Finite Elements in Analysis and Design. - 2009. - 45. - P. 766-781.
10.Якупов Н.М., Серазутдинов М.Н. (1993). Расчет упругих тонкостенных конструкций сложной геометрии. - Казань: ИМН РАН. - 206с.
11.Матвиенко Ю.Г., Чернятин А.С., Разумовский И.А. (2013). Численный анализ несингулярных составляющих трехмерного поля напряжений в вершине трещины смешанного типа// Проблемы машиностроения и надежности машин. - №4. - С. 40 - 48.
12.Бате, К.Ю., Шидловский В.П. (пер. с англ.); Турчак Л.И. (ред.). (2010). Методы конечных элементов. - М.: Физматлит. - 1022 с.
13.Голованов А.И., Тюленева О.Н., Шигабутдинов А.Ф. (2006). Метод конечных элементов в статике и динамике тонкостенных конструкций. - М.: Физматлит. - 391 с.
14.Косицын С.Б., Чан С.Л. (2013). Анализ напряженно-деформированного состояния пересекающихся цилиндрических оболочек при упруго-пластических деформациях с учетом геометрической нелинейности// Строительная механика инженерных конструкций и сооружений, №1. c. 3 - 9.
15. Агапов В.П. ( 2000). Метод конечных элементов в статике, динамике и устойчивости пространственных тонкостенных подкрепленных конструкций. - М.: Издательство АСВ. - 152 с.
16.Николаев А.П., Клочков Ю.В., Киселев А.П., Гуреева Н.А. (2012). Векторная интерполяция полей перемещений в конечно-элементных расчетах оболочек. - Волгоград: ФГБОУ ВПО Волгоградский ГАУ «Нива». - 264 с.
17.НовожиловВ.В. (1962). Теория тонких оболочек. - Л.: Судостроение. - 431с.
References
1. Postnov, VA, Harhur, IJ (1974). The Finite Element Method in the Calculation of Ship Structures. L .: Shipbuilding, 344p.
2. Tang, CB, Liul, SY, Zhou, GX, Yu, JH, Zhang, GD, Bao, YD, Wang, QJ (2012). Nonlinear finite element analysis of three implant- abutment interface designs. Int J Oral Sci, 4, №2, р. 101-108.
3. Mang, S, Vipperman, JS (2013). Mutual validation of muffler performance evaluation methods using field measurements and finite element analysis. J Acoust Soc Am, 134, №5, р. 4217.
4. Li, SY, Xu, BG, Tao, XM, Feng, J (2011). Numerical analysis of the mechanical behavior of a ring-spinning triangle using the Finite Element Method. Text Res J, 81, р. 959-971.
5. Peng Shang, Xin Ye, Xueling Bai, Zhentao, Hou, Zhaobin Xu, and Linan Zhang (2014). Hub bearing fatigue life prediction using finite element method. WIT Transactions on Information and Communication Technologies, 48, p. 1227-1233.
6. Bing, B, Pang, Y, Wang, B, Yang, A (2014). Finite Element Analysis of Strip Steel Jitter For Cooling Tower. Mechanical Research and Application, 1, p. 30-32.
7. Zheng-wei, LU (2014). Static and Dynamic Characteristic Finite Element Analysis of Slide Bearing Block Based on ANSYS. Mechanical Research and Application, 1, p. 33-35.
8. Jing, J (2014). Finite Element Analysis of Wheels. Mechanical Research and Application, 1, p. 102-104.
9. Kyoung-Sik Chun, Samuel Kinde Kassegne, Berhanu Kebede Wondimu (2009). Hybrid/mixed assumed stress element for anisotropic laminated elliptical and parabolic shells. Finite Elements in Analysis and Design. 45, p. 766-781.
10. Yakupov, NM, Serazutdinov, MN (1993). Calculation of Elastic Thin-Walled Structures of Complex Geometry. Kazan: CP RAS, 206 p.
11. Matvienko, YG, Chernyatin, AS, Razumovsky, IA (2013). Numerical analysis of three- dimensional non-singular components of the stress field at the crack tip of mixed type. Problems of Mechanical Engineering and Machine Reliability, №4, p. 40 - 48.
12. Bath, KY, Shydlouski, VP (trans. from English.); Turchak, LI (eds.) (2010). Finite Element Methods. M.: Fizmatlit, 1022 p.
13. Golovanov, AI, Tuleneva, ON, Shigabutdinov, AF (2006). The Finite Element Method in Statics and Dynamics of Thin-Walled Structures. M.: Fizmatlit, 391 p.
14. Kositsyn, SB, Chan, SL (2013). Analysis of the stress-strain state of intersecting cylindrical shells under elastic-plastic deformations considering geometric nonlinearity. Structural Mechanics of Engineering Constructions and Buildings, №1. p. 3 - 9.
15. Agapov, VP (2000). The Finite Element Method in Statics, Dynamics and Spatial Stability of Thin-Walled Reinforced Structures. M.: Publishing DIA, 152p.
16. Nikolaev, AP, Klotchkov, YV, Kiselev, AP, Gureeva, NA (2012). Vector Interpolation Displacement Fields in Finite Element Calculations of Shells. Volgograd: VPO VGAU "Niva", 264 p.
17. Novozhilov, VV (1962). The Theory of Thin Shells. L.: Shipbuilding, 431p.
COMPARATIVE ANALYSIS OF SCALAR AND VECTOR FORMS OF APPROXIMATIONS IN A FEM AFTER THE V.V. NOVOZHILOV'S RELATIONS FOR ELLIPTIC CYLINDER
KLOCHKOV YuV, NIKOLAEV A.P., ISHCHANOV T.R.
A comparative analysis of scalar and vector fields displacement approximations are set out on the example of finite element calculation of an elliptic cylinder with using of the Novozhilov theory of thin shells A curvilinear quadrilateral finite element as a fragment of the middle surface of an elliptical cylinder with eighteen degrees of freedom at the node was used. It was shown on numerical examples that the vector approximation being invariant has fundamental advantages over scalar approximation in the calculation of arbitrary shells.
KEY WORDS: vector approximation, scalar approximation, finite element, an elliptical cylinder.
