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

Библиографический список

1. Хорошилов В.С. Механические модели движения космического аппарата с солнечной батареей // Изв. АН СССР. МТТ. 1978. № 5. С. 18-24.

2. Докучаев Л.В., Климов О.П. Об устойчивости вращения твердого тела с гибкими элементами // Изв. АН СССР. МТТ. 1982. № 5. С. 10-15.

3. Набиуллин М.К. Стационарное движение и устойчивость упругих спутников. Новосибирск.: Наука, Сибир. отд-ние, 1990. 216 с.

4. Мирович Л., Квинн Р.Д. Уравнение движения маневрирующего космического аппарата нежесткой конструкции // Аэрокосмическая техника. 1988. № 6. С. 82-96.

5. Ганиев Р.Ф., Ковальчук П.С. Динамика систем твер-

дых и упругих тел. Резонансные явления при нелинейных колебаниях. М.: Машиностроение, 1980. 208 с.

6. Борисов М.В. Применение метода Релея - Ритца для нахождения собственных частот и форм колебаний сложной упругой системы // Студенческая наука аэрокосмическому комплексу: Сб. тр. студ. и асп. фак. летательных аппаратов. Самара, 2001. Вып. 7. С. 10-16.

7. Ланцош К. Вариационные принципы механики. М.: Мир, 1965. 408 с.

8. Тимошенко С.П., Янг Д.Х., Уивер У. Колебания в инженерном деле. М.: Машиностроение, 1985. 472 с.

9. Круглов Г.Е. Аналитическое проектирование механических систем / Самар. гос. аэрокосм. ун-т им. акад. им. С.П. Королева. Самара, 2001. 132 с.

УДК 539.3

О НЕКОТОРЫХ ОСОБЕННОСТЯХ ПРИМЕНЕНИЯ МЕТОДА КОНЕЧНЫХ ЭЛЕМЕНТОВ К РЕШЕНИЮ КОНТАКТНОЙ ЗАДАЧИ НА БАЗЕ ПРОГРАММНОГО КОМПЛЕКСА ABAQUS

Н.В. Овчинникова, Ю.В. Чеботаревский*

Саратовский государственный технический университет, кафедра теоретической механики E-mail: [email protected], *[email protected]

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

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

Some Aspects of Applying Finite Element Method to Contact Problems Using ABAQUS System

N.V. Ovchinnikova, U.V. Chebotarevsky*

Saratov State Technical University,

Chair of Theoretical Mechanics E-mail: [email protected], *[email protected]

In this paper the selection of preferred computational scheme for numerical research of a contact interaction between a rigid indentor and an elastic-plastic half-space using ABAQUS system is described. Furthermore, the critical values of the load applied to indentor which don’t cause plastic deformation in the material and subsequently don't effect the process of a surface hardening are determined. Some peculiarities of the forming of near-surface plastic region by static loading are investigated.

Keywords: contact interaction, deflected mode, isotropic-kinematic hardening model, plastic deformation, surface hardening, ABAQUS.

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

© Н.В. Овчинникова, Ю.В. Чеботаревский, 2009

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

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

В связи с изложенным целью настоящей работы являются:

а) выбор приемлемой схемы численного решения задачи, сформулированной в [2], на базе программного комплекса ABAQUS путем сравнения результатов расчета с известным аналитическим решением;

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

При выборе схемы численного решения рассматриваемой задачи ограничимся следующим частным случаем. Обрабатываемый объект в силу локальности исследуемого механического процесса представим в виде полупространства, к некоторой точке которого прикладывается с постоянно действующей силой F0 = const абсолютно твердый индентор со сферической рабочей поверхностью радиуса кривизны RKр (рис. 1).

При условии, что поверхность полупространства вне зоны контакта свободна от внешней нагрузки, рассматриваемая осесимметричная задача будет отличаться от сформулированной в работе [2] лишь несколько измененными граничными условиями, а именно

aij ^ 0 (i, j = x,y,z) при x,y,z ^ го,

= -р(х, у, Ь) при г = 0,

О

р(х, у,1) ¿хйу = ^0.

0 0

Здесь Сту — компоненты тензора напряжений, р(х,у,£) — давление в зоне контакта.

В качестве эталонной для рассматриваемого частного случая при упругом поведении материала (до наступления пластических деформаций) примем решенную А. Динником [3] задачу о сжатии двух соприкасающихся, ограниченных идеально гладкими поверхностями тел, в предположении, что одно из них является деформируемым, а другое — абсолютно твердым. Положим в [2] радиусы кривизны абсолютно твердого тела равными радиусу кривизны сферической рабочей поверхности индентора Ли = Л12 = Лкр, а деформируемого тела — равными бесконечности Л21 = Л22 = го. Для случая кругового контура давления в силу осевой симметрии ограничимся записью выражений для напряжений в плоскости х = 0:

Рис. 1. Схема процесса взаимодействия материала с индентором

3Fq Г (1 — 2v) a2

2na2

1 -

VA3

Z Л 1 + V Л 1 + V

+—^ ( 2v +—2—;—г A +------

V A V a + A a

VA ■ arctg —A )

3

2

У

°УУ —

3^о

2па2

Известия Саратовского университета. 2009. Т.9. Сер. Математика. Механика. Информатика, вып.2 (1 - 2v) а2

3

3^о

У2

1

+

а2 г3

л/Л3/ (Л2 + а2г2) л/Л л/Л I а2 + А

. 1 — V л , л/Л а ^

+----^ | 2 . л Л + (1 + V)-------------агс^ ^ - 2

3^о

/Г ^ — -

2п (Л2 + а2г2)л/Л’ У2 2п (Л2 + а2г2) (а2 + Л)’

Здесь — численное значение главного вектора давлений, равное по величине прикладываемой к индентору силе, V и Е — соответственно коэффициент Пуассона и модуль Юнга для материала полупространства, а радиус контура давления а и параметр Л находятся из соотношений:

22

Уг

1л/Л

’ч/Л

— 0.

а—

3(1 - V2) Я

кр

— 1.

4Е а

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

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

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

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

Расчеты проводились при различных значениях коэффициента сгущения для двух типов конечных элементов, предлагаемых программным комплексом ABAQUS: треугольных и четырехугольных. В качестве примера на рис. 2 показаны зависимости численных значений осевого напряжения а22 в точке

контакта от величины коэффициента сгущения сетки при воздействии на индентор силы в 9,81 Н.

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

-з.з

-3.35

-3.4

-3.45

-3.5

-3.55

-3.6

™1" "решение на базе треугольной сетки аналитическое решение А. Динника ““ “решение на базе четырехугольной сетки

чч ^1 ^ _

1 \ \ \ / Ч, 1**, ■ ч : >

\' ...у.. , у/

* V Л

\ \ ъ--‘ / \ / « / / ч ч ч V \ \ \ V \

% % / ( ч ч ч V

>

100

200

300

400

500

600

К,

Рис. 2. Зависимость значений осевого напряжения аг фициента сгущения сетки

от величины коэф-

3

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

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

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

F = f (а, а) — а0 = 0,

где

а0 = а |о + (l - e-b#pi)

/ (а, а) =

\

3 3 3

2 ЕЕS.jSy = аи,

i=1j=1

(l)

п|° и п° — значения предела текучести в момент появления пластической деформации и в процессе дальнейшего нагружения, Б — девиатор тензора напряжений, , Ь, С, 7 — параметры материала, определяемые на основе экспериментальных данных.

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

1.8

1.6

1.4

1.2

1

0.8

0.6

0.4

0.2

0

жении СТи вдоль оси z при =9, 81 Ни Якр = 0,06 м.

Из полученных результатов следует, что в пределах упругости интенсивность напряжений сти достигает максимального значения не в области контакта, а на некотором расстоянии от неё под центром давления. В то время как осевое напряжение имеет наибольшее значение на поверхности обрабатываемого объекта в центре контура давления. На рис. 4 (для тех же значений постоянной составляющей силы и радиуса кривизны рабочей поверхности индентора) наглядно показано, что даже при превышении осевым напряжением предела текучести по модулю в 1,5 раза материал в центре контура давления остается упругим.

Jb г% _ 1 — . ■ численное решение ■^^"аналитическое решение А. Динника

і

6 2 х 101 М

1 2 3 4 5

Рис. 3. Распределение интенсивности напряжений au по оси симметрии

<5 х 10 , П.1

3.5

2.5

1.5

0.5

0

“ ™ ™ интенсивность н —1 —1 ■ осевое напряже ап ряже ний С7И ние — и

1 * 1 / 1 предел 1 еку чес 1

'ц-

! ч V ^ ^ ■■

► X Ч ч 1

Г

О

0.2

0.4

0.6

0.3

1 2Х104,М

Рис. 4. Распределение интенсивности напряжений и осевого напряжения

по оси симметрии

Отсюда следует, что впервые пластические деформации возникают не на поверхности контакта, а на некотором расстоянии от неё под центром давления в точках наибольших значений интенсивности напряжений. На рис. 5 показана зависимость наибольших значений интенсивности напряжений от радиуса кривизны индентора для различных значений прилагаемой к нему постоянно действующей силы. Легко заметить, что приведенные графики не являются гладкими, а имеют точку излома при значениях интенсивности напряжений совпадающих с пределом текучести материала. Её наличие является следствием не только количественного, но и качественного изменения характера поведения напряжений и деформаций при переходе материала из упругого состояния в пластическое. Точка излома соответствует минимальным значениям радиуса кривизны индентора, при которых под действием данной силы в обрабатываемом объекте не возникают пластические деформации.

сих 10, Па

з I

2.8

2.6

2.4 2.2

2

1.8

1.6

1.4 1.2

1

™ ™ ™ решен —1 —1 ■ решен ■«»**** решен ш пластической задачи при нагрузке 9,81 Н ле пластической задачи при нагрузке 19,62 Н ле пластической задачи при нагрузке 29,43 Н текучести

V \ч Г1рнД 0.1

ч

ч > ч*-.

^ ^ ™ 1 М 1

X ^ V

N. (

0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1 1Ц, м

Рис. 5. Зависимость наибольших значений интенсивности напряжений от радиуса кривизны индентора при трех значениях постоянной силы

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

35 Р0„Н

Рис. 6. Зависимости интенсивности напряжений от величины постоянно действующей силы

Результаты исследования напряженно-деформированного состояния материала при наличии пластических деформаций, показали, что независимо от уровня нагружения в пределах ^ = (9,81 — 39, 24) Н в зоне контакта остается приповерхностный тонкий упругий слой толщиной от 12 до 37 микрон. В качестве примера на рис. 7 показано распределение остаточных деформаций по толщине обрабатываемого объекта при полном снятии нагрузки для трех значений силы, прикладываемой к индентору с радиусом кривизны рабочей поверхности 0,02 м.

}Ш0

-0.5

-1.5

-2.5

ч ✓ ✓ Г г1-— 4 * * * к,-'* -гг"-

V 4 / .... ф. / / ✓ / % ■ ■ ■

- > £\ \ч ■■ ч \ : /... 4 4 / Р *

/ у ч +

— — — решение пластической задачи при нагрузке 9,81 Н решение пластической задачи при ’"нагрузке 19,62 Н решение пластической задачи при нагрузке 29,43 Н

Ч-** *•*

О 0.2 0.4 0.6 о.!

1 1.2 1.4 1.6 1.Е

2.2 2.4

Ъ х 10 , М

Рис. 7. Распределение остаточных деформаций по толщине обрабатываемого объекта при полном снятии нагрузки

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

Библиографический список

1. Галин Л.А. Контактные задачи теории упругости и вязкоупругости М.: Наука, 1980. 304 с.

2. Овчинникова Н.В., Павлов Д.Г., Чеботарев-ский Ю.В. Модельная задача для исследования процессов поверхностного упрочнения пластическим деформированием с применением ультразвуковых воздей-

ствий // Вестн. СГТУ. 2007. № 4(28), вып. 1. С. 14-18.

3. Динник А.Н. Избранные труды. Т. I. Киев: Изд-во АН Укр. ССР, 1952. 152 с.

4. ABAQUS Analysis User’s Manual Version 6.4 / Hibbitt: Karlsson Sorensen, Inc. USA, 2002.

УДК 532.5; 532.135

ТОРМОЖЕНИЕ ПЛАСТИНЫ О СЛОЙ «ЗАПАЗДЫВАЮЩЕЙ» ВЯЗКОПЛАСТИЧНОЙ СРЕДЫ С УЧЕТОМ ПРИСТЕННОГО СКОЛЬЖЕНИЯ

М.И. Сафрончик

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

кафедра математической кибернетики и компьютерных наук

E-mail: [email protected]

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

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

Plate Braking Against the Layer of «Delayed» Viscoplastic Fluid with Regard to Wall Sliding

M.I. Safronchik

Saratov State University,

Chair of Mathematical Cybernetics and Computer Sciences E-mail: [email protected]

The paper presents the problem of unstable viscoplastic fluid flow between parallel planes, one of which is fixed, while the other one is put in motion from a standstill under the influence of constant force. Viscoplastic fluid flow develops gradually. The border of the flow is not known in advance and is to be determined in the process of solving the task. The force applied to the upper plate is chosen so as to cause the effect of sliding along the two plates in the course of time. The task definition is given within the limits of five-parameter model, which permits to take up the difference between behavior under stress and without stress as well as possible sliding along the solid walls. Hysteresis of deformation is considered by means of Slibar - Pasly hypothesis. To take the possible sliding along the walls into account, a hypothesis, analogical to the well-known hypothesis for viscous fluid of Prof. N.P. Petrov, is suggested. The offered hypothesis also allows to describe the natural physical condition of smooth transition from «sticking» to «sliding». Moreover, the parameters included into it can be defined empirically. To solve the task with the required border, a modified method of Kolodner is used.

Keywords: viscoplastic fluid, unstable viscoplasticfluid flow, Slibar— — Pasly hypothesis, wall sliding, the method of Kolodner.

1. ПОСТАНОВКА ЗАДАЧИ. ОСНОВНЫЕ УРАВНЕНИЯ И СООТНОШЕНИЯ

Рассмотрим две неограниченные горизонтальные пластины, находящиеся на расстоянии Н друг от друга, пространство между которыми заполнено вязкопластичным материалом. Нижняя пластина остается все время в покое, а верхняя начинает движение из состояния покоя под действием приложенной к ней в момент £ = 0 постоянной силы. Обозначим величину силы, приходящейся на единицу площади пластины через дь а массу единицы площади пластины толщины Н через т. В этом случае уравнение движения пластины можно записать в виде

т— = дг — т (0,£), (1)

© М.И. Сафрончик, 2009