Вычислительные технологии
Том 11, № 6, 2006
ОБ ОДНОМ ИТЕРАЦИОННОМ МЕТОДЕ РАСЧЕТА НАПРЯЖЕНИЙ В НЕОДНОСВЯЗНЫХ ТЕЛАХ
В. В. СтружАнов Институт машиноведения УрО РАН, Екатеринбург, Россия e-mail: [email protected]
An iteration method of stress-strain state calculation in elastic non-simply connected solids with arbitrary number of cavities has been developed for both force loading and kinematic loading. The method is based on the solutions of boundary value problem of linear elasticity theory for simply connected area with given boundary conditions and boundary value problems of residual stress determination in simply connected area with zero boundary conditions and given inelastic residual strains.
Полости в элементах конструкций являются концентраторами напряжений, поэтому определение напряженно-деформированного состояния в их окрестностях имеет большое значение для прогнозирования прочности изделий. Решение же краевых задач теории упругости для неодносвязных тел является трудной математической задачей. Известны аналитические решения только для некоторых конкретных неодносвязных областей [1, 2].
В данной работе предлагается итерационный метод, позволяющий рассчитывать напряжения в неодносвязных телах с любым числом полостей произвольной формы. Показано, что итерации сходятся к точному решению и, следовательно, возможно получение приближенного решения с наперед заданной точностью. Метод основан на использовании решений краевой задачи для односвязной области при заданных условиях на границе и краевых задач об определении собственных напряжений в односвязной же области при нулевых граничных условиях, когда в зонах полостей заданы собственные деформации. Причем данный алгоритм позволяет использовать численные методы решения указанных выше краевых задач.
1. Рассмотрим неодносвязное тело П с полостями V (i = 1, 2,... , m). Внешняя граница тела — кусочно-гладкая поверхность Г, а кусочно-гладкие поверхности Ф^ — границы полостей. В общем случае поверхность Г разбивается на две части — Г1 и Г2 (Г = Г + Г2). На границе Г1 заданы перемещения w(у), y £ Г1, а на границе Г2 — поверхностные силы t(y), у £ Г2. Кроме того, на тело действуют объемные силы /(ж), x £ П. Здесь w,t,f — трехмерные вектор-функции. Свойства материала заданы, вообще говоря, неоднородным анизотропным симметричным тензором четвертого ранга модулей упругости H(ж), ж £ П, или соответствующим тензором модулей податливости W = H-1(W ■ ■H = H ■ ■W = I). Здесь I — единичный тензор четвертого ранга, а двумя точками обозначено двойное скалярное произведение тензоров [3].
Наряду с телом П рассмотрим односвязную область V = П U Пт, где Пт = и™^. Свойства материала, заполняющего эту область, определим тензором модулей упругости
© Институт вычислительных технологий Сибирского отделения Российской академии наук. 2006.
С = Л1Н + Л2Ь. Здесь индикаторные функции А1 = 1, Л2 = 0 при х £ П, Л1 = 0, Л2 = 1 при х £ Пт, а тензор Ь представляет собой непрерывное продолжение тензора Н в области V;. Отметим, что тензор модулей податливости в области V имеет вид Б = А1Ж + Л2Ь-1. Очевидно, что Б ■ -С = (Л1Н + Л2Ь) ■ •(Л1Ж + Л2Ь-1) = Л11 + Л21 = I, так как Л1Л1 = Л1,
Л2Л2 = Л2, Л^2 = 0, Л1 + Л2 = 1.
Напряженно-деформированное состояние в области V при тех же объемных силах и граничных условиях дает решение краевой задачи
V ■ а' = /, е' = def и, а' = С ■ -е', а' ■ п\Г2 = ¿, и\Г1 = /, (1)
где V ■ а' = / — уравнение равновесия, е' = def и — соотношение Коши, а' = С ■ ■е' — обобщенный закон Гука [3], тензоры напряжений и деформаций а' и е' — симметричные тензоры второго ранга; и(х), х £ V, — вектор перемещений; п — вектор внешней нормали к поверхности Г.
Наложим на решение задачи (1) напряженно-деформированное состояние, определяемое краевой задачей
V ■ р'' = 0, е' = defи, р'' = С ■ -(е' - е*), р'' ■ п\г2 = 0, и|г1 = 0, (2)
и потребуем, чтобы выполнялись следующие условия:
Л2(а' + р'') = 0; (3)
(а' + р'') ■ пг[фг = 0, (4)
где пг — векторы внешних нормалей к поверхностям фг относительно области П. Отметим, что уравнения (2) определяют тензор напряжений р'' и тензор совместных деформаций е' при защемленной границе Г1 и свободной от усилий границе Г2, если в области V задан тензор первоначальных несовместных деформаций е*.
Из условия (3) и определяющих соотношений краевой задачи (2) следует, что
Л2(а' + р'') = Л2(а' + д' - а*) = 0
или
Л2а* = Л2(а' + д'), (5)
где д' = С ■ ■е', а* = С ■ ■е*. Условие (4) является следствием условия (3). Действительно, деформации, составляющие тензор е*, не удовлетворяют условиям совместности только в области V, а в областях V имеет место совместность этих деформаций, что следует из равенства (5). Поэтому, если освободить области V от связей, то деформации е* не вызовут появления в них напряжений [4]. Следовательно, напряжения там возникают только в результате воздействия окружающего материала. Тогда для реализации в областях V тензора напряжений р'' = -а' необходимо в силу единственности решений краевых задач теории упругости приложить на границах Фг усилия, равные (-а'■Пг|ф4). Отсюда р''■Пг\ф1 = -а' ■ пг\ф,,, или (а' + р'') ■ П¿\ф, = 0.
Очевидно, что сумма решений краевых задач (1) и (2), а именно а = а' + р'', е = е' + е' при выполнении условий (3) и (4) определяет напряженно-деформированное состояние в неодносвязном теле П для заданных на внешней поверхности Г граничных условиях.
2. Представим теперь равенство (5) в виде операторного уравнения, которому должен удовлетворять тензор а*, чтобы выполнялись условия (3), (4). Имеем
Л2а* = Л2^Л2а* + Л2а'
(6)
где А — оператор, определяющий решение краевой задачи (2), т.е. ставящий в соответствие тензору \2а* тензор д', связанный обобщенным законом Гука с тензором совместных деформаций е'.
Для выяснения конкретного вида оператора А рассмотрим энергетическое гильбертово пространство Т, состоящее из всевозможных симметричных тензоров второго ранга — тензоров напряжений, определенных в области V [5]. Скалярное произведение и норма в нем заданы формулами
(Р1/Р2) = !Р1 • • Б • !И2 = (р,р).
V
Известно [5], что пространство Т является ортогональной суммой следующих подпространств:
Т = {4 : Лик Б • • 4 = 0}, Т2 = {4' : V • 4' = 0, 4' • п|г = 0},
Б1 = {р : Лик Б • •р' = 0, и|Г1 = 0}, В2 = {р'' : V • р'' = 0, р'' • Щ|г2 = 0},
т. е. Т = Т1 + Т2, Т = + 02. Здесь Лик — оператор несовместности [3].
Итак, любой тензор р £ Т можно представить единственным образом — суммой р = 4 + 4' или суммой р = р' + р'', где 4 £ Т1, 4' £ Т2, (4, ('') = 0, р' £ 01, р'' £ Б2, (р',р'') = 0, 4 = Р1р, 4' = Р2р, 4 = ^1р, 4' = Я2р [5]. Здесь Р1, Р2 — операторы ортогонального проектирования соответственно на подпространства Т1, Т2, а Q1, Q2 — ортопроекторы соответственно на Д1, 02, причем Р1 + Р2 = Л, Q1 + Q2 = Л, где Л — тождественный оператор.
Далее очевидно, что Т2 С 02 и С Т1. Тогда [6] Q1P2 = Р^1 = 0 (проекторы Q1 и Р2 ортогональны, так как ортогональны подпространства и Т2) и, кроме того, имеют место равенства
Р^1 = Ql, QlPl = Ql, Q2P2 = Р2, Р2 Q2 = Р2. (7)
Эти равенства проверяются непосредственно. Например,
Р^1 = (Л - Р2^1 = Ql - P2Ql = Ql.
Отметим наконец, что элементами подпространства Т1 являются решения краевой задачи
V • 4 = /, е' = ¿еГи, 4 = С • • е', 4 • п|г = Ь (8)
для всевозможных / и / решения задачи
V • 4' = 0, е,' = аеГи, 4' = С • • (е' - е*)', 4' • Щг = 0 (9)
для всевозможных симметричных тензоров второго ранга е* являются элементами подпространства Т2; решения краевой задачи
V • р' = /, е' = аеГи, р' = С • • е', и|г1 = 0, р' • п|Г2 = Ь
для всевозможных / и Ь входят в подпространство Д1, а решения задачи
V • р'' = 0, е' = аеГи, р'' = С • • (е' - е*), = V, р'' • Щг2 = 0
для всевозможных е* и V входят в подпространство П2. Очевидно, что решение задачи (2) принадлежит подпространству П2.
Найдем теперь общий вид решения системы (2). Задачу будем решать в два этапа. Сначала запишем решение задачи (9). В работе [4] показано, что тензор напряжений здесь определяется выражением д'' = —Р2а* (а* = С • -б*), а совместные деформации и перемещения являются решениями уравнений (8) при / = V • а*, £ = а* • п|г (задача В^), причем д' = Рга*. Итак, получили напряженно-деформированное состояние в теле V со свободной границей при заданном поле первоначальных деформаций, определенном тензором б*. Обозначим через к вектор перемещений, которые получают при этом точки границы (и|г = к).
На втором этапе необходимо решить задачу (1) при / = 0 и граничных условиях и|г1 = —к, а'ь • п|п = 0 (задача В2) и затем наложить это решение на решение задачи (9). В результате получим искомое решение задачи (2).
Рассмотрим дополнительно систему
V • а а = V • а*, б = def и, аа = С • •б, и|Г1 = к, аа • п|Г2 = а* • п|Г2,
где на границе Г1 заданы перемещения, полученные при решении задачи Следовательно, решения этой системы и задачи В1 совпадают, т.е. аа = д' = Р1а*. Сложим решения этой задачи и задачи В2. В результате получим, что а а + а'ь £ Д1, где а а £ Т1, а'ь £ Д2. Отсюда а а = (аа + а[) + (—а[), т. е. имеет место разложение тензора а а на сумму элементов из подпространств и Д2. В силу единственности такого представления находим, что а'ь = — Оа = — ^Аа*.
Итак, искомое решение задачи (2) с учетом равенств (7) можно представить в виде
р'' = —Р2а* — ^Аа* = —(Л — А)а* — ^Аа* = —а* + (Л — ^Аа* =
= —а* + О^а* = —а* + ^1а* = —(Л — ^1)а* = —^2а*.
Наконец, р'' = д' — а* и д' = а* — О2а* = (Л — О2а*) = О1а*. Таким образом, оператор А является ортопроектором О1.
3. Будем искать решение уравнения (6) на множестве Л2Т. Если а* £ Л2Т, то Л2а* = а*(Л2Л2 = Л2). Тогда уравнение (6) принимает вид
а* = Л2^1а* + Л2 а'. (10)
Оценим норму оператора Л2О1. Имеем ||Л2О1|| < ||Л2||||О1||. Известно [6], что норма ортопроектора ||О1|| = 1. Далее
||Л2р||2 = / Л2Р ••£ ^р^ = [ р ••£ • •р^ < [ р ••£ • •р^ = ||р||2.
иу ¿У1 ¿V
Здесь р• •£• •р — положительно-определенная квадратичная форма, р — некоторый элемент из множества М С Т, которое составляют тензоры, определенные в области V', где V' — любая область, входящая в V или совпадающая с V, причем V С V'. Тогда, рассматривая Л2 как оператор, действующий из М в Л2Т, получаем ||Л2|| < 1. Если же Л2 действует из Л2Т в Л2Т, то ||Л21| = 1. Очевидно, что оператор О отображает элементы из Т в М, кроме элементов Л2р' £ Д1. Отсюда следует, что оператор Л2О1, определенный на множестве Л2Т\(Л2ТП ^1), имеет норму ||Л2О1|| < 1. Следовательно, он является оператором сжатия, и решение уравнения (10) представимо сходящимся рядом
а* = ^(Л2О1)гаЛ2а'
П=1
Таким образом, алгоритм решения исходной неодносвязной задачи состоит в следующем:
— для односвязного тела V определяем напряженно-деформированное состояние, отвечающее заданным на внешней границе краевым условиям (краевая задача (1));
— в полостях задаем тензор а* = Л2а', где тензор а' — решение задачи (1), и решаем краевую задачу (2) для е* = Б • • а* = Б • • Л2а', определяя тензор д[;
— задаем в полостях тензор а* = Л2д[ и снова решаем задачу (2) для е* = Б • • а*, определяя тензор д'2;
— и т. д.
4. Проиллюстрируем данную методику на примере равномерного растяжения тонкого кругового кольца с внешним Ь и внутренним а радиусами. Сначала будем задавать радиальное перемещение т точкам внешней границы.
Решение краевой задачи (1) для диска (односвязное тело V) имеет вид
а'г = а" = Ет/Ь(1 — V), бг = б" = т/Ь, (12)
где индексами г, в обозначены соответственно радиальные и тангенциальные напряжения и деформации; Е — модуль Юнга; V — коэффициент Пуассона. В этом примере объемные силы равны нулю, граница Г — внешняя окружность, Г2 = 0. Решение задачи (2), когда в сплошном диске в области, совпадающей с полостью, заданы первоначальные деформации е* = е* = е*, дают выражения
рГ = р" = —[(1 + V )а2/Ь2 + (1 — V )]Ее*[2(1 — V )]-1, дГ = д" = (1 + V )(1 — а2/Ь2)Ее*[2(1 — V)]-1,
еГ = е" = 0.5е*(1 + v)(1 — а2/Ь2), 0 < г < а, (13)
Рг," = —[(1 + V) ± (1 — V)Ь2/г2]Ее*а2[2Ь2(1 — V)]-1,
д'г," = Рг,", е'г" = —(1 ± Ь2/г2)е*а2(1 + V)[2Ь2]-1, а < г < Ь.
Используя это решение, находим значение величины Л2а' и вид оператора Л2^, которые фигурируют в формуле (10). Имеем
Л2а' = Л2Ет[Ь(1 — V)]-1, Л2^1 = 0.5Л2(1 + v)(1 — а2/Ь2),
где Л2 = 1 при г € [0,а), Л2 = 0 при г € [а, Ь]. Теперь, решая уравнение (10) методом последовательных приближений (формула (11)), находим
а* = Л2а* = Л2а* = 2Л2 ЬшЕ [(1 — V )[(1 — V )Ь2 + (1 + V )а2]]-1.
Отсюда е* = (1 — V)Е-1а*. Подставляя это значение в формулы (13), получаем решение задачи (2), сумма которого с решением (12) определяет искомое напряженно-деформированное состояние кольца:
аг = а" = 0, бг = б" = 2тЬ[(1 — v)Ь2 + (1 + v)a2]-1, 0 < г < а,
аг" = ЕЬт[(1 — V)Ь2 + (1 + V)а2]-1(1 т а2/г2),
ег," = шЬ[(1 — ^)62 + (1 + V)а2]-1[(1 — V) ^ (1 + V)а2/г2], а < г < Ь.
Пусть теперь кольцо растягивается посредством приложения к точкам внешней окружности равномерного радиального растягивающего усилия интенсивностью д, т.е. граница Г1 = 0, Г2 — внешняя окружность. В этом случае решение задачи (1) имеет вид
а' = а'д = д, б' = б" = д(1 — v)E-1. (14)
Решение задачи (2) при свободной от напряжений границе с заданными во внутренней области, совпадающей с полостью, первоначальными деформациями б* = б* = б* дают выражения
р' = р"' = 0.5б*Е (а2/Ь2 — 1),
д' = д" = 0.5б*Е (1 — V )-1[(1 + V) + (1 — V )а2/Ь2],
б' = б" = 0.5б*[(1 + V) + (1 — V)а2/Ь2], 0 < г < а, (15)
р"',г = 0.5б*Еа2(1/Ь2 ± 1/г2),
Р," = р',", е",г = 0.5б*а2[(1 — V)1/Ь2 ± (1 + V)1/г2], а < г < Ь.
Используя эти формулы, находим величину Л2а' и вид оператора Л2О1 для подстановки их в уравнение (10). Имеем Л2а' = Л2д, Л2О1 = 0, 5Л2[(1 + V) + (1 — V)а2/Ь2]. Отсюда, решая уравнение (10), находим а* = Л2а* = Л2а" = 2Л2д/[(1 — V)(1 — а2/Ь2)], и тогда б* = 2д/Е (1 — а2/Ь2).
Отметим, что в данной задаче оператор О равняется оператору Р1. Подставляя найденное значение б* в формулы (15) и складывая получившиеся выражения с формулами (14), запишем искомое решение
аг = а" = 0, ег = б" = 2дЬ2[Е(Ь2 — а2)]-1, 0 < г < а,
а">г = дЬ2(1 ± а2/г2)/(Ь2 — а2), б">г = дЬ2[(1 — V) ± (1 + V)а2/Ь2][Е(Ь2 — а2)]-1, а < г < Ь.
Легко видеть, что, как в первом случае кинематического нагружения кольца, так и во втором случае силового нагружения, полученные по изложенной выше методике напряженно-деформированные состояния кольца совпадают с известными решениями задачи Ляме.
5. В заключение отметим, что необходимые для реализации разработанного алгоритма решения краевых задач (1) и (2) можно получать, используя численные методы. И в данном случае методика позволяет находить решения задач для неодносвязных тел. Это утверждение было проверено при определении напряженно-деформированного состояния тонкой однородной изотропной пластины с несколькими прямоугольными отверстиями, растягиваемой равномерно распределенными по торцам усилиями. Пластину высотой 20 мм и шириной 10 мм разбивали на квадраты со стороной 1 мм, которые затем делили диагональю, идущей от правого верхнего угла квадрата в левый нижний угол, на прямоугольные треугольники. Затем из квадратов образовывали несколько отверстий различной конфигурации. Итерационный процесс заканчивали тогда, когда абсолютная
величина разности между величинами деформаций и напряжений, вычисленных на предыдущем и последующем шагах процесса, становилась меньше наперед заданного малого числа (10-9). В результате в областях отверстий возникли напряжения, величины которых имеют порядок 10-4, т.е. практически равные нулю. Затем задачу решали методом конечных элементов с традиционным подходом к неодносвязному телу: явным образом учитывали границы отверстий, разбивая на элементы только области вне отверстий, причем сохраняя вид использованных ранее треугольных элементов. Найденное напряженно-деформированное состояние совпало с состоянием, полученным по изложенной выше методике.
Список литературы
[1] Слвин Г.Н. Распределение напряжений около отверстий. Киев: Наукова думка, 1968.
[2] Мироненко Н.И. Модифицированный метод Д.Н. Шермана // Изв. АН СССР. МТТ. 1987. № 4. С. 143-147.
[3] Лурье А.И. Теория упругости. М.: Наука, 1970.
[4] Соколов А.Г., Стружанов В.В. Об одной задаче оптимизации напряженного состояния в упругом теле // Прикл. мат. и механика. 2001. Т. 65, вып. 2. С. 317-323.
[5] Михлин С.Г. Вариационные методы в математической физике. М.: Наука, 1979.
[6] Канторович Л.В., Акилов Г.П. Функциональный анализ. М.: Наука, 1977.
Поступила в редакцию 10 августа 2005 г.