О напряженно-деформированном состоянии цилиндра из упругопластического материала под действием давления в канале Текст научной статьи по специальности «Физика»

УДК 531 А.В. Горшков

О НАПРЯЖЕННО-ДЕФОРМИРОВАННОМ СОСТОЯНИИ ЦИЛИНДРА ИЗ УПРУГО-ПЛАСТИЧЕСКОГО МАТЕРИАЛА ПОД ДЕЙСТВИЕМ ДАВЛЕНИЯ В КАНАЛЕ

С помощью вариационных принципов исследуется напряженно деформированное состояние цилиндра из упруго-вязкопластического материала при высоком давлении в канале. Эта модель описывает один из завершающих этапов пробоя цилиндра высоковольтным электрическим разрядом. Определены напряжения в материале и скорость течения материала вблизи стенок канала.

При решении задач механики сплошной среды одна из важных проблем - получить конечномерные уравнения, адекватные непрерывным. Один из путей решения этой проблемы - использование вариационных принципов или законов сохранения. Если для непрерывной модели выполняется вариационный принцип или закон сохранения, то это свойство должно выполнятся и для дискретной модели. В качестве примера использования вариационного принципа приведено решение динамической осесимметричной задачи о напряженно деформированном состоянии цилиндра из упруго-вязкопластического материала с тонким каналом по оси (рис. 1). Внешний радиус цилиндра Я0 = 1, радиус канала Я1 = 0.1, высота цилиндра И = 1. В канале действует высокое быстро нарастающее давление Р(/). Внешняя поверхность свободная. Силой тяжести будем пренебрегать. На рис. 1 стрелками А и В отмечены точки закрепления, ООх - ось симметрии. Для решения используется вариационный принцип скоростей и напряжений, развитый в работах В.Л. Колмогорова [1, 2] в сочетании с методом конечных элементов [3].

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

Используется Лагранжева начальная цилиндрическая система координат. Через х' будем

обозначать текущие координаты точек среды, а через у' - их начальные значения, ' = 1,2. Ось

Р и с. 1.

оу1 направлена по радиусу, ось ву2 касательной,

совпадает с осью цилиндра, ось оу направлена по

Яі < у1 < Я0, 0 < у2 < И . (1).

Осесимметричное движение в цилиндрической системе координат описывается уравне-

ниями:

Эа11 Эа1

Эу1 " да12

Эу2 у1 + и

да22

1 - а33 (у1 + и‘ ) = Ро

Эу1 Эу2 у1 + и

1 ро

и

Э і

Здесь ау - компоненты тензора напряжений, иг- - компоненты вектора скорости, и{ - компоненты вектора перемещений. С учетом осевой симметрии компоненты тензора скоростей деформаций ву примут вид:

е11 =

Эи1

э7

1

Эи1 Эи2

эу

е22 =

Эи2

эр

=у + и1 у

е12 =

е

2

Материал цилиндра описывается следующими определяющими соотношениями: при T < оs -ау = leg4 + 2|mey ; при T > оs материал вязкий идеально-пластический - T = os +gH . Объемная

деформация остается упругой. Здесь Т = ^sijS,J I2 - интенсивность касательных напряжений,

= а" _ £} з - компоненты девиатора тензора напряжений, О - первый инвариант тензора напряжений,

I 21

H=

| з [(11 е22 ) +(е11 e33 ) + (33 e22 ) ] + е12 + e13 + e23 ^

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

формации, £ = eijg' - первый инвариант тензора деформации, g,J - компоненты метрического

тензора, в рассматриваемом случае giJ = 0, ' ф j , g11 = g22 = 1, g33 = (у1 + и1)2, 1, т - постоянные Ламе, g - коэффициент трения.

Так как материал цилиндра вязкопластический, то учитывается нелинейная связь перемещений и деформаций.

Функционал вариационного принципа имеет вид [1]:

где I1 = J

V

S {I1 +12 +13 +14 } = 0,

© '

еу S!J e s' f

J s‘J (e)de + J ev (s)ds + Js(e)de + Je(s)ds + p(w! - g*)

(2)

dV; 12 =-J fVdS;

13 = -J Z^’dS; 14 = -J

f!

dS !J Эи1 Эи2 u1

dS e = e gJ =_______________ +________ +__________

’ j ^,.1 ^,.2 ,.1 , ..1

Э y1 Э y2 y1 + u'

/«;_/ ) _}«,(/)#

0 0

- первый инвариант тензора скоростей деформации. На части поверхности тела Б/ заданы силы /', на - скорости, на Б^ - трение скольжения. Звездочкой отмечены заданные величины, штрихом - варьируемые.

В рассматриваемом случае поверхность Б/ состоит из двух участков:

1) поверхность канала, на которой задано давление Р(() - известная функция времени;

2) поверхность торцов цилиндра, на которых давление равно давлению окружающей среды Р0 ; Б^ - внешняя цилиндрическая поверхность, на которой скорость точек равна нулю. Тогда /4 = 0, 1з = 0,

п 140

2pR1 P(t)Ju'ldy2 - 2pP0 J yU/2dy1 .

Принцип скоростей и напряжений для изотропного тела в данном случае принимает вид

i=J

П lev

J T(H)d H + J H(T)d T + J s (e )de + J e (s ) ds + p w‘vi

0 t s 0 0

h R0

-2 p R1 P(t)J u' dy2 - 2pP0 J y1v'1 dy1 .

dV -

(3)

Для построения конечномерной модели используется метод конечных элементов. В области (1) вводится N узлов с координатами (уП,у2), которыми она разбивается на М треугольных

V

0

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

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

На элементе задается система локальных базисных функций

ф, = а у1 + ьг у2 + с,, і =1,2,3.

Коэффициенты аі, Ь, с. зависят от координат вершин элемента и выбираются так, чтобы система функций ф, была нормализованной. В локальной системе координат базисные функции имеют вид

Ф,(г‘Л = -у - +1, ф,(Лг!) = у - ^, фзС-У) = -^,

I щ I щ д

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

Поле перемещений на элементе аппроксимируется локальными базисными функциями:

и(т)і»и Чу).

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

Соответственно, глобальные перемещения аппроксимируются выражениями:

Vі »Хлп гіт)(у)икі =Фк('у)икі,

т

Vі - і-ая компонента вектора перемещения, икг - і-ая компонента перемещения узла с гло-

. М (т) . (т)

бальным номером к, Фк (у) = ^Апк фПт)(у) - глобальная базисная функция. Таблица Лпк задает

т=1

соответствие глобальных номеров узлов сетки и локальных номеров узлов элементов:

к 11, если глобальный узел к связанной модели совпадает с узлом п элемента т;

Л п = -

(т)

[0, в противном случае.

Аналогично аппроксимируются скорости и ускорения:

(т) (т)

V » ХЛ Фпт)(^)ик = фк(^у*, ж1» ^лпк фПт)(^)"й = фк(^)жк.

тт

Обратная процедура - разбиение области на конечные элементы, выполняется с помощью

(т) (т) (т)

таблицы ОП . Таблицы ОП получаются из Лкп транспонированием при фиксированном т .

После подстановки конечно-элементных представлений в (3) и интегрирования, последнее слагаемое в скобке примет вид

1к тт

где

Уя'ФкУ ап ак,„ ...

і J р я "У

т=1

т(п1. = | Ро(У )ёу (У )№ (У ^ (У № - согласованная матрица масс, V" - варьируе-

мая I -ая компонента скорости в узле п, ^ - ] -ая компонента ускорения в узле Сумму

М (т) (т)

Мр„ = £пр цХЗ

т=1

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

Аппроксимация скоростей деформаций на элементе имеет вид

(т) і (т)

е .. =1ак і 2

VI

Эфкт^ + уп Эфкт)

ЭУ.

ду і

Из полученных выражений видно, что на элементе скорость деформаций не зависит от координат, и вариация первого слагаемого в скобке в (3) по скоростям имеет вид

И' м

8| | Т (И )dHdV = У ( + уИ(т) )(т)

Б,

(т)

(т) ЭН(т) Эе 8Н(т) = У а"АН— ^ 8V

,^р " Э. дVpn

т=1

ЭИ(т) р ’ Эе (т)

4е(т)

3И(т)

ЭН

(т)

(т)

Эе

(т)

Н

(т)

-1 * і

Эе

(т)

ЭV"

0, р * і, р * і,

Э^Г

ЭУ1

Эркт)

ЭУ,

р = 1

р = і,

где т) - площадь элемента т, И(т) - интенсивность скоростей деформаций на элементе т, VП - п компонента вектора скорости в узле с глобальным номером р. Функцию ,И(т) будем

1 1

аппроксимировать выражением

И

(т)

Ф

к И (т)к

Аналогично строится вариация интенсивности напряжений Т .

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

dUі

dt

= V.

М.

пку

dt

а,и V),

которая решается численно. Здесь і - номер компоненты вектора перемещений или скорости, к - номер глобального узла. Используется неявная разностная схема:

ип.+, - ип. =—Vп.+, + Vп.),

і+і і2 V і+і і'

м.„, .+, - V* )=— V V+і.+1.VJ+,)+/- V и ,vJ)).

Индекс . - номер временного слоя; и. показывает, что функция может зависеть от значений перемещений и скорости во всех узлах сетки для данного временного слоя. Эта схема для ли-

(е)

V 0

нейной механической задачи обеспечивает сохранение полной энергии системы на шаге времени.

На рис. 2 и 3 представлены графики изменения интенсивности скорости сдвига в среднем сечении цилиндра у2 = 0.5 и интенсивности касательных напряжений. Из графиков видно, что вдоль радиуса распространяется волна сдвига. Расчеты проводились для цилиндра с размерами Я0 = 1, Яг = 0.1, к = 1; модуль упругости Е = 1011; коэффициент Пуассона V = 0.3 ; о= 107;

скорость роста давления = 5' 107.

Р и с. 2.

Р и с. 3.

На рис. 4 и 5 представлены профили интенсивности касательных напряжений для моментов времени t = 0.0003 и t = 0.0007 Увеличение напряжений в точках с координатами к = 0, г = 1. и к = 1, г = 1 на рис. 5 объясняется тем, что это точки закрепления.

ь= 0.0003 Р и с. 4.

ЗхЮ7

Ъ= 0.0007

Р и с. 5.

Предложенный алгоритм позволяет определить скорость выброса материала, график распределения которой представлен на рис. 6.

Рассматривалось несколько вариантов сетки. Первый - сетка равномерная по радиусу и сетка с увеличивающимся шагом по радиусу. Сравнивались также сетки с различным разбиением области на элементы. Сравнение результатов показало, что наилучшие результаты обеспечивает равномерная сетка. Менялось также расположение треугольников. Наилучший вариант соответствует симметричному расположению.

1. Колмогоров В.Л. Вершинин С.В., Спевак Л.Ф., Горшков А.В. Применение метода расчета напряженно- деформированного состояния для некоторых задач обработки металлов давлением // Тез. докл. 10-й зимней школы по механике сплошных сред, Пермь, 1995.

2. Колмогоров В.Л., Федотов В.П., Горшков А.В. Трехмерный анализ напряженно-деформированного состояния // Ковочно-штамповочное производство. №8. 1998. С.23-28

3. Оден Дж. Конечные элементы в нелинейной механике сплошных сред. М.: Мир, 1976. 466 с.

Работа выполнена при поддержке РФФИ, проект № 04-01-00274.

Поступила 26.12.2005 г.