ЖУРНАЛ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ И МАТЕМАТИЧЕСКОЙ ФИЗИКИ, 2015, том 55, № 3, с. 418-428
УДК 519.63
СРАВНЕНИЕ СКАЛЯРНОЙ И ВЕКТОРНОЙ ФОРМ МЕТОДА КОНЕЧНЫХ ЭЛЕМЕНТОВ НА ПРИМЕРЕ ЭЛЛИПТИЧЕСКОГО ЦИЛИНДРА
© 2015 г. Т. А. Киселева, Ю. В. Клочков, А. П. Николаев
(400002 Волгоград, Волгоградская о., пр-т Университетский, 26, ВолГАУ) e-mail: Klotchkov@bk.ru Поступила в редакцию 19.05.2014 г.
Переработанный вариант 28.08.2014 г.
Предлагается и реализуется инвариантная аппроксимация искомых величин в векторной формулировке при формировании матрицы жесткости четырехугольного криволинейного конечного элемента в виде фрагмента срединной поверхности эллиптического цилиндра с восемнадцатью степенями свободы в узле. На численных примерах показано, что векторная аппроксимация обладает принципиальными преимуществами по сравнению со скалярной аппроксимацией при расчете произвольных оболочек со значительными градиентами кривизн линий срединной поверхности. Библ. 14. Фиг. 2. Табл. 3.
Ключевые слова: векторная аппроксимация, скалярная аппроксимация, конечный элемент, эллиптический цилиндр, теория оболочек.
DOI: 10.7868/S0044466915030102
1. ВВЕДЕНИЕ
В настоящее время разработана достаточно подробная теория деформирования оболочек (см., например, [1]—[6]). Однако аналитическое решение систем дифференциальных уравнений, описывающих процесс деформирования оболочек, в большинстве случаев не представляется возможным, поэтому широкое распространение получили численные методы расчета, среди которых наибольшую известность приобрел метод конечных элементов (МКЭ) (см. [7]—[11]).
В подавляющем большинстве современных конечно-элементных алгоритмов, в том числе и в вычислительных комплексах типа ANSYS, ABAQUS, NASTRAN и им подобных, реализована скалярная интерполяционная процедура, суть которой заключается в аппроксимации отдельной компоненты вектора перемещения через узловые значения этой же компоненты. При достаточно плавном изменении геометрии оболочки и отсутствии смещений как жесткого тела скалярная аппроксимация позволяет получать удовлетворительные результаты. Однако при анализе напряженно-деформируемого состояния (НДС) оболочек произвольной формы со значительными значениями кривизн срединной поверхности или при наличии жестких смещений скалярная аппроксимация приводит к существенному росту погрешностей конечно-элементных решений. В этих случаях необходимо применять векторную аппроксимацию, суть которой заключается в использовании аппроксимирующего выражения непосредственно для вектора перемещения произвольной точки конечного элемента через векторы перемещений узлов и их производные.
2. ОПИСАНИЕ МЕТОДА
При компоновке матриц жесткостей криволинейных конечных элементов обычно используется скалярная (неинвариантная) аппроксимация искомых величин (см. [7]—[10]). В [11] для расчета оболочек вращения на основе МКЭ разработана эффективная инвариантная аппроксимация искомых величин в векторной формулировке. В настоящей работе способ аппроксимации из [11] развит и реализован при компоновке матрицы жесткости и столбца внешней нагрузки конечного элемента в виде фрагмента срединной поверхности эллиптического цилиндра с восемнадцатью степенями свободы в узле. Реализованы два варианта аппроксимирующей процедуры:
общепринятая скалярная (см. [7]—[10]) и разработанная векторная аппроксимация искомых величин.
На примере расчета эллиптического цилиндра показана высокая эффективность векторной аппроксимации искомых величин.
2.1. Геометрия эллиптического цилиндра
Радиус-вектор точки срединной поверхности эллиптического цилиндра может быть задан в декартовой системе координат в виде
R0 = xi + yj + с^ 1 - y2/b 2k, (1)
2 2 y z
где b и С — полуоси поперечного сечения эллиптического цилиндра (-- + ^ = 1).
b с
Ковариантные векторы базиса точки срединной поверхности эллиптического цилиндра определяются дифференцированием функции (1) по глобальным координатам x, y
a0 = i, a0 = j - ■ СУ k, (2) b Vi - y /b
где должно выполняться условие
1 - y2/b2 > 0, (3)
что накладывает существенные ограничения на размеры рассчитываемого цилиндра, а именно, исключения из области определения линий пересечения цилиндра с плоскостью Oxy и зон непосредственно примыкающих к данным линиям.
С целью устранения указанных выше ограничений радиус-вектор эллиптического цилиндра предлагается записывать в виде
R0 = xi + r (0) sin 0j + r (0) cos 0k, (4)
где 8 — угол, отсчитываемый от оси Oz против хода часовой стрелки в плоскости, перпендикулярной оси Ox;
r (0) = bc/4 с2 sin2 0 + b2 cos2 0. (5)
Использование (4), (5) позволяет выполнять сплошную параметризацию срединной поверхности и осуществлять расчет эллиптического цилиндра как замкнутой оболочки без каких-либо ограничений на форму и геометрические размеры.
Тангенциальные векторы базиса и орт нормали к срединной поверхности эллиптического цилиндра определяются по формулам
0 Т>0 0т,0 0 0 0 , 0
a! = R = i, a2 = Ra, a = ajXa2/Va , (6)
0 0 0 I 0 \2 где a = ana- - (au)
Положение точки, отстоящей от срединной поверхности на расстоянии в исходном и деформированном состояниях определяется радиус-векторами
И0С = И0 +С а0, Ис = И0С + V. (7)
Входящий в (7) вектор перемещения точки, отстоящей от срединной поверхности на расстоянии с учетом гипотезы прямой нормали может быть представлен в виде
V = V + С (а - а0), (8)
где а — орт нормали в деформированном состоянии, V = vaaa + va — вектор перемещения точки срединной поверхности эллиптического цилиндра (здесь и ниже греческие буквы а, р принимают значения 1, 2).
Производные векторов локального базиса и вектора перемещения произвольной точки срединной поверхности оболочки могут быть представлены компонентами, отнесенными к векторам локального базиса оболочки в исходном состоянии:
0 _ -г0^ 0 . д0 0 0 _ , 0у 0
аа,в - Г авау + Ьава , а,а - -Ьа ау,
_,у 0 . 0 _,у 0 . 0
^а - ' аау + 'а , ¥ ,ав - ' авау + ' ава ,
(9)
где Г Ор — символы Кристоффеля II рода; — ковариантные компоненты тензора кривизны; 1Уа, ¿а, ^р и ?ар — многочлены, содержащие тангенциальные V и нормальную V компоненты вектора перемещения, а также их первые (для г^, ¿а) и вторые (для г^р, ¿ар) производные по глобальным криволинейным координатам.
Ковариантные компоненты тензора деформаций в слое, отстоящем от срединной поверхности на расстоянии определяются соотношением механики сплошной среды (см. [12])
эа|3 5 ав
) /2 •
(10)
Входящие в (10) ковариантные компоненты метрического тензора в исходном g а р и деформированном g ар состояниях определяются скалярными произведениями соответствующих базисных векторов
0 _ 0 0 >ар gа * gр,
х|3 - g а * g в ,
(11)
где gа=R0а, gа=R,
Соотношение (10) можно представить в виде суммы
Л
(12)
еар - еар + СХар ,
где Бар = (аа • V,в + у,а • ав)/2,
Ха° = (аа • а,° + V,а • а,в + а,а ' ав + а,а ' V,° + а 1 ' Ь^а[ + Ь^а• а°) /2.
Входящие в (12) производные орта нормали в деформированном состоянии определяются соотношением
1,а = (а! X а2/4а) ^ = (а1а х а2 + а! х а2,а)/4а + (а! х а2)(1/4а) а,
(13)
где аa = а^ + ^а0 + гаа0, а = апа22 — (а12)2 — детерминант метрического тензора в деформирован-
а ма 1 •'а^у
ном состоянии.
После определения с учетом (9) и (13) входящих в (12) скалярных произведений, деформации и искривления срединной поверхности произвольной оболочки могут быть представлены в развернутом виде:
_ 0 1 0 2 _ _ / 0 1 0 2 0 1 0 п _ 0 1 0 2 6ц — ац/1 + а12г1 , 812 — 821 — (ац/ 2 + а12г 2 + а12'1 + а22г1 )/А 6 22 — а12г2 + a22^2,
Х11 = г2
Г 02 01 т-<01 02 ,
11Й1 -ГПЙ1 +
0
01 02 а,1 02
2а
0"1
+ ^
Г 02 01 , т-<01 02 ,
21а1 +г 2 а +
0
. -,г01 г02 а,1 —0 а1 +1ц 21 —0
2а 2а
01 02 , а,1 01
2 0 01 1 0 02 0 2 0 01 0 1 0 02 + г2 (-Ь11а1 ) + г1 ( Ь21а1 - Ь11) + г1 (Ь21а1 - Ь12) + г2 (Ь11а1 )
, , / 01\ , , / 02\ , р / 0,01 , 0 ,02\ + г11 (-«1 ) + г21 (-«1 ) + ер (( + а12Ь1 ),
г
Xl2 = ¿2
f 0 О
т^02 01 т-01 02 a 2 02 ^02 01 ^01 02 a 1 02
Г12Й1 -Г^ + a1 +Гпа1 -Гпа2 + a2 v 2a 2a
01
+11
0
y^02 01 т-,01 02 , a, 2 01 y^02 01
-Г 22a1 +r22a1 a1 -r21a2 2a
0
y^01 02 , a, 1 01
Г 21a2 + TT a2
2a
01 12
0
+ Г 21
2a j
0
T-02 a ,2
22 —0 2a j
/2/2 +
2 0 01 0 01 0 1 0 02 0 02 0 + 12 (-¿^ - ¿10 - ¿12 )/2 +11 (-¿22a1 - ¿2O2 - ¿21) /2 +
2 0 01 0 01 0 1 0 02 0 02 0 + ¿1 ((1 + ¿2O2 - ¿22)/2 +¿2 (¿^ + ¿10 - ¿11 )/2 +
, , / 01 02 \ , , / 02 \ , , / 02\ , p/ 0 ,01 , 0 ,02 0 ,01 , 0 ,02\ + ¿12 (-0 - a2 )/2 + ¿22 (1 )/2 + ¿11 (-a2 )/2 + sp (an¿2 + au¿2 + au¿1 + a22¿1 ) /2,
(14)
X22 = ¿ 2
/ 0
т-,02 01 у^01 02 a,2 02 т-,01
r12a2 -r12a2 + a2 + Ги
v 2a
2
0 л
^02 0a. 22 0 2a y
/
+ ¿1
0 л
т-,02 01 т-,01 02 a,2 01
-Г 2202 +Г 2202 + 02 v 20 y
, .2/ ,0 01 ,0 \ , Jí ,0 02 \ , Л(,0 01\ , Л/, 0 02 ,0 \ + ¿2 (—¿12a2 - ¿22) + ¿1 (—¿22a2 ) + ¿1 (¿22a2 ) + ¿2 (¿12a2 - ¿12 )
, , / 01\ , , / 02\ , p( 0,01 , 0 ,02 \ + ¿12 (-a2 ) + ¿22 (-a2 ) + e P (a12¿2 + a12¿2 /
где 6Р
Sa|3a
0ав- a0a|\ a0
0a
р, ар — контравариантные, ковариантные и смешанные компоненты мет-
рического тензора исходного состояния.
Соотношение (12) можно представить в матричном виде:
{<£} = [Г] {s} = [r][D]{U},
(15)
где (е)т = (6П6222612X11X222X12}
3Х1 3x6 6x1 3x6 6x3 3x1
[D] — матрица алгебраических и дифференциальных операций
1x6
{и }т = {
1 2 . V v v}.
1x3
2.2. Физические соотношения тонких оболочек
Контравариантные компоненты тензора напряжений в произвольной точке оболочки выражаются через ковариантные компоненты тензора деформаций соотношениями (см. [12])
ав
а ' = ^ (е) ^ав + р*е РУ, (16)
где X, ц — параметры Ламе, gав — контравариантные компоненты метрического тензора,
А (е) = g еи + 2g е 12 + g е 22 + g е33 — первый инвариант тензора деформаций.
Выражения могут быть представлены в матричной форме:
{аав} = [С](е У, (17)
3x1
3x3
3x1
где
КТ={
11 22 12 a a a
; {sÜp} - {sns2226^2}.
2.3. Элемент дискретизации и аппроксимация искомъх величин
В качестве элемента дискретизации выбирается четырехугольный фрагмент срединной поверхности эллиптического цилиндра с узлами I, у, к, I, отображаемый для удобства численного интегрирования на квадрат с локальной системой коо
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.