Математическое моделирование динамического деформирования упруговязкопластических оболочек конечной длины лучевым методом Текст научной статьи по специальности «Физика»

Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки. 2018. Т. 22, № 2. С. 325—343 ISSN: 2310-7081 (online), 1991-8615 (print) cjjhttp://doi.org/10.14498/vsgtu1610

Математическое моделирование, численные методы и комплексы программ

УДК 539.3

Математическое моделирование динамического деформирования упруговязкопластических оболочек конечной длины лучевым методом

Н. Д. Вервейко, М. В. Егоров

Воронежский государственный университет,

Россия, 394018, Воронеж, Университетская пл., 1.

Аннотация

В работе представлено математическое моделирование динамического напряженно-деформированного состояния оболочки вращения из упруговязкопластического материала. Решается модифицированная система уравнений в частных производных типа С. П. Тимошенко путем построения системы уравнений на подвижных поверхностях разрыва с начальными условиями в виде ударной нагрузки на торце, записанной в виде степенного ряда по времени, коэффициенты которого есть начальные условия для дифференциальных уравнений. Решение представляется в виде лучевого ряда Тейлора с точностью до четвертого порядка по координате оболочки. Для моделирования отраженных волн от границ вводятся условия двух типов на границе (жестко защемленной и свободной от напряжений), не зависящие от времени. Разработан комплекс программ, написанных на языке Fortran 90 на платформе Code::Blocks. Реализованы 2 программы для моделирования динамического деформирования оболочки в упругом и упруговязкопластическом состоянии. Использовано разностное представление производных, вычисление интегралов методом трапеций с заданным шагом разбиения отрезка. Результатом работы программ являются сеточные функции коэффициентов рядов Тейлора, которые используются для построения графиков перемещений как функций времени и продольной координаты оболочки. * 3

Научная статья

3 ©® Контент публикуется на условиях лицензии Creative Commons Attribution 4.0 International (https://creativecommons.Org/licenses/by/4.0/deed.ru)

Образец для цитирования

Вервейко Н. Д., Егоров М. В. Математическое моделирование динамического деформирования упруговязкопластических оболочек конечной длины лучевым методом // Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки, 2018. Т. 22, № 2. С. 325-343. doi: 10.14498/vsgtu1610.

Сведения об авторах

Николай Дмитриевич Вервейко © http://orcid.org/0000-0002-6980-9217 доктор технических наук, профессор; профессор; каф. механики и компьютерного моделирования; e-mail: verveyko2017@yandex.ru

Михаил Валерьевич Егоров А © http://orcid.org/0000-0002-8977-7320

аспирант; каф. механики и компьютерного моделирования; e-mail: egorovmv89@mail.ru

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

Получение: 23 февраля 2018 г. / Исправление: 21 мая 2018 г. / Принятие: 11 июня 2018 г. / Публикация онлайн: 1 июля 2018 г.

Введение. Задача динамического напряженно-деформированного состояния осесимметричной оболочки из упругого материала в постановке С. П. Тимошенко [1,2] сводятся к решению квазилинейных систем дифференциальных уравнений в частных производных гиперболического типа. Учет необратимых деформаций при динамическом деформировании возможен с использованием модели упруговязкопластического материала [3]. Модифицированные уравнения динамического деформирования оболочек вращения [4] для упруговязкопластического материала имеют вид

d 2и dt2 д2ф dt2 d2v dt2

ер-

v

,d2u

dv

cf

0 ^

1 dz2 R2(z) dz

d2^ 12c2k2 dv

C2V

4——1—r v----

dv 2p ddzz

R(z) dz 12c2k2

p dz

dz2

,d2v

h2 dz

У2 // + <2 k2 ^ -

h2

c2v du

. 2cfv d2v

^ - Щг) VW°;

c2v (dv\2

dz

+

(OV\ 2

\dz)

c2v сЪ

1 -V + ^т- V

R(z) dz 2R(z)\dz,

d2v 3^k2v 2 2p depz

R2(z)

dep-

а^з _ fp _ _ dt ij

d p

--P ■ ■

dt 13

6ep-

St

R(z) dz2 2R3(z)

1 f d2(a') - КУ2 p\ l2(a')

h(P) =

v —

p dz

) aij;

— c

dep-

dn

<з <i ■

(1)

Здесь u — продольные перемещения, м; v — поперечные перемещения, м; ф — угол поворота срединного сечения, рад; z — координата оболочки, отсчитываемая от левого торца; t — время, с; с2 — скорость распространения продольного возмущения в неограниченной среде, м/с; с2 —скорость распространения поперечного возмущения в неограниченной среде, м/с; R(z) — радиус срединной поверхности оболочки в зависимости от координаты z, м; р — модуль сдвига, Па; р — плотность, кг/м3; е?- —тензор пластических деформаций; h — толщина оболочки, м; к — параметр, характеризующий распределение касательных напряжений в сечении; v — коэффициент Пуассона; —скорость остаточных пластических деформаций, 1/с; р — коэффициент динамической вязкости материала, Па-с; ( ■) — полная производная по времени; I2(ad) — второй инвариант девиатора тензора упругих напряжений, Па; К — предел пластичности, Па; ad — девиатор тензора упругих напряжений; (■) — частная производная по времени; ^ (■) —производная по времени на фронте возмущения; с — скорость движения фронта возмущения.

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

го деформирования £?• при выполнении условия неупругого деформирования [5].

Наряду с аналитическими методами, пригодными для частных случаев, широкое применение имеют численные методы [6-8], существенной особенностью которых является выбор корректного типа конечно-разностной схемы, пространственного и временного шага. Зачастую реализация таких задач требует привлечения больших вычислительных мощностей. Лучевой метод [5,8-12] предоставляет возможность построить степенные ряды повышенного порядка точности решения системы нестационарных дифференциальных уравнений в частных производных и соотношения на подвижных границах, которые движутся со скоростью распространения возмущений. Задача определения коэффициентов ряда сводится к решению систем обыкновенных дифференциальных уравнений на подвижных границах. Численное решение таких уравнений дает выражение для перемещений, деформаций, их градиентов, градиентов второго порядка на подвижных границах. Решение в окрестности волнового фронта строится в виде продолжения степенного ряда Тейлора.

1. Постановка задачи численного моделирования динамического деформирования оболочки вращения лучевым методом. Рассмотрим динамическое деформированное состояние оболочки вращения (рис. 1), порожденное ударным продольным нагружением на левом торце х = 0, с граничными условиями на правом торце х = L, соответствующими свободному от напряжений торцу или жестко защемленному.

Решение задачи для перемещений и и v и угла поворота срединного сечения ф за фронтами волн £i, £2 может быть представлено в виде степенного ряда Тейлора по расстоянию п, отсчитываемого от соответствующей волны, с коэффициентами, зависящими от времени t:

Коэффициенты лучевого ряда (2) находятся из последовательно решаемой системы обыкновенных дифференциальных уравнений переноса скачков функций на фронтах волн £1, £2. Система уравнений (1) для и, V, ф допускает разрывные решения на подвижных поверхностях, распространяющихся вдоль оболочки со скоростями [11]

Е

Е

Р(1 - v2)

,2)’

2р(1 + и)

х

Рис. 1. Модельное представление оболочки вращения [Figure 1. Model representation of a shell of revolution]

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

Отметим, что из условия непрерывности перемещений разрывы перемещений [и], [г>] отсутствуют: [u] |s= 0, [v] |s = 0, что ведет к непрерывности угла поворота срединной поверхности ф на £? и £2: [ф] = 0, так что представление (2) для перемещений и угла поворота начинается с к = 1, соответствующего градиенту перемещений.

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

ж - 2^ + с2['“”“] =

=с?[,1“‘] + ты + - ^ [eL’”];

сЪ

2^

^ - 2с^Ц ^„,] =

St

^ - 2сФф + Х[,,„„] =

2 12с2к2 12с^к2 2с\и

= с-Жп] —Ы —[1\ - МЫ;

= c2k2[vnn] + с^к2[фп\ - ДД[ип\ + Д [vn\2-

„2,, г2,.

„ ,______, , C?U [v 12

п] ад[“га] + 2R(n)[n

c?v , r 1Г , 3c2k2v r l2 2u

1 + MK«] - 2 r-12

R(n)2 R(n)

[v]2 _ ^[eP ]

2R3(n)[U] pc [Crz’n]'

(3)

Скачки градиентов пластических деформаций [e^z,n ], [Zrz,n] определяются из выражения для скорости пластической деформации:

Г dep-1

[4] = Г * 1

8[% ]

\9Д 1 \9Д ]

- С = - с

- dt J 5t L dn \ L dn -Вследствие непрерывности пластических деформаций на фронтах волн

(4)

К, ] = 0.

На фронте разрыва скоростей скачки градиентов пластических деформаций [e%j,n] определяются из скачка скорости пластической деформации:

г del 1

1

1

= -7 4 ] = 4

h((Tr) - кV2 ,

-а,

ч

. dn J cL~гр cpl 12(^г)

где скачки напряжений определяются из закона Гука через скачки градиентов скоростей:

[ер

PzzJ cq(1

к Д2

л/2/3 (Х([и,п } + 1/r[v}) + 2p[u,n})2 + 2p?[v,n]2 -

х (Л([и,п} + 1/r[v}) + 2^[и,п\);

[еР } =______111 -

\Утгчп\

Ку/2

Л/2/3(Х([и,п] + 1/r[v}) + 2 p[u,n})2 + 2^2[v,ra]2

)Ф,п}

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

Е п

Рис. 2. Представление разрыва функции [Figure 2. Representation of the jump of a function]

Система ОДУ (3) для [на £1 ([u,n], [ф,п], [v,n] = 0) записывается в виде

на £2 (К„]

S[u,n]

5N

+

р(Х + 2р.) рс\р

[u,n]

5[Ф,п]

5п

у/зкр,

рС3р

0, [ф,п] = 0) —

ф,п] + £У +

4Е(п)с2&2

\v,n}2 +

^2

рс2к3р

[V,n]

0;

рс2к3р’

(5)

(6)

(7)

где У = ct — расстояние, пройденное фронтом £1 со скоростью с.

Дифференциальные уравнения для членов ряда (2) при к = 2, 3,... получаются путем дифференцирования системы (1) по п необходимое число раз и представления ее в виде разности на волновом фронте.

Решения уравнений для скачков градиентов скоростей на £1, £2 для (5)— (7) имеют вид

[v,n]

где

-КУ+адп ^зК м(А+ЗДп ^зж ^

[“'"} =е (x+y +[“-}»-x+d;

.А + 2^

[фуп] = [Ф 5«.}0;

^^^tg((C2i + n)yj621C21 - 0.25а21 yj621С21 - 0.25а21 - 0.5&2i),

&21

с2^

4R(n)c2k2 ’

^21

^2

рс2 к3р' С21 pc2k3p,

С21

1

Д21С21 - 0-25&21

arctg ^

&2i[^,n]o + 0-5а21 Vb2lC2l - 0.25&21

.

2. Граничные условия на левом торце оболочки. В качестве начального условия (при t = 0) примем состояние покоя по всей длине оболочки вращения.

Дифференциальные уравнения переноса для [u,n], [u,nn], [и,ппп], [ф,п], [ф,пп], [ф,ппп], [v,n], [v,nn], [v,nnn] решаются при начальных условиях [и,п](0) = [и,п]0 и т. д., которые определяются из разложения в ряд Тейлора по времени граничного задания нагрузки на левом торце оболочки:

&ZZ

z=0

<Х

£

к=0

d к azz (z,t)

dtk

- Ф

,=o\ Ы,

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

&rz (0) = ffrVzz (0),

ffr — коэффициент трения.

Из реологического уравнения для упруговязкопластического материала

&ij — Л&кк bij + 2f^^ij 2f^^ij (8)

можно получить выражения для скачков градиентов перемещений и градиентов более высокого порядка при z = 0, t = 0:

[^zz ]0 — ^[^,га]0 + 2^[^,га]0; [^,га]0

\&zz ]0

Л + 2^.

Дифференцирование соотношения (8) по времени t дает выражение для градиентов второго порядка при z = 0, t = 0:

' dazz -. dt -I0

Л[-^ ,^]0 + 2M^ ,n ]0 2^[ezz ]0;

ii0 =

[^,rara]0

du dt o, 2^ сщ

[ik,ra]0 -

[-^,ra]0

C['U,rara]o,

1

Лс1 + 2/лс1

[^j nn\0

dazz -dt J0

C [-U,rarara]0;

2ф3ц,К ( Л c 1 + 2/j.c 1)^'

Двойное дифференцирование выражения (8) по времени t позволяет определить

\d2azz 1 ... 1,0 г..-, 0 \деф 1

— ['^',га]0 + 2№Щ ,п]0 2^

- dt2 J 0

1 dt V

1 2^zz 1 .2/1

+-----[и, пп] 0.

[и,ппп]0 = с?( Л + 2ц)У dt2 J0 1 с1Г}

По аналогии с вышеизложенным запишем выражения для градиентов поперечного перемещения [v,n]0, [v,nn]0 [v,nnn]0 при z = 0, t = 0:

[^,ra]0 = ~2^[^rZ ]0;

[F,rara]0 —

1

~dar

[F,rarara]0 —

2ц,с2к1 dt Jo 1 Г d2a,

+

^ К

[V,n\<0 —

C2tq

+

С2^’

2ц,с2к21 dt2 Jo c2tq

[F,rara]0.

3. Граничные условия для перемещений и их градиентов за отраженными Si, S2 в случае жестко закрепленной границы. Представим в общем случае жестко защемленной или свободной границы схематичное положение отраженной волны следующим образом (см. рис. 3).

/+ J rej fref f. f+ J vnc J me

Ci At d At |

V' V- Г V-

Si Si L Si

Рис. 3. Схематичное изображение границы, падающей волны, отраженной волны и их гипотетического положения [Figure 3. Schematic representation of the boundary, the incident wave, the reflected wave, and their hypothetical position]

Условия закрепления оболочки состоят в отсутствии перемещений на конце z — L:

uz (z,t)\z=L — 0, vr (z,t)\z=L — 0. (9)

Поскольку условия (9) сохраняются в любой момент времени t, их можно дифференцировать частным образом по t несколько раз (достаточное число, обусловленное требованиями внешних условий):

дк дк

Ofiu* СМ) —0, ^vr (М) —0, к — 1, 2, 3. (10)

Условия (9), (10) преобразуем применительно к решению за падающей Sjrac и отраженной Eref волнами от закрепленного конца z — L.

Рассмотрим представление решения для u, v и ф в виде лучевого ряда Тейлора:

ГО

u(z,t) —

к=0

п

¥

dku(z, t) - dk u(z,t)'

dzk z=ct . dnk .

(11)

где п — расстояние до геометрической точки z, отсчитываемое от положения фронта волны Е (рис. 3), которая координатой z отсчитывается от левого конца оболочки z — 0; ( ■ )| —производная по геометрической координате,

вычисляемая за фронтом падающей волны слева, которая является значением этой величины перед фронтом отражающей волны (рис. 3).

Условие равенства нулю перемещения правого торца z — L оболочки запишем с учетом скачков функций:

и (z,t)Sref — u(z,t)^inc - [u(z,t)]Sref

Производная (частная) от перемещений по времени в подвижной точке z — ct за фронтом волны для первого члена ряда (11) примет вид

du(z, t) dt

Su ди S [и]

St °дп St

+ с

ди' дп. ’

В случае падающей прямой волны £пс впереди нее предполагается состо-

+ I du + I

яние покоя, т. е. иК^ и -q^-|s отсутствуют:

и+

du+

дп

s

0.

За падающей волной £пс производная по времени t от перемещения за £ref запишется следующим образом:

ди

~dt

5и

~st - с

ди' дп.

Для отраженной волны £ref решение и производная по времени t от перемещения за £ref примет вид

du(z, t) dt

6и+

~W

ди+ С дп

6 [и]

~W

+ с

ди

дп

•}

где и+ = 0 - [u], 1ST = - [ш].

Условие равенства нулю скорости правого конца оболочки z приводит к соотношению

т ди Ь, dt

0

д[и]

dt inc

+ С

-dw ад + с ди -

ldn\ гпс ut ref ldn\

0.

(12)

Так как при z = L скачок перемещений [и] = 0, из (12) следует граничное условие для скачка градиента перемещения:

гди,л ди I

.дп. ref .дп.

vac

(13)

Для скачков градиентов более высокого порядка от перемещений условия на границе z = L будут аналогичны (13):

\дк ил лдк и I

-дпк. ref -дпк-

гпс

(14)

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

[,di,n]ref

Узк

2д

4 Дп

_|_ g 3рс13^

[и,п] (L)inc

Узк -

2р .

(15)

Выражение (15) представляет собой скачок продольной деформации за отраженной волной от границы L.

Скачки градиентов более высокого порядка на отраженной волне могут быть получены из (14).

Решение за отраженной волной в случае жесткого закрепления правого торца оболочки будет представлено следующим образом:

fref (п,п) = П{ДДД - \fref j) + у - [^f ,пп^) +

. «3 ( 9Sfinc

6 V дп3

[fref

уППП

О.

4. Граничные условия и решение для перемещений за отраженными Si, S2 в случае свободной границы. Условие на свободной границе состоит в отсутствии напряжений на ней:

vzz(L, t) = orz(L, t) = 0. (16)

Запишем условие на свободной границе в скачках напряжений с учетом состояния покоя перед волновым фронтом £1 а также того, что напряжения за £ при z = L складываются из напряжений за £1 и £1:

Vine + Vref — 0. (17)

Для получения граничных условий для градиентов напряжений продифференцируем (17) по времени t на фронтах £inc и £ref:

9<Jinc . dVref dt + dt

St

(Vinc

+ Vref )

1A

c dn

(Vinc

+ Vref) — 0.

(18)

Условие (18), представленное в точке z = L в окрестности фронтов £inc и £ref, позволяет получить условие при z = L для градиентов напряжений с учетом ([oinc] + [oref]):

-дог

- дп -

Сами перемещения непрерывны:

.

до ref дп -

0.

Мтс — Мге/ — °.

Граничные условия для градиентов перемещений при условии (16) следуют из закона Гука [ozz] = у—д [u,n], [vrz] = [v,n\-

'dw ди ■

.дп. ref .дп.

vac

(19)

Для построения выражений скачков градиентов более высокого порядка для отраженной волны воспользуемся уравнением движения элемента оболочки на границе z = L:

,du'

a = РЦ,

где ш = и' — продольная скорость; р' = m/S, m — масса элемента оболочки, S — площадь элемента оболочки.

Запишем выражение для ускорения на фронте возмущений в окрестности свободной границы:

о

<> (

5и' ди''

— с-

St

дп.

0.

(20)

Преобразуем выражение (20), представив скорость и1 через локальную

А?» till

производную ^ и скорость переноса с^

на фронте возмущения:

St

Su ди St С дп

)

д /5и ди С дп\ St Сдп

0.

(21)

В терминах скачков (21) примет вид

ё2[и] 8 г ди

St2 С8t .дп

+ с

2

д2и

дп2

0.

(22)

Выражение (22) дает соотношения для градиентов второго порядка падающей и отраженной волн:

'д 2w _2 — 4 ди ■

Удп2\ ref С St 1дп\

гас

(23)

- с—:

с дп'

Для получения выражений для скачков градиентов третьего порядка на

Л _ s_ dt St

отраженной волне продифференцируем выражение (21) по времени 44 _ 44 —

S2 ди 5 д2и St2 дп + С St дп2 + 2 '

S д2и

— с

д 3 и

0.

St дп2 дп3

Выражение для скачка градиента третьего порядка отраженной волны принимает вид

г£>3„. 1 х2 гя„,-,

(24)

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

3

до третьего порядка по п запишется в виде

■д3и■ _ 4—42 ди ■

1дп3\ ref С2 St2 1дп\

гас

U

■ ди'

-дп\ ref

(N — п) +

(N — п)2 гд2и] (N — п)3 гд3и'

.дп2\ ref

+

6 1.дп3) ref'

В случае пластического деформирования материала за фронтом волны имеем соотношение для напряжений для случая упруговязкопластического динамического деформирования в скачках:

[Oij] _ \Щ8гз + 2^[eij] - 2^[epj].

Скачок пластической деформации на волновом фронте [ер-] _ 0. Следовательно, отраженная волна для упруговязкопластической падающей волны будет только упругой. На границе при этом будут выполняться соотношения (19), (23), (24).

5. Оценка пластической деформации элемента оболочки в окрестности фронта волны. На самом фронте волны пластического деформирования при условии [6] неупругого деформирования пластическая деформация отсутствует при упругом состоянии материала перед волной, т. е. фронт волны не добавляет пластическую деформацию непосредственно за фронтом.

Наличие скорости пластического деформирования при превышении предела пластичности (1) дает градиент пластической деформации за фронтом:

г dep-1

£р = _ с _v ч VdnY

Выражения для скоростей пластических деформаций на £1 имеют вид

de^z dt Si

= [^Lsi ] = -

de^z

(Л + 2Д), , л/3К

[U,n\ + •

V

V

= fei ] = 0.

dt s1

Выражения для скоростей пластических деформаций на £2 имеют вид

depzz

dt S2

— — о-

— bzzS2 — 0-

depz dt s2

— £rzS2 — - “ [V,™\ + —.

К

V

Знание скоростей пластического деформирования на фронте волны позволяет оценить пластические деформации за фронтом с точностью до величин 3-го порядка по расстоянию от фронта:

е-7 — И^ - *) - К,,-]

(N - z)2

1,

+ [ер- ]

L IJ )ППП1

(N - z)3

6

\ер ] — --[ер--

LS ni с [ttj-

1

хгР

\ер 1 —~Ц[ер- ] + \—^} V

LSпп\ Д[cij,n] + |_£wJ )-

Л £ГР я2рР

\еР 1 — ] + Г ч,п\ + Г____Щ)

\У,ппп\ Д} [ §]у YdN2\)

(25)

(26) (27)

Выражения (25)-(27) получены по аналогии с (4) для [е^-,п], [^,пп]-Выражения для остаточных перемещений после снятия нагрузки могут быть получены путем разложения остаточных перемещений в степенной ряд за продольной и сдвиговой волнами:

u(N,z)— epzz (N - z) - [epzz J(N о z) + [epzz ,nn](N z) ■

v(N, z) — eprz(N - z) - [eprz,n]^—^ + [eprz,nn]

6

(N - z)3 6 '

Второй инвариант тензора скоростей пластических деформаций вычисляется исходя из выражения [3]:

ы4) — ДД4,

где efj = — 35ijepkk — тензор-девиатор скоростей пластической деформации.

Тензор-девиатор скоростей пластической деформации представим в виде матрицы

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

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

6. Описание программы. Программа [13] написана на языке Fortran 90 в кроссплатформенной среде Code::Blocks. Программа реализует 2 случая динамического деформирования: упругий и упруговязкопластический. Результатом работы программы является массив данных коэффициентов рядов Тейлора для построения решения перемещений за фронтами волн в зависимости от времени. Далее из массива по средствам интерполяции создаются непрерывные функции (в Mathcad 15). Выбор случая расчета (упругий или упруговязкопластический) осуществляется путем сравнения максимальных эквивалентных напряжений по Мизесу на торце оболочки с пределом текучести. Блок-схема алгоритма программы представлена на рис. 4.

7. Численный пример. В качестве примера рассмотрим удар по левому торцу оболочки с постоянным усилием при наличии трения на торце: azz(0) = 250 МПа, arz(0) = ffrazz(0), ffr = 0.001. Геометрические параметры оболочки выберем равными L = 10 м, R = 1 м, h = 0.1 м, что задает безразмерные характеристики оболочки L = L/R, h = h/R. Безразмерный параметр к характеризует распределение касательных напряжений в сечении. Физические параметры выбраны следующим образом: с\ = 5422 м/с, с2 = 3208 м/с, ц = 1010 Па ■ с, К = 200 МПа, v = 0.3, что определяет fj = r\ci/(KV), с = ci/c2. Задание конкретных геометрических и физических параметров примера определяет единственное значение безразмерных комплексов L, h, к, fj, с, которые определяют подобие рассматриваемого явления при других геометрических и физических константах, сохраняя критерии подобия. Численные расчеты приведены для двух вариантов граничных условий на правом торце: жестко защемленном и свободном в виде u = u/L, v = v/L, N = N/L, V = z/L.

На рис. 5 представлены 3Б-графики изменения продольного и поперечного перемещения до момента возвращения волны Si к левому торцу оболочки для случая жесткого закрепления правого торца оболочки. В рассматриваемом примере на левой границе заданы положительные деформации arz (0, 0). На рис. 5 видны моменты отражения волны Si в точке (V = 1, IV = cit/L = 1).

На рис. 6 представлены формы мгновенного изменения поверхности оболочки в момент времени, когда волна Si дошла до правого торца оболочки,

Рис. 4. Блок-схема алгоритма программы [Figure 4. Flowchart of the program algorithm]

и в момент времени, когда волна £i вернулась к левому торце оболочки, для случая жесткого закрепления правого торца оболочки.

На рис. 7 представлены изменения остаточной формы поверхности оболочки вращения при снятии нагрузки в момент времени, когда волна £i дошла до правого торца оболочки, и в момент времени, когда волна £i вернулась к левому торцу оболочки, для случая жесткого закрепления правого торца оболочки.

На рис. 8 представлены ЭБ-графики изменения продольного и поперечного перемещений до момента возвращения волны £i к левому торцу оболочки для случая свободного правого торца оболочки. Здесь от точки (z = 1, TV = 1) начинается движение правой границы.

На рис. 9 представлены формы мгновенного изменения поверхности оболочки в момент времени, когда волна £i дошла до правого торца оболочки, и в момент времени, когда волна £i вернулась к левому торцу оболочки, для

й

0.018

0.016

0.014

0.012

0.010

0.008

0.006

0.004

0.002

0

v

0.018

0.016

0.014

0.012

0.010

0.008

0.06

0.004

0.002

О

-0.002

Рис. 5. Распределение продольных (сверху) и поперечных (снизу) перемещений в случае жесткой правой границы, рассчитанное до момента времени возвращения отраженного волнового фронта Ei к левой границе (онлайн в цвете)

[Figure 5. The distribution of longitudinal (top) and transverse (bottom) displacements in the case of a rigid right-hand boundary calculated until the moment of return of the reflected wave front Ei to the left boundary (color online)]

0.3

0.2

v{z)

0.1 0

8 • 10-5 4 • 10-5 0 ” 0.012 0.008 0.004 0

«Ф й(2)

Рис. 6. Изменение мгновенной формы поверхности оболочки за счет продольных и поперечных перемещений в случае жесткого закрепления правой границы при IV = 1 (слева)

и N = 2 (справа)

[Figure 6. Change in the instantaneous shape of the shell surface due to longitudinal and transverse displacements in the case of a rigid right-hand boundary when IS! = 1 (left)

and N = 2 (right)]

й(2)

Рис. 7. Изменение остаточной формы поверхности оболочки за счет продольных и поперечных перемещений в случае жесткого закрепления правой границы при снятии нагрузки в момент времени JV = 1 (слева) и JV = 2 (справа)

[Figure 7. The change in the residual shape of the shell surface due to longitudinal and transverse displacements in the case a rigid right-hand boundary when the load is removed at the time JV = 1 (left) and JV = 2 (right)]

u 0.018 0.016 0.014 0.012 0.010 0.008 0.006 0.004 0.002 0

v

0.018

0.016

0.014

0.012

0.010

0.008

0.06

0.004

0.002

0

-0.002

Рис. 8. Распределение продольных (сверху) и поперечных (снизу) перемещений в случае свободной правой границы, рассчитанное до момента времени возвращения отраженного волнового фронта Ei к левой границе (онайн в цвете)

[Figure 8. The distribution of longitudinal (top) and transverse (bottom) displacements in the case of a free right-hand boundary calculated until the moment of return of the reflected wave front Ei to the left boundary (color online)]

8 • 1(Г5 4 • 1(Г5

u(z)

Рис. 9. Изменение мгновенной формы поверхности оболочки за счет продольных и поперечных перемещений в случае свободной правой границы при Ж = 1 (слева)

и Ж = 2 (справа)

[Figure 9. Change in the instantaneous shape of the shell surface due to longitudinal and transverse displacements in the case of a free right-hand boundary when 14 = 1 (left)

-4•10“6

-2■10“6 u{z)

Рис. 10. Изменение остаточной формы поверхности оболочки за счет продольных и поперечных перемещений в случае свободной правой границы при снятии нагрузки в момент времени 14 = 1 (слева) и 14 = 2 (справа)

[Figure 10. The change in the residual shape of the shell surface due to longitudinal and transverse displacements in the case a free right-hand boundary when the load is removed at the time N = 1 (left) and N = 2 (right)]

случая свободного правого торца оболочки.

На рис. 10 представлены изменения остаточной формы поверхности оболочки вращения при снятии нагрузки в момент времени, когда волна Si дошла до правого торца оболочки, и в момент времени, когда волна Si вернулась к левому торцу оболочки, для случая свободного правого торца оболоч-

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

фронтом волны, использовано для построения численного алгоритма расчета деформированного состояния оболочки за фронтами продольной и сдвиговой волн в процессе их распространения. Условия на свободном и жестко закрепленном торце оболочки приведены в терминах скачков градиентов первого, второго и третьего порядков от перемещений на торце. Волны, отраженные от свободной границы волны являются волнами разгрузки, а отраженные от защемленного торца волны — волнами пластического нагружения. Сами поля деформаций и перемещений носят кусочный характер в плоскости (N,z), что связано с распространением продольных и сдвиговых волн в оболочке при ударе. Программа для ЭВМ требует времени порядка одной минуты для расчета кинематических параметров оболочки.

Конкурирующие интересы. Конкурирующих интересов не имеем.

Авторский вклад и ответственность. Все авторы принимали участие в разработке концепции статьи и в написании рукописи. Авторы несут полную ответственность за предоставление окончательной рукописи в печать. Окончательная версия рукописи была одобрена всеми авторами.

Финансирование. Исследование выполнялось без финансирования. Библиографический список

1. Timoshenko S., Woinowsky-Krieger S. Theory of Plates and Shells. New York: McGraw-Hill Book Comp., 1959. 580+xiv pp.

2. Сагомонян А. Я. Волны напряжения в сплошных средах. М.: МГУ, 1985. 416 с.

3. Ивлев Д. Д. Теория идеальной пластичности. М.: Наука, 1966. 332 с.

4. Егоров М. В. Динамическое деформирование осесимметричной оболочки вращения из упруговязкопластического материала вблизи ударных волн // Вестник Чувашского государственного педагогического университета им. И. Я. Яковлева. Сер. Механика предельного состояния, 2016. №2. С. 132-143.

5. Вервейко Н. Д. Лучевая теория упруговязкопластических волн и волн гидроудара. Воронеж: ВГУ, 1997. 204 с.

6. Марчук Г. И. Методы вычислительной математики. М.: Наука, 1989. 456 с.

7. Седов Л. И. Методы подобия и размерности в механике. М.: Наука, 1977. 440 с.

8. Самарский А. А., Гулин А. В. Численные методы. М.: Наука, 1989. 429 с.

9. Буренин А. А., Севастьянов Г. М., Штука В. И. Лучевой метод в приближённом решении задачи об ударном нагружении несжимаемого цилиндрического слоя // Вычислительная механика сплошных сред, 2016. Т. 9, №4. С. 400-411. doi: 10.7242/1999-6691/ 2016.9.4.33.

10. Герасименко Е. А., Завертан А. В., Рагозина В. Е. Об использовании прифронтовых лучевых разложений в динамике деформирования // ПММ, 2009. Т. 73, №2. С. 282-288.

11. Штука В. И. Применение лучевого метода в задаче определения напряженнодеформированного состояния предварительно продеформированного упругого слоя // Ученые записки КнАГТУ, 2017. №2. С. 40-44.

12. Бабичева Л. А., Быковцев Г. И., Вервейко Н. Д. Лучевой метод решения динамических задач в упруговязкопластических телах// ПММ, 1973. Т. 37, №1. С. 145-155.

13. Егоров М. В. Программа расчета динамического напряженно-деформируемого состояния оболочки вращения из упруговязкопластического материала при ударе ее по торцу: Свидетельство о государственной регистрации программы для ЭВМ. Зарегистрировано в Реестре программ для ЭВМ № 2017660609. 22.09.2017, 2017.

Vestn. Samar. Gos. Tekhn. Univ., Ser. Fiz.-Mat. Nauki

[J. Samara State Tech. Univ., Ser. Phys. Math. Sci.], 2018, vol. 22, no. 2, pp.325—343 ISSN: 2310-7081 (online), 1991-8615 (print) ^)http://doi.org/10.14498/vsgtu1610

MSC: 74H15, 65Z05

Mathematical modeling of dynamic deformation of elasto-viscoplastic shells of finite length by a ray method

N. D. Verveiko, M. V. Egorov

Voronezh State University,

1, Universitetskaya pl., Voronezh, 394018, Russian Federation.

Abstract

The paper presents a mathematical modeling of dynamic stress-strain state of the rotation of the shell elasto-viscoplastic material. We solve the modified system of S. P. Timoschenko’s partial differential equations by constructing a system of equations on moving surfaces of the gap with the initial conditions in the form of a shock at the end, written in the form of a power series in time, whose coefficients have initial conditions for the differential equations. The solution is presented in the form of the Taylor’s row series up to the fourth order in the shell coordinate. To simulate the waves reflected from the boundaries, the conditions at the boundary of two types (rigidly restrained and stress-free), independent of time, are introduced. A set of programs written in Fortran 90 on the Code::Blocks platform is developed. Two programs for the simulation of dynamic deformation of shell in elastic and elasto-viscoplastic state are implemented. We use the difference representation of the derivatives, the calculation of the integrals by the trapezoid method with a given step of partitioning the segment. The result of the programs is the grid functions of the coefficients of the Taylor rows, which are used to construct the displacement graphs as functions of time and the longitudinal coordinate of the shell.

Keywords: dynamic deformation, rotating shell, ray method, reflected wave, modeling, elasticity, viscosity, plasticity.

Received: 23rd February, 2018 / Revised: 21st May, 2018 /

Accepted: 11th June, 2018 / First online: 1st July, 2018

Competing interests. We have no competing interests.

Authors’ contributions and responsibilities. Each author has participated in the article concept development and in the manuscript writing. The authors are absolutely

Research Article

3 ©® The content is published under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.Org/licenses/by/4.0/)

Please cite this article in press as:

Verveiko N. D., Egorov M. V. Mathematical modeling of dynamic deformation of elasto-viscoplastic shells of finite length by a ray method, Vestn. Samar. Gos. Tekhn. Univ., Ser. Fiz.-Mat. Nauki [J. Samara State Tech. Univ., Ser. Phys. Math. Sci.], 2018, vol. 22, no. 2, pp. 325-343. doi: 10.14498/vsgtu1610 (In Russian).

Authors’ Details:

Nikolai D. Verveiko © http://orcid.org/0000- 0002- 6980- 9217 Dr. Techn. Sci.; Professor; Dept. of Mechanics & Computer Modeling; e-mail: verveyko2017@yandex.ru

Mikhail V. Egorov & © http://orcid.org/0000-0002-8977-7320

Postgraduate Student; Dept. of Mechanics & Computer Modeling; e-mail: egorovmv89@mail.ru

responsible for submitting the final manuscript in print. Each author has approved the

final version of manuscript.

Funding. The research has not had any funding.

References

1. Timoshenko S., Woinowsky-Krieger S. Theory of Plates and Shells. New York, McGraw-Hill Book Comp., 1959, 580+xiv pp.

2. Sagomonian A. Ya. Volny napriazheniia v sploshnykh sredakh [Stress Waves in Continuous Media]. Moscow, Moscow State Univ., 1985, 416 pp. (In Russian)

3. Ivlev D. D. Teoriia ideal’noi plastichnosti [Theory of Ideal Plasticity]. Moscow, Nauka, 1966, 332 pp. (In Russian)

4. Egorov M. V. Dynamic deformation of an axisymmetric shell of rotation from an elasto-viscoplastic material near shock waves, Vestnik Chuvashskogo gosudarstvennogo pedagogich-eskogo universiteta im. I. Ia. Iakovleva. Ser. Mekhanika predel’nogo sostoianiia, 2016, no. 2, pp. 132-143 (In Russian).

5. Verveiko N. D. Luchevaia teoriia uprugoviazkoplasticheskikh voln i voln gidroudara [Ray Theory of Elasto-Viscoplastic Waves and Hydrodynamic Waves]. Voronezh, Voronezh State Univ., 1997, 204 pp. (In Russian)

6. Marchuk G. I. Methods of Numerical Mathematics, Application of Mathematics, vol. 2. Springer-Verlag, New York, 1982, xiii+510 pp. doi: 10.1007/978-1-4613-8150-1.

7. Sedov L. I. Similarity and Dimensional Methods in Mechanics. New York, London, Academic Press, 1959, 363 pp. doi: 10.1016/c2013-0-08173-x.

8. Samarskii A. A., Gulin A. V. Chislennye metody [Numerical Methods]. Moscow, Nauka, 1989, 429 pp. (In Russian)

9. Burenin A. A., Sevastyanov G. M., Shtuka V. I. On the localization of discontinuities in calculations of incompressible elastic media, Computational Continuum Mechanics, 2016, vol. 9, no. 4, pp. 400-411 (In Russian). doi: 10.7242/1999-6691/2016.9.4.33.

10. Gerasimenko Ye. A., Zavertan A. V., Ragozina V. Ye. The use of near-front ray expansions in deformation dynamics, J. App. Math. Mech., 2009, vol. 73, no. 2, pp. 204-208. doi: 10. 1016/j.jappmathmech.2009.04.010.

11. Shtuka V. I. The ray metod application for stain-stress state of elastic cylindrical layer with preliminary deformations definition, Uchenye zapiski KnAGTU, 2017, no. 2, pp. 40-44 (In Russian).

12. Babicheva L. A., Bykovtsev G. I., Verveiko N. D. Ray method of solving dynamic problems in elastic-viscoplastic media, J. Appl. Math. Mech., 1973, vol. 37, no. 1, pp. 132-141. doi: 10. 1016/0021-8928(73)90142-1.

13. Egorov M. V. A program for numerical calculation of dynamic stress-strain state of the shell of revolution made of elasto-viscoplastic material when strikes the butt end: Certificate of state registration of the computer program, Registered no. 2017660609. 22.09.2017, 2017 (In Russian).