Разработка алгоритма контактного взаимодействия на основе альтернирующего метода Шварца Текст научной статьи по специальности «Физика»

УДК 539.3

И. В. Станкевич, М. Е. Яковлев, Си Ту Хтет

РАЗРАБОТКА АЛГОРИТМА КОНТАКТНОГО ВЗАИМОДЕЙСТВИЯ НА ОСНОВЕ АЛЬТЕРНИРУЮЩЕГО МЕТОДА ШВАРЦА

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

E-mail: me-yakovlev@rambler.ru

Ключевые слова: контактное взаимодействие упругопластических тел, метод Шварца, задача Герца, метод конечных элементов, итерационное решение.

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

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

Актуальной является также проблема создания новых эффективных алгоритмов и на их основе современного прикладного программного обеспечения для решения контактных задач вычислительной механики. Ниже рассматривается алгоритм решения контактной задачи механики в геометрически сложных двумерных областях, базирующийся на альтернирующем методе Шварца [1, 2].

Математическая постановка контактной задачи теории упругости. Рассмотрим два двумерных однородных и изотропных линейно-упругих контактирующих тела А и В, занимающих на плоскости (в пространстве ^2) области 01 и 02 и ограниченных кусочно-гладкими границами 801 и 8С2. Введем на плоскости декартову систему координат 0х\х2. Математическая формулировка контактной задачи теории упругости в этом случае будет включать:

— уравнения равновесия

а]К] (и(х)) + Qг(x) = 0, х е 0а,г,3 = 1,2, а € {А, В}; (1)

— граничные условия (кинематические и силовые соответственно)

и(х) = <(х), х е Б* С 80а,а е {А, В}, (2)

а^(и)п = (х), х е Бра С 80а,г,з = ТДа е{А,В}; (3)

— соотношения Коши 1 2

— определяющие соотношения в форме закона Гука

а = Б : в; (5)

М

здесь Б — матрица Гука, а = {а} = < а22 > — вектор напряжений,

м

£ = {е} = < е22 > — вектор деформации, Qi — компоненты массовой (объемной) силы, и^(х) = | 1 — вектор заданных перемещений

ег](x) = -(Ui,j (x)+ Uj(x)) x G Ga,i,j = 1, 2,a G {A, B}; (4)

точек поверхностей Б^ и Б^, р«(х) — компоненты заданной распределенной нагрузки р«(х) = {на поверхностях Бр1 и БВ.

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

и1(х) - (х) = 8п (х) (6)

и по напряжениям (силовое условие)

аА(х} = -аВ(х) < 0, (7)

А В

где иАА и «В — проекции перемещении граничных точек на внешнюю нормаль к границе тела А; 5п — начальное расстояние по нормали

А В

между граничными точками; аА и аВ — составляющие напряжении по внешней нормали к границе тела А. Соотношения (6) и (7) соответствуют случаю, когда трение не учитывается.

Совокупность соотношений (1)-(7) составляет математическую формулировку контактной задачи теории упругости. Для ее решения в настоящей работе был использован алгоритм, основанный на альтернирующем методе Шварца.

Основные процедуры альтернирующего метода Шварца. Альтернирующий метод Шварца является итерационным методом. Его суть состоит в следующем. На первом шаге на контактных поверхностях тел А и В соответственно Б^ и Б^ задают начальные перемещения и(х) | а = и^(х) и и(х) \3в = и^ (х), которые имеют смысл дополнительных кинематических условий. Значения компонент векторов и^(х) и и^ (х) выбирают из диапазона ожидаемых значений для зоны контактного взаимодействия.

В начале решения контактные поверхности Бк1 и Б^ могут совпадать, но могут и не совпадать. Это зависит от особенности постановки решаемой контактной задачи. В конце решения контактной задачи

геометрически скорректированные поверхности контакта совпадают

Б А _ б В

Бк = Бк .

Далее решают независимо две подобные задачи теории упругости для тел А и В. Затем вычисляют поверхностные силы РА(х) и рВ (х) на контактных поверхностях БкА и БкВ и их корректируют так, чтобы выполнялись силовые контактные условия (7).

На втором шаге на контактных поверхностях Бк1 и БВ задают силовые контактные условия, в качестве которых используют скорректированные поверхностные силы РА(х) и рВ (х),

OK \sb = pf(x), x e Sa С ¿Ga, i,j = 1, 2,a e{A,B};

и вновь решают независимо задачи теории упругости отдельно для тел А и В. По результатам полученных решений выполняют коррекцию компонент векторов перемещений иА(х) и иВ (х) соответственно точек контактных поверхностей Бк1 и БкВ с тем, чтобы не было взаимного проникания контактирующих тел, т.е. выполнялись кинематические условия контактного взаимодействия (6). Следует отметить, что коррекция векторов перемещений иА(х) и иВ (х) влечет геометрическое

изменение поверхностей контакта БА и Б;, например, из рассмотрения исключаются те точки поверхностей контакта БА1 и , в которых отсутствуют сжимающие нормальные напряжения.

Скорректированные перемещения контактных поверхностей иА(х) и и; (х) рассматривают в качестве новых кинематических граничных условий на геометрически измененных поверхностях контакта Б^1 и

и вновь решают независимо задачи теории упругости отдельно для тел А и В.

Вопросы сходимости подобного типа алгоритмов рассмотрены в работах [1-3].

Для численного решения контактной задачи был использован метод конечных элементов. Пронумеруем узлы контактной поверхности Бк и введем в рассмотрение два вектора — {ик} и {Як}, первый из которых составлен из компонент перемещений и и V узлов, расположенных на контактной поверхности Бк, а второй — из компонент / и д узловых сил тех же узлов.

Коррекция компонент вектора {ик} а выполняется в соответствии с соотношением

>, 2 п

и

(A),m

=

2 n-1

+ а

(A),m

2n-1 (A),m

n = 0,

(A),m 2n 1

(B),s

2n 1

(A),m.

(8)

n = 1, 2 ...,

где а

2n-1

(A),m

— итерационный параметр, m (1 ^ m ^ MA) — узел, лежа-

щий на контактной поверхности SA тела A,

2n 1

(B),s

вектор пере-

мещений сходственной точки в, лежащей на контактной поверхности Б В тела В.

Для коррекции компонент вектора контактных узловых сил { Як } а используется формула

f

2n+1

(A),m

2n

f

2n

(A),m

2n f 2n

+

/

2n

(A),m

, n = 0,1, 2

(9)

где а(А) т — итерационный параметр, т (1 ^ т ^ МА) — узел, лежащий на контактной поверхности БА1 тела А, '

9

2n

вектор кон-

(B),s

тактных узловых сил сходственной точки в, лежащей на контактной поверхности тела В.

v

9

9

принадлежащего контактной поверхности БА тела А, с отрезком (стороной конечного элемента), лежащем на контактной поверхности $В

Чтобы воспользоваться формулами (8) и (9) необходимо знать ком-

{ч 2п—1 2п

> и вектора сил < > в сход-

(В),. ^ В (В),.

ственной точке в, лежащей на контактной поверхности Б^ тела В.

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

А

к

,октной поверхности Ск

тела В.

Для эффективной работы рассмотренного алгоритма важным является вычисление численных значений итерационных параметров аой и а^.т, а е {А, В}, 1 < т < Ма, здесь Ма - число контактных узлов тела а.

Пусть для определенности а = А и 1 ^ т ^ Ма. На четных итерациях, т.е. когда выполняется коррекция компонент вектора перемеще-

{ч 2п—1

и „ 2и — 1

> , итерационный параметр а^) т определяется с

помощью выражения

КА"1 II

а2п-1 =_11 (А)'т"__(10)

а(А),т = и, 2п-1 || + ||,,2п-1|| • V1";

||и(А),тН + |и(В),в 1

На нечетных итерациях, т.е. когда выполняетмся коррекция компо-

2п Ш 2П

нент вектора контактных узловых сил р^) т = < ^ > , итерационный параметр а^) т определяется с помощью выражения

11Р2п I

а2п = ИН(А)'тИ (11)

а(А)'т = IIр2п || + ||р2п || • (11)

Как показали результаты численных исследований, использование для вычислений итерационных коэффициентов соотношений (10) и (11) является достаточно эффективным.

В данной работе трение скольжения учитывалось по закону Кулона [4]. Условие кулоновского трения вызывает следующее дополнение к постановке задачи (1)-(7):

условия сцепления

uf - uB = 0;

И < ei^Ai; kBi < вКвi;

условия скольжения

iaAi = kB i = в kAi

(13)

где в — коэффициент трения; а^4, — компоненты вектора напряжений в направлении касательной к границе тела А.

Вначале вычислений с учетом трения все узлы контактной поверхности считаются сцепленными. С этими граничными (кинематическими) условиями выполняют начальную (нулевую) итерацию. После усреднения усилий в зоне контакта проверяют выполнение силовых условий (12). В случае их невыполнения переводят соответствующие узлы в число скользящих и ограничивают в них касательные усилия по (13). Далее итерационно уточняют контактные условия по формулам (8) и (9) с учетом соотношений (12) и (13), причем на нечетных итерациях условия трения удовлетворяются точно, а на четных итерациях приближенно. Решение задач показало, что подобный учет трения не добавляет контактных итераций. Кроме того, опыт показывает, что сходимость достигается после нечетных итераций, т.е. когда условия трения удовлетворены точно.

Результаты численных исследований. На основе разработанного алгоритма был создан комплекс прикладных программ для решения упругопластических контактных задач. Расчеты упруго пластического деформирования проводили по схеме метода переменных параметров упругости. Для проверки работоспособности программного комплекса было выполнено численное решение ряда задач, имеющих известное аналитическое решение, в том числе и задачи Герца. Рассматривалось контактирование длинного упругого цилиндра радиусом Я = 25 мм с упругим полупространством (рис. 1). Модули Юнга и коэффициенты Пуассона контактирующих тел были приняты равными соответственно Е1 = Е2 = 2,1 ■ 105 МПа и = ^2 = 0,3. Распределение давления по поверхности контакта определяли с помощью формулы

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

Рис. 1. Задача Герца

B = , Q (v1 + v2)

_ п (К 1 + К2)'' 4(1 - К 1

где . = , К = д -

Для численного решения были построены конечно-элементные модели контактирующих тел, представленные на рис. 1, а. Цилиндр состоял из 1546 конечных элементов, а пластина — из 1081 конечного элемента. Полученная область контакта изображена на рис. 1, б. Максимальная интенсивность напряжений 370 МПа.

В результате упругого расчета было получено значение полуширины поверхности контакта В = 0,189 мм, при этом аналитическое решение дает значение В = 0,185 мм. Сходимость достигнута за 11 итераций, а относительная погрешность выполнения контактных условий составила е = 0,01.

Был проведен анализ напряженно-деформированного состояния замковых соединений типа "ласточкин хвост" лопаток и дисков компрессоров ГТД с различными углами раствора (45°, 60°, 70°) [3]. На рис. 2 представлены результаты (перемещения увеличены). Наименьшая концентрация напряжений достигается при наибольшем угле раствора — 70° (рис. 2, в).

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

а б в

Рис. 2. Замковые соединения

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

СПИСОК ЛИТЕРАТУРЫ

1. Ц в и к Л. Б. Принцип поочередности в задачах о сопряжении и контакте твердых деформируемых тел // Прикладная механика. - 1980. - Т. 16. - Ч. I. - С. 1318.

2. МожаровскийН. С. Качаловская Н. Е. Приложение методов теории пластичности и ползучести к решению инженерных задач машиностроения: В 2 ч. Ч. 2: Методы и алгоритмы решения краевых задач. - Киев: Вища шк., 1991.-287 с.

3. МавлютовР. Р. Концентрация напряжений в элементах конструкций. - М.: Наука, 1996. - 240 с.

4. Б а р л а м Д. М. Решение контактной задачи теории упругости методом конечных элементов // Проблемы прочности. - 1983. - № 4. - С. 39-43.

Статья поступила в редакцию 03.10.2011