Трехмерный конечный элемент для расчета оболочек средней толщины

Д. В. Бережной, М. К. Сагдатуллин

ТРЕХМЕРНЫЙ КОНЕЧНЫЙ ЭЛЕМЕНТ ДЛЯ РАСЧЕТА ОБОЛОЧЕК СРЕДНЕЙ ТОЛЩИНЫ

Ключевые слова: оболочечный конечный элемент, упругие деформации, метод двойной аппроксимации, гипотеза малости

напряжений обжатия.

Работа посвящена построению конечного элемента для расчета оболочек средней толщины на основе модификации трехмерного изопараметрического восьмиузлового элемента сплошной среды. Вводится гипотеза малости напряжений обжатия и используется техника понижения порядка аппроксимации по толщине элемента. Используется метод двойной аппроксимации по точкам суперсходимости. На численных примерах показана эффективность данного подхода.

Key words: shell finite element, elastic strains, method of double approximation, hypothesis of infinitesimality of compression stress.

The article regards the construction of a finite element for calculating the middle thickness of shells on the basis of a modification of three-dimensional isoparametric 8-node element through introduction of a hypothesis of infinitesimality of compression stress and usage of the technique of approximation order reduction. The method of double approximations on superconvergence points has been applied. Effectiveness of the given approach is shown on numerical examples.

Введение

Целью настоящей работы является разработка методики построения изопараметрического трехмерного конечного элемента (КЭ) сплошной среды, приспособленного для расчета оболочек средней толщины при однослойной аппроксимации по толщине. Обычно в подобных элементах узлы располагаются на срединной поверхности и в качестве степеней свободы кроме проекций вектора перемещений используют и углы поворота нормального волокна. Некоторый обзор таких КЭ приведен в [1-3]. Обычно эти элементы демонстрируют хорошую точность и достаточно эффективны в расчетах оболочек малой и средней толщины, однако их использование весьма затруднительно при моделировании сопряжений оболочек с массивными телами, так как необходимо выражать узловые перемещения трехмерных элементов через углы поворота оболочечных КЭ. В этом случае приходится использовать специальные процедуры (введение специальных множителей Лагранжа, использование метода штрафа и построение семейства специальных переходных конечных элементов) для стыковки трехмерных конечных элементов и элементов оболочек в единую расчетную модель или строить специальные трехмерные конечные элементы сплошной среды, включающие в качестве степеней свободы еще и углы поворота. Поэтому получили распространение специальные элементы оболочек, которые имеют в качестве узловых степеней свободы перемещения узлов, расположенных на лицевых поверхностях. Примеры подобных КЭ и их использования описаны в работах [3-6].

В предлагаемом конечном элементе вводятся изопараметрнческие аппроксимации геометрии и неизвестные компоненты перемещений. Узловые точки располагаются в вершинах элемента, в локальных координатах представляющего собой куб. По толщине используется линейная аппроксимация, что эквивалентно введению кинематической гипотезы

о линеином распределении перемещении по толщине с учетом обжатия оболочки.

Вариация потенциальной энергии деформации (работа внутренних напряжений на возможных деформациях) вычисляется численно по квадратурной формуле Гаусса-Лежандра (два узла по каждому координатному направлению). В каждом квадратурном узле вводится статическая гипотеза, обычно используемая для оболочек средней толщины, об отсутствии влияния поперечного напряжения на мембранные, изгибные и деформации поперечного сдвига. Фактически задача поперечного обжатия отделяется от классической задачи теории оболочек с учетом поперечного сдвига.

Главный источник ошибок, возникающих при использовании подобных КЭ (особенно для лагранжевой аппроксимации перемещений и геометрии в плане КЭ), которые делают их практически непригодными, - появление «ложных» деформаций поперечного сдвига. Для устранения этого недостатка предполагается использование метода двойной аппроксимации по точкам суперсходимости.

1. Трехмерный изопараметрический элемент оболочки средней толщины

В глобальной декартовой системе координат х, у, 2 предлагаемый КЭ представляет

собой искривленный параллелепипед, у которого верхняя и нижняя поверхности существенно искривлены, а четыре боковые грани являются линейчатыми поверхностями (рис. 1).

Для задания функций формы такого элемента вводится локальная система координат

7, С (в которой представленный объемный элемент будет прямолинейным кубом с гранями -1 <%, 7, С< 1). В общем случае связь между координатами х, у, 2 и £, 7, С имеет вид

Г (^((х І + ((ЄРрх 1 + гвеРхк )1^ +

І = 1 2

+^І

І=1 2

N 1 + — N 1 _ —

= У N.. (#, ц)гв(рх —-+у N.. (#, ц)>иж —-

Рис. 1

Характерной особенностью этого элемента является линейная аппроксимация по толщине оболочки. Функции формы N. (#,ц) в плане

конечного элемента, приведенного на рис. 1,

определяют стандартную квадратичную

аппроксимацию лагранжева типа.

Изопараметрическая концепция предполагает определение вектора перемещений в виде,

аналогичном представлению геометрии, т. е. в виде

N 1 +— N 1 _-

V = У N. (#,ц))>('}--------+У N. (#,ц))>('} --,

1 ’ / / верх ^ 1 ’ ' / нижи ГК ’

І =1 2 1 = 1 2

откуда видно, что искомые функции суть значения

декартовых проекций точек, лежащих на лицевых

поверхностях оболочки. Матрица жесткости строится

следующим образом. Первоначально определяется

матрица Якоби преобразования координат

д дх ду дх ' д ' д '

д# д# д# д# дх дх

д дх ду дх д = [з ] д

дц дц дц дц ду ду

д дх ду дг д д

[д- _д- д- И-\ дг. дг,

Компоненты этой матрицы вычисляются прямым дифференцированием выражения (1), при этом видно, что они будут зависеть от координат текущей точки #, ц, -, т.е. их можно вычислять только поточечно. Приведенный вектор деформаций

ієУ =

строится с помощью обычных соотношений теории упругости

= дух = дуу = дУх дУу

Єхх = ~дх , Єуу = ~ду , Уху ="д7 ~дх ,

= ду = дуу ду2 = ду2 дух

ЄЄ ^ , 7 уг ^ ^ , 7гх ^ ^ ,

дг дг ду дх дг

где производные по х, у, г вычисляются на основе матрицы Якоби [ 3 ] следующим образом

дух ~дУх дУу дУ; '

дх дх дх д# д#

дУх дУу =Г з11 дух д

1ду ду ду і і дц дц дц

дух дух дУу дУ;

дг дг ді г д- д- д-і

Предполагается справедливость закона Гука, т.е. приведенный вектор напряжений

7Т = {7 , 7 , 7 , г , т , г }

( хх ’ уу’ 22’ ху ’ уг^ 2Х у

связан с (2) при помощи линейных соотношений

7} = [ ^]{^},

где [ Б] — матрица упругих постоянных материала.

Составляется выражение для потенциальной энергии деформации

V = {№г}Т И{ФП .

Аналитически вычислить все встречающиеся в (3) интегралы практически невозможно, поэтому обычно используется численное интегрирование. Так как в области изменения параметров £, 7, С элемент представляет собой куб, а подынтегральные функции — полиномы, то наиболее точные значения дают квадратурные формулы Гаусса-Лежандра по каждой из координат £,7, С . При этом

следует учесть, что элемент объема ёП в этих координатах тоже выражается через матрицу Якоби [Д] в виде

ёП = det [3] с1^с17с1С .

В результате фактически вычисляется сумма

у=2 "} 71 ^ * > ) }

2 • 1 *

• { \ 71), С( *))} [ Т ), 7( 1), С( *) )] ,

где 1,1, * — индексы суммирования по

квадратурным точкам соответственно в направлении осей £,7, С ; £('),7< 1 ),С<*} —

координаты квадратурных точек; а1, а), ак —

весовые множители. Порядок квадратурной формулы определяется исходя из необходимости точного вычисления главных членов в выражении энергии. Для квадратичной аппроксимации в плане конечного элемента этот порядок равен — 3 х 3 х 2.

Особенностью представленной группы элементов является тот факт, что использование классических трехмерных соотношений упругости (закон Гука) приводит к неверным решениям при отличном от нуля коэффициенте Пуассона.

Так как удельная потенциальная энергия является инвариантной величиной по отношению к системе координат, то, если ввести декартову систему координат, ориентированную относительно нормали, можно учесть гипотезу о малости напряжений обжатия, характерную для теории оболочек с учетом деформаций поперечного сдвига. В соответствии с этим в каждой квадратурной точке вводится местная декартова система координат

х , у , г с ортами рі, р2, р3, ориентированная таким образом, чтобы орт р3 совпадал с нормалью. Для определенности будем ее выбирать в виде

р3 = дг/д#х дг/д/дг/д#хдг/д/,

Р = дт/д^\дт/д#, р2 = рз хр1.

„ дг дг дг

Разложим векторы —, —, —

д# д/ дС

по ортам

р1, р2, р3. Получим дг

дг

= 4 = Ві рі + В2 ^

д# д/

дг

"ТТГ = Сі р1 + С2 р2 + С3 рз ,

дС

где

. дг _ л дг _ дг _

4і = тт-Рр в =т— Р,,С, =77'Рг.

д# д/ дС

Тогда можно получить явные выражения производных по местным координатам

' д д

дХ д#

д . = [*]• д

ду' д/

д д

.дг\ [дС\

где

[*]=

-ч

Ві

ЛА

Сі ВС

і -+- і 2

0 0

і

0

В

С2 і

В2С3

энергия

АС АВС

В этих осях потенциальная деформации может быть представлена в виде

V = {|[{е'}Т {а'}dQ.,

где

{е} = {ехх, еуу, е.., Уху, Уу., Ух },

{°'}Т = {, СТуу , стй,, Тху , V, г,., }.

В практической реализации используются модифицированные физические соотношения, в которых напряжения обжатия а.,., зависят только от деформаций е.,., в виде упрощенного закона Гука:

О,,

= Е3Є,’,’

напряжения о,, х,, Оу

соответствующими деформациями е

соотношениями, справедливыми для

напряженного состояния, т.е.

Е

связаны с

,,, є ,,, у ,,

с х ? у у ? / х у

плоско-

О, , =-

і-л

Е

(єхх + мєуу),

ієуу + мє,), 7ху = Суху,

уу і-м2

а в зависимостях напряжений поперечного сдвига т.,,, тл, от соответствующих деформаций у,,,, у,х,

вводится коэффициент Файлона к, учитывающий различие между реальным распределением этих напряжений по квадратичному закону и распределением, следующим из соотношений упругости и линейных аппроксимаций перемещений по толщине (постоянные значения). Здесь введены оси X, у', г ', которые ориентированы относительно нормали таким образом, что г' направлена вдоль нее, а X, у' определяют касательную плоскость к эквидистантной поверхности. Вводятся они

самостоятельно в каждой квадратурной точке по схеме, описанной ранее. Таким образом,

соотношения упругости принимаются в виде

а }=[ о' ]{е' },

где

" Е * /лЕ * 0 0

мЕ* Е* 0 0

= = 0 0 Ез 0

0 0 0 О

0 0 0 0

_ 0 0 0 0

Здесь Е* = Е

і-м2

0

0

0

0

О*к

0

Е3

0

0

0

0

0

О*к

— модуль

поперечного обжатия, О — модуль поперечного сдвига, к — коэффициент Файлона.

Переход к деформациям и напряжениям в осях х, у, г осуществляется с помощью

направляющих косинусов ортов рі, р2, р3,

определяющих локальные оси х’ , у’, г ’. Если ввести матрицу

где

[ рв ] =

[РС ] = 2

[ ро ] =

РххР: РухР Р,’хР>

Рх’хРу’ Ру'хР,'-,

Р,’хРх’

РххРуу

РухР,'у

[р] " рл = рс 1 1 ,

Р2 Гх'у РІ’,

РуХ р у г уу Ру, ,

РІХ р2 Г , у Р2 ’ г , ,

х’у Рх’ ,Рх’х

уу Ру’уРу’, Р у, Рух

1 у Р,’уР,1 Р,’, Р,’х

Рх’уРу’у

РууР,у

Р,’уРх’у

Рх,Ру’, Ру,Р,1 Р*Рхг

РхуРу, Рх,РуХ

РууР,1

Ру,Р,х

Р,’хРх’у Р . Р , Р . Р , Р , Р ,

г ухг ху г ууг X, г у,г хх

+ Р,'хРу 'у Р,'уРу ’, Р,',Ру 'X

Р, Р, Р, Р, Р, Р,

г ххг ,у г хуг 21 ГХ2Г ,х

то матрица упругих постоянных будет иметь вид

[в ] = [р ]т [о ’ ][р ]

Р,’уРх’, Р,’,Рх’

и

и в каждой квадратурной точке вычисляется самостоятельно.

2. Трехмерный изопараметрический элемент оболочки средней толщины с двойной аппроксимацией деформаций по точкам суперсходимости

Для удобства построения подобного КЭ вводится тройная нумерация узлов в соответствии с рис. 1 . Радиус-вектор аппроксимируется в виде биквадратичного полинома по координатам #,7 и линейного — в поперечном направлении (по толщине), т . е .

Г (#,7,С) = £ £ г ( т),7( п),С(к) )•

т,п=1 к=1

•Н„ (£)Ит (ч)1к (С), где нп (#) нт (л) и ^ (#) задаются

соотношениями

Н (#) = 1 #(#-1), н 2 (#) = 1 -#2,

Н (#) = 2#(# +1),

А (0 = 1 (1 -8) ^ (0 = 1 (1+0-

Вектор перемещений определяется в виде

V(#,7,С) = X XV(т)Уп),С(к))•

т, п=1 к=1

•Н„ (£)Ит (7)4 (С),

что обеспечивает линейное изменение перемещений по толщине оболочки. Вводятся ковариантные компоненты де формаций:

дг ду дг ду = дг дУ

>СС =~дС"д£’

дГ ду дГ дУ

В## д#' д#’ 87,71 дц дц'

дг ду дг ду

#ц д# дц дц д#’ 88 дц д8 д8 дц

є#8

дг дУ дг дУ

+ -

д# д8 д8 д#

Если ввести в рассмотрение тензоры деформаций

1

1

є XX 2 Ух'у' 2 ух^'

(8) = 1 2 у єУУ 1 2 Уу

1 Iу 1 2 Угу є

Ґ є## 1 2 У#ц 1 ^ 2 у##

8) = 1 2 Уг># є 1 2 Уц8

1 і і У 88 1 2 У8ц є88 )

то они связаны между собой следующим образом

и=[* ]даг.

При вычислении деформаций (8) через перемещения используется метод двойной

аппроксимации. В соответствии с этой методикой деформация є## вычисляется в точках

# = ±1/73; ц = -1,0,1; 0 = ±1, и по этим значениям строится квадратичная аппроксимация вдоль координаты ц и линейные вдоль # и 8 -

Деформация єцц вычисляется в точках

# = -1,0,1; ц = ±1/л/3; 8 = ±1, и для нее определяются линейные полиномы по ц и 8 , и квадратичные — по #. Деформация обжатия є0 вычисляется в точках #,ц = -1,0,1; 8 = 0, что дает квадратичные аппроксимации по # и ц, и постоянное значение по толщине. Деформация уц

вычисляется в точках #, ц = ±1/л/3; 8 = ±1 и определяется в виде линейных аппроксимаций по всем координатам. Деформация поперечного сдвига

у вычисляется в точках

#8

# = ±1/л/з; ц = -1,0,1; 8 = 0 , по которым строится квадратичная аппроксимация по ц, линейная — по

# и постоянная — по 8. Сдвиговая деформация у вычисляется при # = -1,0,1; ц=±1/-\/3; 8 = 0, и это дает квадратичный полином по # , линейный по ц и постоянный вдоль 8 .

Результатом таких построений будут

соотношения

є

= X X Нт](#Н, ц (8)-

і, і ,т=1 к=1

-г (, ц( т), 8( к)) - У (, ц{ і), 8( к)),

8ц = X X Ні (#Нт. ц (8)-

і, і ,т=1 к=1

-г (#( т), ц{ і), 8(к))- У (#(0, ц( і), 8(к)),

Є88 = Х X Н. (#Ні Ш(°)Ь'к(°)-

і, і =1І ,к =1

-г (ці \8(1))- V (V і' ),8( к)),

У#ц = 'X X ( ( (ц+ РЛ#Р.т ц)-

і, і,т=1 к=1

-ьк (8) г (,ц(”\8(к))- V (), ц( і), 8(к)),

у#8 = X X Ні(ц-

і, і,п=1 к=1

\Рп] (#К(0)4 (0) + Р (# (0)4(0)}-

-г (#(і),ц(п\8(І))-V(#('),ц(і'),8(к)),

Уц8 = 'X 'XНі (#-

і, і ,п=1 к=1

-{Рч,(ч)Ц(0) 4 (0) + Рп Ш, (0Ж (0)}-

-г (#(и),ц( і'), 8(І))- V (#(0,ц( і'), 8(к)), где Нті (#) определяется соотношениями

и

н„(0=2 (1 -03 )н;(0’ )н; (О )+ +2 (1+03) н;(т) н; (02>),

а Рп] (О). — соотношениями

р„ (О = 2 (1 ) Нп (О) н; (О)+

+2 (10) н (о2>) н; (О2’) .

Далее эти аппроксимации подставляются в выражение (9), позволяющее вычислить деформации (5), которые, в свою очередь, входят в закон Гука (6) с матрицей упругих характеристик (7).

Потенциальная энергия деформаций вычисляется путем численного интегрирования по всем трем координатам О, 7, С, в виде, аналогичном

(4), но в местных координатах х', у', 2', т.е.

V = 2 7 ',С‘к ’) [»']•

2 , ] к

• {г' (, 7( 7), С( к))} [ J ((0 , 7( 7), С( к))] ЩЮ]Щ.

В качестве квадратур используются формулы Гаусса-Лежандра порядка 3 х 3 х 2 соответственно по координатам 0, 7, С .

Подобные подходы к построению конечноэлементной методики расчета оболочек средней толщины также развиваются в работах [6-17].

3. Тестовые задачи

В настоящем разделе рассматриваются примеры численной реализации решения линейных задач. Рассмотрены две тестовые задачи, приводится сравнение с известными результатами.

3.1. Сферический купол. Рассматривается свободный сферический купол под действием системы самоуравновешенных сил (рис. 2). Эта задача

рассматривается в качестве тестовой в работах [16], где она решается с помощью 9-узлового КЭ с сокращенным интегрированием и специальной процедурой повышения ранга матрицы жесткости элемента и [1]. Параметры принимались следующими: Я = 304.8 см, к = 1.2192 см,

E = 4.В -106 кГ/см2, P = 0.907 кГ, ц = 0.3.

Рис. 2

Из соображений симметрии рассматривается 1/4 часть, разбитая на три КЭ. Точное решение для прогиба под силами дано в [16] и равно

м?т = 0.0563 см . В таблице 1 приведены значения погрешности

w - w

5 =-----т- -100%

wт

для КЭ, построенного по приведенной ранее методике в сравнении с результатами, полученными в [1, 16].

Таблица l

Prog 1 Prog 2 [1] [16]

S В. 9% 0.31% -0.4% 7.9%

«Prog 1» соответствует КЭ, техника построения которого приведена в разделе 1 настоящей работы, «Prog 2» - этому же элементу с двойной аппроксимацией деформаций (раздел 2 статьи).

3.2. Сжатый цилиндр. Рассматривается деформирование замкнутой круговой

цилиндрической оболочки со свободными торцами под действием самоуравновешенной системы двух сосредоточенных сил (рис. 3). Особенностью этой задачи является близость напряженного состояния чистому изгибу, что, как указывалось выше, представляет определенные трудности при использовании МКЭ. Наиболее обширные численные данные накоплены для оболочки со следующими параметрами:

L = 26.2 см, R = 12.5 см, h = 0.238 см,

E = 7.4 • 105 кГ/см2, P = 45.3 кГ, ц = 0.3125.

Аналитическое решение, полученное в предположении нерастяжимости оболочки в [7], дает значение прогиба под силой:

Wmax = -0.275 см ,

тогда как конечно-элементное решение дает несколько большее значение [17]

Wmax = -0.289 см .

В таблице 2 представлены значения погрешности (10) при вычислении прогиба под силой на различных сетках для 1/8 части оболочки с наложением условий симметрии на трех гранях для различных типов КЭ: «8-RI» — 8-узловой элемент с сокращенным интегрированием по работе [17]; «9-DA» — 9-узловой элемент с двойной

аппроксимацией деформаций, подобной описанной ранее [17]; а также КЭ из работ [1, 15]. Отметим, что 9-узловой элемент с сокращенным интегрированием в этой задаче приводит к сингулярной матрице жесткости и не дает никакого результата [16].

Таблица 2

сетка 2Х2 5Х5 10Х10

Prog 1 -52.З% -22.7% -8.6%

Prog 2 -10.1 % -2.В7 % -0.6%

[15] - -27 % -0.8%

[i] -1.7 % - -

8-RI [17] -15 % - -

9-DA [17] -7.З % - -

Заключение

Вычислительный эксперимент показал, что построенный KЭ без двойной аппроксимации деформаций (раздел 1 статьи) дает приемлемые результаты только при расчете пластин. При повышении кривизны рассчитываемой области точность расчетов начинает резко снижаться. KЭ с двойной аппроксимацией деформаций (раздел 2 статьи) достаточно эффективен для пластин и оболочек средней толщины, не «заклинивает» с уменьшением толщины, а также не приводит к сингулярной матрице жесткости в задачах со слабо связанной границей.

Работа выполнена при частичной финансовой поддержке РФФИ в рамках научных проектов № 12-01-009ЗЗ, № 12-01-97026

Литература

1. Голованов А.И. Современные конечно-элементные модели и методы исследовании тонкостенных конструкции / A.K Голованов, AB. Песошин, O.H. Тюленева. - ^зань: Е^зан. гос. ун-т, 2005. - 442 с.

2. Yang H.T.Y. A survey of recent shell finite elements / H.T.Y. Yang, S. Saigal, A. Masud, R.K. Kapania // Int. J. for numerical methods in engineering. - 2000. - V. 47. - P. 101-127.

3. Сахаров А.С. Метод конечных элементов в механике

твердых тел / A.G. Сaхaров, В.К ^слоок™, В.В. Kиричевский, И. Aльтенбaх, У. Габберт, Ю. Данкерт, X. ^пилер, З. - Киев: Вища шк., 1982. - 480 с.

4. Бережной Д.В. Искривленный конечный элемент пластин и оболочек средней толщины с учетом обжатия / Д.В. Бережной // Труды XVII междунар. конф. но теории оболочек и пластин. - ^зань: ^зан. гос. ун-т, 1996. - Т. 2.

- С. 94-99.

З. Баженов В.А. Mоментнaя схема метода конечных элементов в задачах нелинейной механики сплошной среды /

6. Sze K.Y. Three-dimensional continuum finite element models for plate/shell analysis / K.Y. Sze // Prog. Struct. Engng Mater. - 2002. - V. 4. - P. 400-407.

7. Тимошенко С.П. Пластинки и оболочки / С.П. Тимошенко, С. Войновский-Кригер - М.:Наука, 1966. -636 с.

8. Голованов А.И. «Трехмерный конечный элемент для расчета тонкостенных конструкций» / А.И. Голованов, М.К. Сагдатуллин // Ученые записки Казанского государственного университета. Серия физ. - мат. наук.

- Казань 2009. - Т. 151, кн. 3, с. 121-129.

9. Перелыгин О.А.Исследование прочности цилиндрических оболочек с вмятинами в области радиальных соединений / О.А. Перелыгин, Ш.Ш. Галявиев, Р.Х. Зайнуллин, Д.В. Бережной // Вестник КГТУ, Казань, изд-во КГТУ, 2001, №1-2, С.75-77.

10. Перелыгин О.А.Исследование прочности

цилиндрических оболочек при наличии увода или смещения кромок сварных швов / О.А. Перелыгин, Н.М. Туйкин, Д.В. Бережной, М. Н. Серазутдинов // Вестник КГТУ, Казань, изд-во КГТУ, 2001, №1-2, С.77-79.

11 Бережной Д.В. Многослойный ортотропный конечный элемент оболочек средней толщины. / Д. В. Бережной, М.К. Сагдатуллин, А.И. Голованов // Вестник Саратовского государственного технического университета. - 2011. -№3 (57). Выпуск 1. - С. 9-19.

12. Голованов А.И. Метод конечных элементов в механике деформируемых твердых тел / А.И. Голованов, Д.В. Бережной - Казань: Издательство «ДАС», 2001. - 301с.

13. Бережной Д.В. Расчет комбинированных конструкций

методом конечных элементов / Д.В. Бережной, М.К. Сагдатуллин, А.А. Саченков // Научно-технический вестник Поволжья. - Казань: Научно-технический

вестник Поволжья, 2012. - № 4.- С.13-16.

14 Бережной Д.В. Универсальный конечный элемент для расчета комбинированных конструкций / Д.В. Бережной, М.К. Сагдатуллин, А.А. Саченков // Вестник Казанского государственного технологического университета. - 2012. - №17. - С.150-157.

15. Сагдатуллин М.К. Статический расчет простых и комбинированных оболочечных конструкций МКЭ: дис. канд. физ.-мат. наук. Казан. фед. ун-т, Казань, 2012.

16. Belytschko T. Implementation and application of a 9-node Lagrange shell element with spurious mode control / T. Belytschko, W.K. Liu, J.S.J. Ong, D. Lam // Computers and Structures. — 1985. — V. 20, №1-3. — P. 121-128.

17. Park K.C. A curved C(0) shell element based on assumed natural-coordinate strain / К.С. Park, G.M. Stanley // J. Appl. Mech. — 1986. — V. 53, №2. - P. 278-290.

© Д. В. Бережной - канд. физ.-мат. наук, доц. каф. теоретической механики К(П)ФУ, berezhnoi.dmitri@mail.ru; М. К. Сагдатуллин - канд. физ.-мат. наук, доц. каф. теоретической механики и сопротивления материалов КНИТУ, ssmarat@mail.ru.