CC BY
103
Механика
Вестник Нижегородского университета им. Н.И. Лобачевского, 2013, № 1 (3), с. 19-24
УДК 539.3
РЕШЕНИЕ ОСЕСИММЕТРИЧНОЙ ЗАДАЧИ УСТОЙЧИВОСТИ УПРУГОПЛАСТИЧЕСКИХ ОБОЛОЧЕК ВРАЩЕНИЯ ПРИ КОМБИНИРОВАННОМ НАГРУЖЕНИИ МКЭ
© 2013 г. В.Г. Баженов1, Е.В. Павленкова1, А.А. Артемьева1, В.А. Иванов2, М.Н. Жесткое2
'НИИМ Нижегородского госуниверситета им. Н.И. Лобачевского 2Чувашский государственный университет им. И.Н. Ульянова, Чебоксары
Поступила в редакцию 10.12.2012
Излагается методика численного решения нелинейных нестационарных задач осесимметричного упругопластического деформирования оболочек вращения с учетом кручения при заданных кинематических и силовых нагружениях. Она основывается на геометрически нелинейной теории оболочек типа Тимошенко и теории пластичности с изотропным упрочнением. Решение задачи осуществляется вариационно-разностным методом в сочетании с явной схемой интегрирования уравнений движения по времени. Для иллюстрации эффективности данной методики проведены исследования предельных состояний упругопластического процесса деформирования цилиндрической металлической оболочки под действием внутреннего давления и последующего кручения.
Ключевые слова: упругопластичность, большие деформации, осесимметричные оболочки,
кручение, численное моделирование, вариационно-разностный метод.
Введение
Обзор работ, посвященных численным решениям осесимметричных геометрически нелинейных задач упругопластического деформирования и устойчивости оболочек при статических нагружениях, представлен в [1], а при динамических нагружениях - в [2]. Эти публикации в основном посвящены решению задач с малыми деформациями и немалыми перемещениями металлических оболочек. Численное решение задач о больших безмоментных деформациях гипер-упругих оболочек рассматривалось в монографии [3]. Постановка и решение обобщенных осесимметричных упругопластических задач с кручением в литературе отсутствуют. Такие задачи возникают при действии внутреннего давления и (или) растяжения с кручением. Они характерны для экспериментальных исследований поведения металлических трубчатых образцов при больших деформациях и неоднородном сложном напряженно-деформированном состоянии. Наиболее простым способом численного решения нестационарных упругопластических задач с большими деформациями является применение метода конечных элементов и явной схемы интегрирования по времени, которая позволяет осуществлять пошаговое перестроение геометрии оболочки. Современные компьютеры позволяют с приемлемыми вычислительными затратами осуществлять моделирование процессов дефор-
мирования в широком диапазоне скоростей нагружения от квазистатических до динамических. Отметим, что обзор математических моделей, постановок задач и численных методов решения нелинейных задач динамики оболочек представлен в работах [4 - 6].
Кинематические соотношения
Обобщенная осесимметричная задача (с кручением) формулируется в цилиндрической системе эйлеровых координат Оть , Oz - ось вращения. Для каждого оболочечного элемента введем местную сопутствующую (лагранжеву) систему координат 5, % [4]. Здесь 5 - длина дуги меридиана оболочечного элемента (0 < 5 < Ь); % - координата, отсчитываемая от срединной поверхности и нормальная к ней (—Ь /2 <%<Ь/2); Ь = Ь(5,?) - толщина оболочки; / - время.
Общая и местная системы координат связаны соотношениями:
5 = г^2 - , % = гуг + гу2,
где \у2 = —г,э, \уг = —z,s - направляющие косинусы нормали к срединной поверхности, символ буква 5 после запятой в нижнем индексе означает дифференцирование по пространственной переменной 5 .
Оболочечные элементы полагаем тонкими, изменением метрики по толщине пренебрегаем. Деформации поперечного сдвига и изгиба счи-
Исследование выполнено при поддержке Министерства образования и науки Российской Федерации, соглашение 14.B37.21.2019 «Проблемы динамического состояния сложных сред и конструкций».
таем малыми, а деформации срединной поверхности - большими. Деформациями сдвига е ^
пренебрегаем (р - окружная координата).
Распределение скоростей перемещений по толщине оболочки представим в виде:
й, О, £ 0 = 0 + (5,0 ,
(5, £ 0 = 0,0 , йр (5, 0 = Г<9(>, 0 .
Здесь мД5,?), и Л ■''./) - скорости перемещений срединной поверхности в направлении касательной и нормали; мД.?,?) - угловая скорость поворота поперечных сечений в плоскости меридионального сечения; 0(.\.1) - угловая скорость поворота относительно оси вращения; символ точка над переменной означает дифференцирование по времени.
Кинематические соотношения формулируются в скоростях и строятся в метрике текущего состояния, что позволяет учитывать большие деформации и формоизменения оболочки. Распределение компонент скоростей деформаций (симметричной части градиента скорости перемещений) по толщине примем в виде:
ёи=ё1 + &*,0 = *,Р),
ё,г=ё°г(1-(2 $1Ь)2),ё,р=ё%,
8 =и Щ -и Щ ,
^ т г т г '
ё°РР = йгГ\ё* =мЛ„
Л
2
1 л
ставляющих:
ё = ёе + ёр, ёр = 0, /,/ = 5,Д£.
г/ 1/ 1/ ’ в ’ ’ ^
Компоненты тензора напряжений Коши, нормальные к срединной поверхности оболочки, принимаются равными нулю (о^ = 0). Из
этого условия и обобщенного закона Гука опре-
• е
деляются компонента , компоненты скоростей напряжений а и напряжения ст.:
ё:е=^(ё>ё;д
/и-1
Е
1^7
£
(3)
аРР 1 ,,2
1-^
(ё® + /лё^),
& я = -
" ■> °'
(1)
где ё° - компоненты скорости деформации
*-» £ • И __
срединнои поверхности, деи — компоненты скорости деформации изгиба.
Компоненты тензора скоростей деформации и скорость вращения элемента выразим
через компоненты обобщенных скоростей перемещений йг, й2, 11^ , в :
(2)
ё ,= — (й +й ш +й ш ), б) „ = — гв, .
2 Г г^Г г'5 эр 2
Физические соотношения
Учет упругопластических свойств материала оболочки (металлы, сплавы) осуществляется в рамках теории течения с нелинейным изотропным упрочнением. Предполагается, что упругие деформации малы, а пластические деформации могут быть большими. Компоненты тензора скоростей деформаций ё можно представить в
виде суммы упругой ё* и пластической ёр со-
1^ * \ + /и ^
В^у = + <7*©^ , и ], к = 5, Р,
г г
= I, /, 7 = 5,/? ; сг^ = |,
о о
где Е - модуль Юнга, ^ - коэффициент Пуассона, символ ^ обозначает производную по Яуманну [7 - 9], которая учитывает поворот элемента оболочки за счет сдвиговой деформации как квазижесткого целого относительно нормали к ее срединной поверхности при кручении. При этом местная сопутствующая система координат 5, % определяется формоизменениями оболочки без учета сдвиговых деформаций кручения, то есть она является ла-гранжевой только для осесимметричной деформации. Заметим, что деформации сдвига при кручении тонких оболочек ограничены из соображений устойчивости осесимметричного процесса деформирования. Металлические трубки при чистом кручении, как показывают эксперименты, теряют устойчивость при углах закручивания порядка ж /12 с образованием неосесимметричных форм. Поэтому введение производной Яуманна для учета квазижесткого вращения при кручении вполне обоснованно при использовании явной схемы интегрирования с малыми шагами по времени. В итоге существенно упрощается запись основных уравнений обобщенной осесимметричной задачи с кручением, так как за базовую систему взята лагранжева система координат осесимметричной задачи.
Скорости пластических составляющих деформации определяются ассоциированным законом течения:
. ?
ёр = Ха'. , сг'сг' = — сг2(лЛ а'.. = сг -сгд ,
У У’ У У 3 4 4 ' (4)
0 = (055 +0РР )/3,
___г
к = ЗПъ\^<ь,
О
здесь параметр X, тождественно равный нулю при упругом деформировании, определяется при пластическом деформировании из условия прохождения мгновенной поверхности текучести через конец вектора догрузки; о, (к) - радиус поверхности текучести; к - параметр Од-квиста.
Уравнение движения
Рассмотрим элемент оболочки, ограниченный контурными сечениями 5 = 0, Ь . Пусть к контурам приложены внешние усилия [Р (?), Р (?), р (0]*=оь и изгибающий момент
[МИ (?)]=0 ь. Кроме контурных нагрузок, на
оболочечный элемент могут действовать нагрузки рг (5, ?), рг (5, ?), распределенные по поверхности.
Уравнения движения оболочки следуют из уравнения баланса виртуальных мощностей работы механики сплошных сред. Учитывая принятые гипотезы теории оболочек, выполняя интегрирование по толщине и пренебрегая силами Кориолиса, получим вариационное уравнение:
N*^1 + Nр8ё°рр +
+Мр5ё^р + 2Т5ё°!;/д + М р{(йг — г6г)8йг +
эР /А\ г
+йг5йг + г2656) +
+-1 „иФ5К ~ Рг5К ~ Рг§йг^ -
схеме второго порядка точности. Скорости и перемещения вычисляются по рекуррентным формулам:
• ■ Агк+У2
гк+1/2 = гк-И2 .ф ч* АГ--
■> І ■> І у ’
/к+1 = /* + /к+У2Мк+1, Аґк+У2 = - (Д^+1 + Дґ*), (7)
і і і 2
Г = и ,и ,и ,0, а = г,г,ф,0,
С= М,/, 7 , і = Щ .
Алгоритм решения задачи
При решении задачи расчетная область разбивается на конические конечные элементы. Узлы сопряжения конечных элементов с толщинами Аі+1/2(5) определяются координатами
(5)
~[г(РгЗйг + Рг5йг +МК5й(р + Рр86)\__^ = 0, где N, ^ , Т , Q , М , Мр - внутренние усилия и моменты; Мр, У - масса и момент инерции, определенные интегралами:
Ь/2 Ь/2
Г (5), (5), г = 0,Т — 1. Задаются начальные,
граничные условия и поверхностные нагрузки. Начальный временной шаг Д?1 выбирается согласно критерию устойчивости разностной схемы [10].
По известному на момент времени ? = ?к распределению параметров напряженно-деформированного состояния из уравнения движения определяются компоненты скоростей перемещений (м..й_.йгр. 6)к+У2, перемещений
(иг, иг, м , ^)к 1 и новая геометрия оболочки:
т*+1 = т° + ^+1, т = г, г, р, / = м, и2,6,
7 = 1Т
Дальнейшим этапом расчета на временном слое является вычисление компонент тензора скоростей деформаций по формулам (1), (2).
При для нахождения при-
ращений пластических деформаций (Дер)к+1У22 организуется итерационный процесс «посадки» девиатора напряжений на поверхность текучести (4). При этом выполняются следующие операции.
1. Вычисляются:
Ні = [ , Мі = [ , і = 5,Р , (6)
-А/2 -А/2
А/2 А/2
0 = | ^(1 -(2£/А)2Ж , Т = | ^,
- параметр ДЛ =
Мр=\ рй% , У, = | р%2й% .
—Ь/2 —Ь/2
Решение полученной системы уравнений (1)
- (6) строится по явной конечно-разностной
2 (О + q| 3) ’
О = о; + 2ОДер (в первом приближении
Дер =0, о =О);
- приращения пластических деформаций (Де' ) =ДХ-°'г;
- приращение параметра Одквиста
Дк = -ЛТз -^ДєрД£р ;
- по диаграмме деформирования о, (к) находится q = [О, (к* + Дк) — о, (кк ^ / Дк ;
- упругие составляющие приращений деформаций Дее = Дек+1 — (Дер)к;
- напряжения о по формулам (3);
- шаровая составляющая тензора напряжений о = (о55 +орр)/3;
- компоненты девиатора напряжений о; =о. —о.
V V
2. Проверяется условие текучести
1 4°°
расчетных и аналитических данных.
л/27зо,
(Де%%):
%% //+1/2
(Де55 + ДеРр)
к+1/2
рр//+1/2
Ь
—(Де55 + Де^р )к+1;?,
£+ = 1/2[1 + +А+)к+ ] .
< е = 10 3, где о, = о, (кк +1). При
невыполнении условия текучести осуществляется переход к пункту 1.
3. Находится приращение деформации Де^
и изменение толщины элементов оболочки:
\к+1/2 №
Далее по формулам численного интегрирования находятся усилия и моменты (6) и формируются коэффициенты уравнений движения для следующего временного слоя, где процесс вычислений повторяется.
Решение тестовой задачи о поперечном колебании упругой балки
Для тестирования методики решалась задача о синусоидальных поперечных колебаниях упругой балки, шарнирно опертой на торцах. Балка имеет следующие геометрические и физические параметры: длина Ь\ = 0.08 м, высота и ширина Ь2 = Ь3 = 0.01 м, модуль упругости Е = =210 ГПа, коэффициент Пуассона V = 0.3, плотность р = 7.8 г/см3. В начальный момент времени задано распределение скорости перемещений следующего вида:
й2 =А§т(тг г /Ц),
А = 10 м/с.
В силу симметрии геометрии, начальных и краевых условий рассматривается 1/2 часть балки. На рис. 1 представлено сопоставление численного решения задачи (сплошная линия) с аналитическим решением [11] (штриховая линия) в виде зависимости прогиба и = и2 / Ц в центре балки от времени т = ?с / Я
(с=<е р ). Получено хорошее соответствие
Рис. 1. Колебания упругой балки
Цилиндрическая оболочка под действием внутреннего давления и кручения
Проводилось численное моделирование деформирования оболочки с неподвижно защемленными торцами под действием внутреннего давления. Образцы выполнены из стали марки 12Х18Н10Т. Геометрические размеры оболочки: длина рабочей части 92 мм, внутренний диаметр 28 мм, наружный диаметр 30 мм. При численном моделировании использовалась диаграмма деформирования, полученная из эксперимента на растяжение аналогичной оболочки (рис. 2) [12].
0,2 0,4 0,6 0,1
Рис. 2. Диаграмма деформирования
На рис. 3 приведено сопоставление изменения формы оболочки в процессе нагружения, полученное в результате численного решения задачи (сплошная линия), с данными, полученными в работе [13] (штриховая линия). Введены следующие обозначения: Кй - начальный радиус оболочки, Ь - начальная толщина оболочки, Ц - начальная длина оболочки, г - расстояние от торца оболочки, q - внутреннее давление, и - радиальное перемещение, от - предел упругости, q* = qR / (Ь0от). Различие по полученным прогибам не превышает 1%.
О 0,2 0,4 0,6 0,8 г/Іо
Рис. 3. Изменение формы оболочки в процессе нагружения
I 1
/
1 5/ а)
1 1,3 1,6 1,9 2,2 9* 0 0,8 1,6 2,4 3,2 Р
Рис. 4. Максимальные прогибы оболочки
Осуществлено исследование раздутия цилиндрической оболочки при сложном нагружении. На оболочку действовало внутреннее давление до достижения параметром q * величин, отмеченных маркерами 1-6 (рис. 4, а), затем при неизменном давлении один торец оболочки закручивался с постоянной угловой скоростью, а второй оставался неподвижно защемленным. На рис. 4 представлены максимальные прогибы иг / ^ на плоскости симметрии оболочки при нагружении давлением в зависимости от параметра q * (а) и при нагружении давлением-кручением в зависимости от угла закручивания Р (б). Численное решение задачи показывает увеличение прогиба оболочки при кручении и неизменном давлении. Таким образом, в зависимости от предварительно достигнутого уровня давления при наложении кручения могут быть достигнуты предельные состояния оболочки.
Заключение
Разработанная методика и программа позволяют проводить численное моделирование осесимметричных упругопластических процессов деформирования и предельных состояний оболочек вращения в широком диапазоне скоро-
стей нагружения от квазистатических до динамических. Методика проста в численной реализации и позволяет рассматривать оболочки с любым очертанием меридиана. Приведенный пример решения задачи о выпучивании цилиндрической оболочки под действием внутреннего давления и последующего кручения демонстрирует ее работоспособность и эффективность при больших деформациях.
Работа выполнена в рамках ФЦП «Научные и научно-педагогические кадры инновационной России» на 2009-2013 годы (соглашения 14.132.21.1812, 14.В37.21.2019), при поддержке Совета по грантам Президента РФ для ведущих научных школ (грант НШ-2843.2012.8), при финансовой поддержке РФФИ (проекты 11-08-00565а, 12-08-31190-мол_а, 12-08-33106-
мол_а_вед, 12-08-12044-офи_м).
Список литературы
1. Коробейников С.Н. Численное решение уравнений с особенностями деформирования упругопластических оболочек вращения // Вычислительные технологии. 2001. Т. 6, № 5. С. 39-59.
2. Выпучивание упругопластических цилиндрических и конических оболочек при осевом ударном нагружении / В.Г. Баженов, М.С. Баранова, А.И. Ки-бец, В.К. Ломунов, Е.В. Павленкова // Ученые запис-
ки казанского государственного университета. Физико-математические науки. 2010. Т. 152, кн. 4. - С. 86-105.
3. Колпак Е.П. Устойчивость безмоментных оболочек при больших деформациях. СПб.: СПбГУ,
2000. 248 с.
4. Баженов В.Г., Ломунов В.К. Устойчивость и закритическое состояние оболочек вращения при осевом ударе // Прикладная механика. 1986. Т. 22, №
9. С. 28-33.
5. Абросимов Н.А., Баженов В.Г. Нелинейные задачи динамики композитных конструкций. Н.Новгород: ННГУ, 2002. 400 с.
6. Баженов В.Г., Чекмарев Д.Т. Численные методы решения задач нестационарной динамики тонкостенных конструкций // Механика твердого тела.
2001. № 5. С. 156-173.
7. Поздеев А.А., Трусов П.В., Няшин Ю.И. Большие упругопластические деформации: Теория, алгоритмы, приложения. М.: Наука, 1986. 232 с.
8. Аннин Б.Д., Коробейников С.Н. Допустимые формы упругих законов деформирования в определяющих соотношениях упруго-пластичности // Сиб. ж. индустр. математики. 1998. Т. 1, № 1. С. 21-34.
9. Баженов В.Г., Павлёнкова Е.В., Артемьева А.А. Численное решение обобщенных осесимметричных задач динамики упругопластических оболочек вращения при больших деформациях // Вычислительная механика сплошных сред. 2012. Т. 5, №4. С. 427-434.
10. Баженов В.Г., Чекмарев Д.Т. Решение задач
нестационарной динамики пластин и оболочек вариационно-разностным методом: Учебное пособие.
Н.Новгород: ННГУ, 2000. 107 с.
11. Тимошенко С.П. Колебания в инженерном деле. М.: Наука, 1967.
12. Баженов В.Г., Ломунов В.К. Экспериментально-теоретическое исследование процесса образования шейки при растяжении стального трубчатого образца до разрыва // Проблемы прочности и пластичности. Межвуз. сб. Н.Новгород: ННГУ, 2001. С. 35-41.
13. Баженов В.Г., Ломунов В.К. Влияние статического давления на устойчивость упругопластических цилиндрических оболочек при продольном ударном нагружении // Прикладные проблемы прочности и пластичности: Всесоюз. межвуз. сб. / Горьковский ун-т, Горький. 1979. Вып. 12. С. 39-42.
SOLUTION OF THE AXISYMMETRIC PROBLEM OF STABILITY OF ELASTOPLASTIC SHELLS OF REVOLUTION UNDER COMBINED LOADING BY FINITE ELEMENT METHOD
V.G. Bazhenov, E. V. Pavlyonkova, A.A. Artemyeva, V.A. Ivanov, M.N. Zhestkov
We present a technique for solving numerically the nonlinear nonstationary problems of axisymmetric elastoplastic deformation of shells of revolution taking into consideration torsion under the prescribed kinematics and loads. This method is based on the geometrically nonlinear timoshenko theory of shells and the plasticity theory with isotropic hardening. A variational difference method in combination with a scheme for explicit time integration of the equations of motion is used for solving the problem of interest. The efficiency of the proposed method is verified by studying the limit states of the elasto-plastic deformation of a cylindrical metal shell under internal pressure and subsequent torsion.
Keywords: elastoplasticity, large deformation, axisymmetric shells, torsion, numerical modeling, variational difference method.
