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

УДК 539.3

М. Е. Яковлев

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ПОЛИКОНТАКТНОГО ВЗАИМОДЕЙСТВИЯ

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

E-mail: me-yakovlev@rambler.ru

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

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

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

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

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

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

j(u(x)) + Qi(x) = 0, x G Ga,i,j = 1, 2, а = 1,n; (1) — граничные условия (кинематические и силовые соответственно)

u(x) = <(x), x G S^ С 6Ga, а = 1, и; (2)

ап(и)и, = pf(x), x G Spa С 6Ga, i, j = ТД а = T~U; (3)

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

ец (x) = 2 (ui,j (x) + Uj,i (x)) x G Ga, i, j = ТТ2, а = T~U; (4)

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

а = D : (е - е0). (5)

Здесь D — матрица Гука, а = {а} = — вектор напряже

Гп

ний, £ = {е} = <е2Л — вектор деформации, £0 = |е0} — век-1^12

тор начальной деформации, Qi — компоненты массовой (объемной)

Г«1о(х)1

силы, и^(х) = < > — вектор заданных перемещений точек по-

I и20(х) I

верхности Р«(х) — компоненты заданной распределенной нагрузки

|р?(х)1

ра (х) = < > на поверхности Б"

[Р2 (Х^

В случае термоупругости вектор начальной деформации имеет вид (аАТ ^

£0 = < аАТ >, где ДТ — изменение температуры относительно неде-1 0

формированного состояния, а — коэффициент температурного расширения.

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

иЩ(х) - <(х) = ¿Щ7(х) (6)

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

аЩ(х) = -аЩ(х) < 0, (7)

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

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

Основные процедуры альтернирующего метода Шварца. Альтернирующий метод Шварца является итерационным методом. Рассмотрим его работу в случае двух тел А и В. Для численного решения контактной задачи используется МКЭ. Пронумеруем узлы контактной поверхности Бк и введем в рассмотрение два вектора — {ик} и {Дк}, первый из которых составлен из компонент перемещений и и V узлов, расположенных на контактной поверхности Бк, а второй — из компонент / и д узловых сил тех же узлов.

На первом шаге на контактных поверхностях тел А и В соответственно Б^ и Б^ задают начальные перемещения и(х) |5л = и^(х)

1Станкевич И. В., Яковлев М. Е., Си Ту Хтет Разработка алгоритма контактного взаимодействия на основе альтернирующего метода Шварца // Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки. Спец. выпуск "Прикладная математика". — 2011. - С. 134-141.

и и(х) | 3б = и^ (х), которые имеют смысл дополнительных кинематических условий. Далее решают независимо две подобные задачи теории упругости для тел А и В. Затем вычисляют поверхностные силы р^(х) и р^ (х) на контактных поверхностях Б^ и Б^ и их корректируют так, чтобы выполнялись силовые контактные условия (7). Для коррекции компонент вектора контактных узловых сил {Дк} используется формула:

2n+1 f J (A),m

2n

-JA

f J (A),m

а2 n

— а(А),т

+

(B),s

2n

(A),m

и = 0,1, 2 ..., (8)

где а

2n

(A),m

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

2n

вектор кон-

л2

•fJ (B),s

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

щий на контактной поверхности БАА тела А;

тактных узловых сил сх поверхности БВ тела В.

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

ра {Uk}

(A)

выполняется в соответствии с соотношением

2n

(A),m

u

=

2n 1

+ а

(A),m

2n-1 (A),m

, и = 0;

(А),т

}2n-1

(B),s и = 1, 2 ...,

2n 1

(A),m

(9)

где а

2n-1

(A),m

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

щий на контактной поверхности S^A тела A;

2n 1

(B),s

вектор пере-

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

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

о

v

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

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

поверхностях тел задаются перемещения и(х) = и^ (х) и реша-

к

ются п отдельных задач теории упругости. На втором шаге часть этих кинематических граничных условий (на подмножестве К1 множества контактных поверхностей К) заменяется скорректированными силовыми аналогично случаю двух тел, а часть (на подмножестве К2 = К \ К!) остается прежними, таким образом, чтобы для любого тела соответствующая задача теории упругости имела единственное решение, и снова решаются п отдельных задач. На третьем шаге кинематические граничные условия на контактных поверхностях из множества К1 заменяются на скорректированные силовые по формуле (10), а силовые граничные условия на поверхностях из множества К2 — на скорректированные кинематические по формуле (9). Далее корректировка кинематических и силовых граничных условий чередуется таким же образом до сходимости на всех поверхностях.

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

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

расположенных нагретых пластин. Температура во всех пяти телах имеет одинаковое распределение и линейно убывает вдоль оси абсцисс. Сверху приложена распределенная нагрузка р = -150 МПа. На рисунке приведена полученная геометрия с учетом деформирования (масштаб увеличен).

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

Работа выполнена в рамках гранта поддержки ведущих научных школ № НШ-255.2012.8.

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