CC BY
35
УДК 539.3
Вестник СибГАУ Том 17, № 3. С. 587-594
РАСЧЕТ КОМПОЗИТНЫХ ЦИЛИНДРИЧЕСКИХ ОБОЛОЧЕК С ПРИМЕНЕНИЕМ МНОГОСЕТОЧНЫХ ЭЛЕМЕНТОВ
А. Д. Матвеев1 , А. Н. Гришанов2
1 Институт вычислительного моделирования СО РАН Российская Федерация, 660036, г. Красноярск, Академгородок, 50/44
2Новосибирский государственный технический университет Российская Федерация, 630073, г. Новосибирск, просп. К. Маркса, 20 *E-mail: mtv@icm.krasn.ru
Предложена процедура расчета трехмерных упругих композитных цилиндрических оболочек с различными коэффициентами наполнения, которая сводится к построению дискретных моделей, состоящих из криволинейных сложных многосеточных конечных элементов. В основе построения таких элементов лежат криволинейные двухсеточные конечные элементы. Двухсеточные и сложные многосеточные элементы проектируются на основе базовых конечно-элементных моделей композитных оболочек, которые учитывают их неоднородную структуру и имеют высокую размерность.
Показаны процедуры построения в локальных декартовых системах координат криволинейных двухсеточ-ных и сложных многосеточных элементов. Поля перемещений аппроксимируются известными степенными полиномами различных порядков, напряженное состояние описывается уравнениями трехмерной задачи теории упругости (без введения упрощающих гипотез о характере распределения полей перемещений, деформаций и напряжении). Аппроксимирующие полиномы и уравнения трехмерной задачи упругости записываются в локальных декартовых системах координат.
Достоинства предлагаемых элементов состоят в том, что они описывают трехмерное напряженное состояние в композитных оболочках, учитывают их неоднородные структуры, сложное закрепление и порождают многосеточные дискретные модели с малым числом узловых неизвестных. Размерности многосеточных дискретных моделей оболочек на несколько порядков меньше размерностей базовых моделей. Временные затраты реализации метода конечных элементов (МКЭ) на ЭВМ для многосеточных дискретных моделей композитных оболочек существенно меньше, чем для базовых моделей.
Предложен сложный многосеточный элемент 3-го порядка для расчета композитных цилиндрических оболочек. Приведен пример расчета по МКЭ консольной трехслойной оболочки с использованием сложных многосеточных элементов 3-го порядка. Результаты расчетов оболочки показывают высокую эффективность применения предложенных сложных элементов.
Ключевые слова: композиты, упругость, цилиндрические оболочки, метод конечных элементов, криволинейные сложные многосеточные конечные элементы
Sibirskii Gosudarstvennyi Aerokosmicheskii Universitet imeni Akademika M. F. Reshetneva. Vestnik Vol. 17, No. 3, P. 587-594
CALCULATION OF COMPOSITE CYLINDRICAL SHELLS USING MULTIGRID ELEMENTS
А. D. Matveev1*, A. N. Grishanov2
institute of Computational Modeling SB RAS 50/44, Akademgorodok, Krasnoyarsk, 660036, Russian Federation 2Novosibirsk State Technical University 20, K. Marksa Ау., Novosibirsk, 630073, Russian Federation E-mail: mtv@icm.krasn.ru
Calculating the three-dimensional elastic composite cylindrical shells with different coefficients of fullness that is reduced to the construction of discrete models consisting of complex curvilinear multi-grid finite elements has been proposed. The basis of such elements construction is curvilinear double-grid finite elements. Double-grid and complex multi-grid elements are designed based on the basic finite element models of composite shells which take into account their heterogeneous structure and have high dimension.
Constructing the curvilinear double-grid and complex multi-grid elements in the local Cartesian reference systems has been shown. Displacement fields are interpolated by known degree polynomials of various orders, the stress state is
described by the three-dimensional elasticity problem (without introduction of the simplifying hypotheses on the nature of the displacement fields, strain and stress distribution). Approximating polynomials and the equations of three-dimensional elasticity problem are recorded in the local Cartesian reference systems.
Advantages of the proposed elements are that they describe the three-dimensional stress state in composite shells, take into account their heterogeneous structures, complex fixing and generate multi-grid discrete models with a small number of nodal unknowns. The dimensions of multi-grid discrete shell models are by several orders less than the dimensions of the basic ones. Time spending of the finite element method (FEM) realization on a computer for discrete models of multi-grid composite shells is significantly less than for the basic models.
A complex multi-grid element of the 3rd order to calculate the composite cylindrical shells is proposed. The FEM calculation example of cantilever sandwich shell using complex multi-grid elements of the 3rd order has been given. Shell calculation results show the high efficiency of the proposed complex elements.
Keywords: composites, elasticity, cylindrical shell, finite element method, complex multigrid curvilinear finite elements.
Введение. Цилиндрические композитные (многослойные) оболочки широко применяются на практике, особенно в авиационной и космической технике. При расчетах композитных оболочек на прочность используют различные теории [1-4]. В основе существующих теорий оболочек лежат гипотезы, которые накладывают определенные ограничения на поля перемещений, деформаций и напряжений, существенно искажают реальные распределения напряжений в окрестностях закреплений оболочек или при действии локальных нагружений и порождают решения с неустранимой погрешностью. Чтобы избежать этих трудностей, необходимо решать проблему расчета оболочек в рамках трехмерной задачи теории упругости.
Как известно [5], расчет на прочность конструкции сводится к нахождению ее коэффициента запаса, который определяется с помощью максимального эквивалентного напряжения. Использование метода конечных элементов (МКЭ) при расчетах на прочность трехмерных композитных оболочек связано с большими трудностями [6-8]. При нахождении максимальных эквивалентных напряжений с малой погрешностью необходимо использовать достаточно мелкие разбиения, которые учитывают неоднородную (микронеоднородную) структуру оболочек, сложный вид их крепления и нагружения. Это порождает дискретные модели композитных оболочек очень высокой размерности и, как следствие, большую вычислительную погрешность. Таким образом, возникает необходимость в разработке новых конечных элементов (КЭ) и алгоритмов расчетов, более эффективных, чем существующие, которые сочетают в себе достоинства МКЭ и просто реализуются в расчетах напряженного деформированного состояния (НДС) конструкций со сложной неоднородной структурой.
В данной работе показано многосеточное моделирование трехмерных линейно-упругих композитных цилиндрических оболочек с различными коэффициентами наполнения. Предлагаемое моделирование сводится к построению дискретных моделей оболочек, которые состоят из криволинейных сложных многосеточных конечных элементов (МнКЭ). При проектировании сложного МнКЭ используется конечное число двухсеточных конечных элементов (ДвКЭ), которые представляют область сложного МнКЭ. Показаны процедуры построения в локальных декартовых системах координат криволинейных
ДвКЭ и сложных МнКЭ. Двухсеточные и сложные многосеточные элементы формы прямоугольного параллелепипеда (прямоугольника, треугольника) для анализа деформирования трехмерных (двумерных) упругих тел неоднородной структуры рассмотрены в работах [9-11]. В ДвКЭ и сложных МнКЭ функции перемещений интерполируются многочленами в форме степенных полиномов различных порядков [6-8], напряженное деформированное состояние описывается уравнениями трехмерной задачи теории упругости [12], которые записываются в локальных декартовых системах координат. В процедурах построения криволинейных ДвКЭ и сложных МнКЭ не требуется определять функции, которые связывают декартовые координаты с локальными криволинейными, и вычислять матрицу Якоби, что необходимо выполнять в существующих процедурах построения криволинейных КЭ [6-8].
При построении сложного МнКЭ используется криволинейная мелкая сетка (высокой размерности) и крупные сетки (малой размерности), вложенные в мелкую, на которых определяются аппроксимирующие функции перемещений. Базовое разбиение сложного элемента, которое учитывает его неоднородную структуру, порождает мелкую сетку. На базовом разбиении строится матричное выражение функционала полной потенциальной энергии сложного МнКЭ, в котором узловые перемещения мелкой сетки представляются через узловые перемещения крупных сеток. В результате полученное матричное выражение функционала энергии имеет малую размерность. Минимизируя построенный функционал, получаем формулы вычисления матрицы жесткости и вектора узловых сил сложного МнКЭ.
Приведен пример расчета трехслойной цилиндрической оболочки с применением сложных МнКЭ 3-го порядка. При построении ДвКЭ используются криволинейные однородные односеточные КЭ 1-го порядка, построенные в локальных декартовых системах координат.
1. Однородные криволинейные односеточные КЭ
Кратко рассмотрим процедуру построения криволинейных КЭ 1-го порядка в локальных декартовых системах координат (более подробно изложена в работах [13-15]). На рис. 1 показан криволинейный
односеточный КЭ У у 1-го порядка с характерными размерами Иу х Иеу х Иу, где ау - угол раствора КЭ У у; 01 х1 у12Х - локальная декартова система координат; 2101 У - плоскость симметрии; с! - ось цилиндрической оболочки; Я] , Щу - радиусы нижней и верхней поверхностей конечного элемента У у; Иу - толщина, Неу - длина КЭ ¥е, Иу = ауЩ ; узлы отмечены точками и пронумерованы (8 узлов). Итак, форма КЭ Уу есть прямая призма высотой Иеу, основанием которой
является криволинейный четырехугольник 1-3-4-2. Конечный элемент Уу рассматриваем как тело, которое испытывает трехмерное напряженное состояние и для которого перемещения, деформации и напряжения удовлетворяют соотношениям Коши и закону Гука трехмерной задачи теории упругости [12]. Для мелких разбиений угол раствора ау КЭ Уу мал. Тогда формы КЭ Уу и прямоугольного параллелепипеда мало отличаются. Поэтому для аппроксимации функций перемещений иу, V у, м>у элемента Уу используем полиномы 1-го порядка, записанные в локальной декартовой системе координат 01 х1 у12Х (рис. 1). Функцию перемещений иу (V у, м>у) представим в виде
иу = а1 + а2 х1 + а3 у1 + а4 х1 + а5 х У\ + + а6 2Х Х1 + а7 ^ у + а8 х У1
Используя представления вида (1), по алгоритмам МКЭ строим аппроксимирующие функции перемещений иу, V у, м>у для КЭ Уу, которые применяем в виде [7; 8]
иу =£ Ыги), V у = £ М1у1, Wу = £ ММ, (2)
1=1 1=1 1=1
где и], у) , м1 - перемещения 1-го узла; и) (у1, м1) -перемещение 1-го узла КЭ Уу в направлении оси 01 х1 (01 у1, 012Х); N = N (х1, у1, г1) - функция формы 1-го узла КЭ Уу.
Рис. 1. КЭ Уу 1-го порядка
Полную потенциальную энергию Wу (ву) для КЭ Уу запишем в матричной форме [7; 8]:
Wу (ву) = -2| (ву )т [ Ву ]т [ Бу ] [ Ву ]в)СУ -
2 Уу
-1 (в))т [му]т ¥уёУ - | (в))т [N у ]т qydS, (3)
Уу Яу
где [Ву ], [Бу ] - матрицы деформаций и модулей упругости КЭ У у; в] - вектор узловых перемещений КЭ У у, имеющий (в системе координат 01 х1 у1 г1) следующую структуру:
=1 г 1 11 11 1 ,т
где N - число узлов КЭ Уу (для рис. 1 имеем N = 8).
Применяя принцип минимума потенциальной энергии [8] для КЭ У у, т. е. выполняя условие
дWу (в]) / д6] = 0, получаем соотношение [К] ] в] = Р], где
[К]] = | [Ву]т[Бу][Ву]!У ,
Уу
Ру1 = | [N у ]т ¥у с!У +| [N у ]т qуСЯ , (5)
Уу Яу
где [К]] - матрица жесткости; Р] - вектор узловых сил КЭ У , которые определяются в локальной декартовой системе координат 01 х1 у1 хх КЭ У у.
Учитывая (5), выражение полной потенциальной энергии W у (в]) представим в виде
Wу (в))=2(в])т [ к)] в) - (в])т р]. (6)
Криволинейность области КЭ Уу учитывается в формулах (5). Далее выражения полных потенциальных энергий криволинейных композитных ДвКЭ и сложных МнКЭ для удобства изложения будем записывать в форме, аналогичной (6). Процедуры построения однородных криволинейных односеточных КЭ 2-го, 3-го порядка, которые геометрически подобны форме КЭ У (рис. 1), аналогичны вышеописанной. Здесь рассматриваем композитные криволинейные ДвКЭ и сложные МнКЭ, геометрически подобные форме КЭ У (рис. 1).
2. Композитные криволинейные ДвКЭ
Рассмотрим две процедуры построения трехмерных криволинейных ДвКЭ с неоднородной структурой на примере ДвКЭ Уа (Уь) 3-го порядка с характерными размерами Иу х Иу х Иу (рис. 2, где 0хух -локальная декартова система координат; Иу - толщина и Иуа - длина ДвКЭ Уа ).
Пусть между компонентами неоднородной структуры ДвКЭ Va связи идеальны. Компоненты есть линейно-упругие однородные тела, перемещения, деформации и напряжения которых удовлетворяют соотношениям Коши и закону Гука трехмерной задачи теории упругости [12]. Базовая дискретная модель (базовое разбиение) Ra ДвКЭ Va состоит из однородных криволинейных КЭ Ve 1-го порядка с характерными размерами hex х hey х hi (см. рис. 1). Дискретная модель Ra учитывает неоднородную структуру ДвКЭ Va и порождает мелкую сетку ha. Отметим, что дискретные модели Ra всех ДвКЭ Va образуют базовую дискретную модель оболочки. Считаем, что hey, hi, ae = const, где e = 1, ..., M; M - общее число КЭ Ve (см. рис. 1). На мелкой сетке ha определяем крупную сетку Ha ДвКЭ Va. Узлы сетки Ha (32 узла), отмеченные на рис. 2 точками, равномерно расположены на ребрах ДвКЭ Va (по 4 узла на каждом ребре). Порядок ДвКЭ определяется порядком интерполяционного полинома, который построен на его крупной сетке. Считаем, что угол раствора aa ДвКЭ Va мал. Тогда криволинейная форма ДвКЭ Va мало отличается от формы прямоугольного параллелепипеда. Поэтому при построении функции перемещений на крупной сетке Ha ДвКЭ Va используем известные интерполяционные полиномы 1-го, 2-го и 3-го порядка [7; 8], записанные в локальной декартовой системе координат Oxyz. В связи с этим, ДвКЭ Va будем называть ДвКЭ соответственно 1-го, 2-го и 3-го порядка.
Рис. 2. ДвКЭ Уа (Уь )
2.1. Первая процедура построения ДвКЭ. Первую процедуру (не теряя общности суждений) рассмотрим на примере построения ДвКЭ Уа 3-го порядка (рис. 2). На крупной сетке На строим аппроксимирующие функции перемещений иа, уа, wa ДвКЭ ¥а, используя полиномы 3-го порядка, записанные в системе координат Охуг . Для функции иа имеем
иа = а1 + а2 х + а3 у + а4 г + а5 ху + а6 гх + а7 гу + а8 хуг +
+ а9 х2 + а10у2 + а11 г 2 + а12 х2 у + а13 х2 г + а14 у2 х +
2 2 2 + а15у г + а16х х + а17г у +
3 3 3 3 3 3
+ a18x + a19y + a20z + a21x y + a22xz + a23yx + + a24 y3 z + a25 z3 x + a26 z3 y + a27 x2 yz + a28 y2 xz + + a29 z 2 xy + a30 x3 yz + a31 y3 xz + a32 z3 xy. Функции перемещений ua, va, wa, построенные по МКЭ, запишем в форме [7]
"0 -О -О
Ua = Ё Nq , Va = Ё Nqv , wa = Ё Npqpw , (7)
z=1
"О
Ё
i=1
"О
Ё
i =1
где д'и, дУу, д^, - перемещения и функция формы Р-го узла сетки На, р = 1, ...,п0; п0 - общее число узлов сетки На (для ДвКЭ 3-го порядка п0 = 32).
Обозначим 6а = {д1, ...,g3u2,д1, ...,д^д?, ...,д32}Т -вектор узловых перемещений крупной сетки На, т. е. ДвКЭ Уа. Пусть ось О1 у1 КЭ Уе (см. рис. 1), который расположен в базовой модели ДвКЭ ¥а, параллельна оси Оу (см. рис. 2) и между осями О1 х1 и Ох угол
равен Ре. Вектор 6е узловых перемещений КЭ ¥е, отвечающий системе координат Охуг , имеет вид
eT
0( e e e e e e v
e = {u1 V1 v8 , w1 w8}
(8)
где иг- , у , м>е - перемещения узла / КЭ Уе, I = 1, ..., 8.
Между векторами , 6е КЭ Уе (см. п. 1) имеем связь = [Те ]6е, где [Те ] - матрица вращений (размерности 24 х 24), которая в силу (4), (8) имеет вид [8]
Щ] [Е0] [ Е2]
[Te ] =
[ Eo] [ E] [ Eo]
-[ E2] [ Eo] [ E1]
(9)
где [E0], ...,[E2] - матрицы размерности 8x8; [E0] -нулевая и [ E] - единичная матрицы; [EJ = cos Pe [£]; [E2] = sin Pe [E].
Подставляя = [Te ]5e в (6) и выполняя условие dWe (5e) / 55 e = 0 , получаем соотношение [Ke ]5e = Fe, где [ Ke ] - матрица жесткости и Pe - вектор узловых сил КЭ Ve, определяемые в декартовой системе координат Oxyz (рис. 2) по формулам
[Ke ] = [Te ]T [K]][Te ], Pe = [Te í Pj . (10)
Выражение полной потенциальной энергии Wa для ДвКЭ Va запишем в матричной форме
M f 1 ч
Wa = £ [2 5T [Ke ] 5e - 5^Pe J.
(11)
Используя функции (7), узловые перемещения вектора 6е КЭ ¥е выражаем через узловые перемещения крупной сетки На, т. е. вектора 6а ДвКЭ Уа. В результате получим равенство
8е = [4й ]8в
(12)
где [Лае ] - прямоугольная матрица размерностью 24 х 96, элементы которой (в силу (7)) определяются значениями функций Жр в узлах КЭ Уе, е = 1, ...,М.
Подставляя (12) в (11) и выполняя условие дШа (8а)/ дЪа = 0, получаем матричное соотношение [Ка ]8а = Ра, где [Ка ] - матрица жесткости; Ра - вектор узловых сил ДвКЭ Уа, которые определяем по
М М
формулам [Ка ] = X [А-е ] [Ке ] [А ], ¥а = £ [А ]Т Ре .
е=1 е=1
2.1.1. Процедура определения напряжений в ДвКЭ Уа. Пусть найден вектор 8а узловых перемещений ДвКЭ Уа в локальной декартовой системе координат Охух (см. рис. 2). По формуле (12) находим (в декартовой системе координат Охух) вектор 8е узловых перемещений КЭ Уе, е = 1, ...,М . Определяем вектор 8^ узловых перемещений КЭ Уе в локальной декартовой системе координат О1 х1 у12Х (см. рис. 1) по формуле 8:е = [Те]8е, е = 1, ...,М. Используя вектор , функции перемещений КЭ Уе, по алгоритмам МКЭ находим функции напряжений в каждом КЭ Уе базового разбиения ДвКЭ Уа. Это позволяет определять напряжения в любом компоненте неоднородной структуры ДвКЭ Уа.
2.2. Вторая процедура построения ДвКЭ. Основные положения 2-й процедуры рассмотрим на примере построения криволинейного ДвКЭ Уь 3-го порядка, который имеет размеры, геометрическую форму, неоднородную структуру и расположен в локальной декартовой системе координат Охух, как
ДвКЭ Уа (рис. 3). В данной процедуре используем мелкую На и крупную На сетки, базовое разбиение Яа и функции перемещений иа, уа, wa (см. формулу (9)) ДвКЭ Уа, т. е. 8а есть вектор узловых перемещений ДвКЭ Уь. На базовом разбиении Яа ДвКЭ Уь с помощью метода конденсации [8] строим суперэлемент О/, полную потенциальную энергию IV, которого запишем в матричной форме:
W = |(86')T[KS] 8b -(8b)TPb,
(13)
где [К, ], Р/ , 8Ь - матрица жесткости, векторы узловых сил и перемещений суперэлемента О/, определяемые в декартовой системе координат Охух (рис. 3).
Используя функции (7), узловые перемещения суперэлемента О, выражаем через узловые перемещения крупной сетки На ДвКЭ У , т. е. для векторов 8/, 8а построим связь вида
8 b = [ AS ]8а
(14)
где [ AS ] - прямоугольная матрица, элементы которой (в силу (7)) определяются значениями функций Np
в узлах суперэлемента GS .
Подставляя (14) в (13) и выполняя условие dWs (8a) / d8a = 0, получаем матричное соотношение [Kb ]8a = Fb, где [Kb ] - матрица жесткости; Fb - вектор узловых сил ДвКЭ Vb, которые определяем
по формулам [Kb ] = [4 ]T [K£ ][4 ], Fb = [4 ]T P/ .
Процедура определения напряжений в базовых КЭ ДвКЭ Vb аналогична п. 2.1.1. Как показывают расчеты, ДвКЭ Vb (построенные по 2-й процедуре) порождают более точные решения, чем ДвКЭ Va (построенные по 1-й процедуре). Однако 2-я процедура включает построение суперэлементов, что связано с обращением матриц высокого порядка. Это усложняет алгоритм и увеличивает временные затраты на реализацию 2-й процедуры на ЭВМ. Криволинейные ДвКЭ 1-го, 2-го порядка с неоднородной (микронеоднородной) структурой, строим по процедурам, аналогичным процедурам п. 2.
3. Композитные криволинейные сложные МнКЭ
Рассмотрим процедуру построения композитных криволинейных сложных МнКЭ (не теряя общности суждений) на примере сложного МнКЭ Vf 3-го порядка с характерными размерами hf х hf х hf,
который имеет неоднородную структуру и расположен в локальной декартовой системе координат Oxyz
(рис. 3). Область МнКЭ Vf представляем криволинейными ДвКЭ УП 3-го порядка, построенными по алгоритмам, аналогичным п. 2. Базовые дискретные модели ДвКЭ УП состоят из однородных криволинейных КЭ Уе 1-го порядка с характерными размерами he х hy х he (см. рис. 1) и учитывают неоднородную структуру ДвКЭ УП, т. е. сложного МнКЭ Vf. Пусть hey,hi,ae = const, где e = 1, ...,Ma, Ma - общее число КЭ Уе базовой дискретной модели МнКЭ Vf . ДвКЭ УП с характерными размерами bxП х b" х bnz расположен в локальной декартовой системе координат Oj xУ1 Zj, причем ось Oy МнКЭ Vf параллельна оси O1 y1 ДвКЭ V^b, где n = 1, ...,N, N - общее число ДвКЭ Vb . Для рис. 3 имеем N = 27, hf = 3bb , hf = 3bb - толщина, hf = 3^ - длина МнКЭ Vf . Пусть b"n, bb, у n = const, n = 1, ..., N.
На мелких сетках ДвКЭ Vb определяем крупную
сетку Hm МнКЭ Vf. Узлы сетки Hm равномерно
расположены на ребрах МнКЭ Vf , по 4 узла на каждом ребре (рис. 3), узлы отмечены точками (32 узла).
Отметим, что МнКЭ VМ включает конечное множество криволинейных различных мелких и крупных вложенных сеток ДвКЭ ^ и крупную сетку Hm.
При мелких разбиениях угол раствора ат МнКЭ VМ мал. Тогда формы МнКЭ VМ и прямоугольного параллелепипеда мало отличаются. В связи с этим при построении аппроксимирующих функций перемещений um, vm, wm для сложного МнКЭ УМ (рис. 3), определяемых по МКЭ на крупной сетке Hm, используем полиномы 3-го порядка [7; 8] и поэтому МнКЭ Vm называем сложным МнКЭ 3-го порядка. Функции um , vm, wm запишем в форме
32 32 32
ит = Ь Ni ^ , vm = Ь Ni V , wm = Ь Ni ^ , (15)
1 = 1 I = 1 ¿ = 1
где ^ - базисная функция ¿-го узла крупной сетки
Ят т т ^
т; и , vi , wi - значения перемещений в ¿-м узле сетки Нт, определяемые в декартовой системе координат Охух , I = 1, ..., 32 .
о х
Рис. 3. Сложный МнКЭ vm
Обозначим через qm вектор узловых перемещений (размерностью 96) крупной сетки Нт в декартовой системе координат Охух, т. е. qМ есть вектор узловых перемещений сложного МнКЭ УМ. Пусть ось О1 у1 параллельна оси Оу и между осями О1 х1 и Ох угол равен фП. Пусть векторы бП , бП узловых перемещений ДвКЭ VП, отвечающие соответственно системам координат О1 х1 у1 х1 и Охух , имеют структуру вида (4). Между векторами бП, бь установим
связь бП = [ТПЬ ]бьП, где [ТПЬ ] - матрица вращений, имеющая такую же структуру, как и матрица вращений [Те ] (см. формулу (9)). Учитывая, что
бП = [ТПь ]бЬ , получим соотношения [кЬ ] = [ТПь ]Т х
х [КП][Ть], р = [Тпь]Р,1, где [Кп], [Кьп] - матрицы
жесткости; РП , РП - векторы узловых сил ДвКЭ VП ,
отвечающие соответственно декартовым системам координат O1x1 y1 z1 и Oxyz , n = 1, ..., N. Полную
потенциальную энергию сложного МнКЭ V^ представляем матричным выражением
wm = £ )T K ] 6b - (6b )T p). (16)
n=1 2
Используя функции (15), узловые перемещения
с-Ь
вектора 6n выражаем через узловые перемещения крупной сетки Hm сложного МнКЭ vm, т. е. между векторами 6b , qm установим связь вида
6bn=[ Am ]qm, (17)
где [Am ] - квадратная матрица размерности, элементы которой (в силу (15)) определяются значениями функций N™ в узлах ДвКЭ Vt, n = 1, ..., N.
Подставляя (17) в (16) и выполняя условие dwm (qm) / dqm = 0 , получим матричное соотношение [Km ]qm = Fm, где [Km ], Fm - матрица жесткости и вектор узловых сил сложного МнКЭ Vm, которые
N
определяем по формулам [Km ] = £А í [Kbn ][Anm ],
n=1
N
F^ = £[Am]TPt . Криволинейные сложные МнКЭ
n=1
1-го и 2-го порядка с неоднородной (микронеоднородной) структурой, которые геометрически подобны форме сложного МнКЭ vm, строим по процедуре, аналогичной процедуре п. 3.
3.1. Процедура определения напряжений в сложном МнКЭ Уат. Пусть найден вектор q^ узловых
перемещений МнКЭ vm , отвечающий системе координат Oxyz (рис. 3). По формуле (17) находим вектор 6b узловых перемещений ДвКЭ V¡b, n = 1, ..., N. С помощью матрицы вращений [Tnb ] вектор 6n узловых перемещений ДвКЭ V¡b определяем (в системе координат Oj xj y Zj ) по формуле 6^ = [Tnb ]6^ . Используя вектор 6n, напряжения в КЭ Ve базового
разбиения ДвКЭ V¡b вычисляем по процедуре, аналогичной процедуре п. 2.1.1.
4. Результаты расчетов
В качестве модельной задачи рассмотрим расчет в безразмерных параметрах цилиндрической трехслойной оболочки V0 постоянной толщины h, которая расположена в декартовой системе координат Ox0y0z0 (рис. 4). При y0 = 0 имеем u = v = w = 0, т. е. левый торец оболочки жестко закреплен. Жесткие (верхний и нижний) слои составляют 1/6 толщины оболочки и параллельны оси Oy0. Модуль Юнга
связующего материала равен 1, жестких слоев - 10, коэффициент Пуассона для жестких слоев и связующего материала равен 0,3. Радиус внутренней поверхности оболочки равен 140, радиус внешней поверхности - 150, толщина оболочки h = 10, длина L = 350 . На верхней поверхности S0 оболочки, симметричной относительно плоскости Oy0z0, приложена равномерно распределенная радиальная нагрузка q = 0,01. Область S0 на рис. 4 заштрихована, угол раствора а области S0 равен % /2. В расчетах, силу симметричного вида нагружения, используем половину оболочки относительно плоскости Oy0z0. Базовая дискретная модель R0 оболочки V0 состоит из однородных криволинейных КЭ Ve 1-го порядка с характерными размерами hex х hey х hZ, (см. рис. 1). Базовая модель R0 оболочки V0 учитывает ее неоднородную структуру и порождает мелкую сетку h0 размерностью 145 х109 х 7, h = 6hZ . Имеем hey, h,, ае = const, где е = 1, ..., M ; M - общее число КЭ Ve.
Рис. 4. Схема нагружения оболочки
Рис. 5. ДвКЭ Vd
Многосеточная дискретная модель Rh оболочки V0 состоит из криволинейных сложных МнКЭ V^ 3-го порядка с характерными размерами \Шех х1Шу х6Ие2, мелкие сетки которых имеют размерность 19 х19 х 7 . Сложные МнКЭ V^ построены по процедуре, анало-
гичной процедуре п. 3, а = 1, ..., N0, N0 - общее число МнКЭ V^ , N0 = 48.
При построении сложных МнКЭ Vm используем девять криволинейных ДвКЭ VП 3-го порядка с характерными размерами 6h, х 6hey х 6h, (рис. 5), которые построены по процедуре п. 2.2 и мелкие сетки которых имеют размерность 7 х 7 х 7 . Результаты расчетов оболочки V0 показывают, что максимальное эквивалентное напряжение am = 2,206 (вертикальное перемещение wm = 156,188) многосеточной модели Rm отличается от максимального эквивалентного напряжения ст0 = 2,270 (перемещения w0 = 158,153) базовой модели R0 оболочки на 2,81 % (на 1,24 %). Размерность базовой модели R0 равна 327348, ширина ленты системы линейных алгебраических уравнений (СЛАУ) МКЭ - 2295. Многосеточная модель R m
оболочки V0 имеет 1776 узловых неизвестных, ширина ленты СЛАУ МКЭ равна 252. Реализация МКЭ для многосеточной модели Rm требует в 1678 раз меньше объема памяти ЭВМ, чем для базовой модели R0. Построение решения для многосеточной дискретной модели Rm требует в 420 раз меньше временных затрат, чем для базовой модели R0 . Эквивалентные напряжения определяем по 4-й теории прочности в центрах тяжести КЭ Ve базовой модели R0 оболочки.
Заключение. В данной работе показаны процедуры построения композитных криволинейных сложных МнКЭ 1-го, 2-го и 3-го порядка. Достоинства предлагаемых конечных элементов состоят в том, что они:
- описывают трехмерное напряженное деформированное состояние в композитных цилиндрических оболочках;
- учитывают неоднородную (микронеоднородную) структуру оболочек;
- порождают многосеточные дискретные модели оболочек, размерности которых на несколько порядков меньше размерностей базовых дискретных моделей;
- учитывают сложный характер крепления и на-гружения оболочек.
Реализация МКЭ для многосеточных дискретных моделей композитных цилиндрических оболочек требует мало временных затрат и значительно меньше объема памяти ЭВМ, чем для базовых моделей.
Благодарности. Работа выполнена при финансовой поддержке РФФИ (код проекта 14-01-00130).
Acknowledgments. This work was financially supported by RFBR (code project 14-01-00130).
Библиографические ссылки
1. Болотин В. В., Новиков Ю. Н. Механика многослойных конструкций. М. : Машиностроение, 1980. 375 с.
2. Голушко С. К., Немировский Ю. В. Прямые и обратные задачи механики упругих композитных пластин и оболочек вращения. М. : ФИЗМАТЛИТ, 2008. 432 с.
3. Андреев А. Н., Немировский Ю. В. Многослойные анизотропные оболочки и пластины: изгиб, устойчивость, колебания. Новосибирск : Наука, 2001. 288 с.
4. Киселев А. П. Расчет тонких оболочек на прочность в трехмерной постановке без упрощающих гипотез // Изв. вузов. Сер. «Строительство». 2008. № 1. С. 18-23.
5. Москвичев В. В. Основы конструкционной прочности технических систем и инженерных сооружений. Ч. 1. Постановка задач и анализ предельных состояний. Новосибирск : Наука, 2002. 106 с.
6. Сегерлинд Л. Применение метода конечных элементов в технике. М. : Мир, 1975. 541 с.
7. Зенкевич О. Метод конечных элементов в технике. М. : Мир, 1975. 541 с.
8. Норри Д., де Фриз Ж. Введение в метод конечных элементов. М. : Мир, 1981. 304 с.
9. Матвеев А. Д. Некоторые подходы проектирования упругих многосеточных конечных элементов / Ин-т вычислительного моделирования СО РАН. Красноярск, 2000. 30 с. Деп. в ВИНИТИ № 2990-В00.
10. Матвеев А. Д. Многосеточное моделирование композитов нерегулярной структуры с малым коэффициентом наполнения // ПМТФ. 2004. № 3. С. 161-171.
11. Матвеев А. Д. Построение сложных многосеточных элементов с неоднородной и микронеоднородной структурой // Известия АлтГУ. Сер. «Математика и механика». 2014. № 1/1. С. 80-83.
12. Самуль В. И. Основы теории упругости и пластичности. М. : Высш. шк., 1982. 264 с.
13. Матвеев А. Д., Гришанов А. Н. Смешанные многосеточные дискретные модели трехмерных цилиндрических композитных панелей и оболочек сложной формы // Сб. ст. XIX зимней школы по механике сплошных сред. Пермь, 2015. С. 198-211.
14. Матвеев А. Д., Гришанов А. Н. Многосеточные криволинейные элементы в трехмерном анализе цилиндрических композитных панелей с полостями и отверстиями // Ученые записки Казанского университета. Сер. «Физико-математические науки». 2014. Т. 156, кн. 4. С. 47-59.
15. Матвеев А. Д., Гришанов А. Н. Расчет композитных цилиндрических панелей и оболочек с учетом их структуры на основе свыше миллиарда уравнений МКЭ с малыми временными затратами // Численные методы решения задач теории упругости и пластичности : материалы ХХ1У Всерос. конф. Омск : Изд-во ОмГТУ, 2015. С. 141-150.
References
1. Bolotin V. V., Novikov Yu. N. Mekhanika mnogosloynykh konstruktsiy. [Mechanics of multilayer structures]. Moscow, Mashinostroenie Publ., 1980, 375 p.
2. Golushko S. K., Nemirovskiy Yu. V. Pryamye i obratnye zadachi mekhaniki uprugikh kompozitnykh plastin i obolochek vrashcheniya. [Direct and inverse
problems of mechanics of elastic composite plates and shells of revolution]. Moscow, FIZMATLIT Publ., 2008, 432 p.
3. Andreev A. N., Nemirovskiy Yu. V. Mnogosloynye anizotropnye obolochki i plastiny: Izgib, ustoychivost', kolebaniya. [Multilayer anisotropic shells and plates: bend, resistance, hesitation], Novosibirsk, Nauka Publ., 2001, 288 p.
4. Kiselev A. P. [Calculation of thin shells for strength to resist three-dimensional setting, without simplifying hypotheses]. Izv. vuzov, ser. Stroitel'stvo. 2008, No. 1, P. 18-23 (In Russ.).
5. Moskvichev V. V. Osnovy konstruktsionnoy proch-nosti tekhnicheskikh sistem i inzhenernykh sooruzhenij. Ch. 1: Postanovka zadach i analiz predel'nykh sostoyaniy. [Basics of structural integrity of engineering systems and civil engineering. Part 1: Setting targets and analysis of limit states]. Novosibirsk, Nauka Publ., 2002, 106 p.
6. Segerlind L. Primenenie metoda konechnykh elementov v tekhnike. [Application of the finite element method in engineering]. Moscow, Mir Publ., 1975, 541 p.
7. Zenkevich O. Metod konechnykh elementov v tekh-nike. [Finite Element Method in engineering]. Moscow, Mir Publ., 1975, 541 p.
8. Norri D., Zh. de Friz. Vvedenie v metod konechnykh elementov. [An Introduction to Finite Element Analysis]. Moscow, Mir Publ., 1981, 304 p.
9. Matveev A. D. [Some approaches of designing elastic multigrid finite elements]. VINITI Proceedings no. 2990-B00. Krasnoyarsk, 2000, 30 p. (In Russ.).
10. Matveev A. D. [Multigrid modeling of composites of irregular structure with a small filling ratio] PMTF, 2004, No. 3, P. 161-171 (In Russ.).
11. Matveev A. D. [Construction of the complex multigrid elements of inhomogeneous and microinhomo-geneous structure]. Izvestiya AltGU, seriya: Matematika i mekhanika, 2014, No. 1/1, P. 80-83 (In Russ.).
12. Samul' V. I. Osnovy teorii uprugosti i plastich-nosti. [Fundamentals of the elasticity and plasticity theory]. Moscow, Vysshaya shkola Publ., 1982, 264 p.
13. Matveev A. D., Grishanov A. N. [Mixed multiugrid discrete models of three - dimensional cylindrical composite shells and panels of complex shape]. Sbornik statey XIX zimney shkoly po mekhanike sploshnykh sred. Perm'. 2015, P. 198-211 (In Russ.).
14. Matveev A. D., Grishanov A. N. [Multigrid curvilinear elements in three-dimensional analysis of cylindrical composite panels with cavities and holes]. Uchenye zapiski Kazanskogo universiteta. Seriya: Fiziko-matematicheskie nauki. 2014, Vol. 156, book 4, P. 47-59 (In Russ.).
15. Matveev A. D., Grishanov A. N. [Calculation of composite cylindrical panels and shells based on their structure on the basis of more than a billion FEM equations with small time consuming]. Chislennye metody resheniya zadach teorii uprugosti i plastichnosti. Materialy XXIV Vserossiyskoy konferentsii [Proceedings of the IV All-Russian Conference]. Omsk, Izd-vo OmGTU Publ., 2015, P.141-150 (In Russ.).
© Матвеев А. Д., Гришанов А. Н., 2016
