Software & Systems

Поведение элементов сложных систем определяется их физико-механическим откликом на силовые, температурные, диффузионные и другие воздействия. Это приводит к постановке связных, в основном нелинейных задач. Перекрестные эффекты определяются феноменологическими составляющими соответствующих систем уравнений, в которых на основе подхода Онзагера для связных физико-механических задач устанавливается связь [1] между обобщенными потоками Jk и обобщенными силами Xk: Jk=LkmXm.

Для деформационных задач эта связь, как правило, устанавливается линейным постулатом Гука. Однако на его основе невозможно описать фазовый переход второго рода от упругости к диссипативному пластическому течению [2], а непрерывные физические величины, такие как температура, концентрация примеси и т.п., связать с условным пределом текучести, поскольку последняя величина является экспериментально-точеч­ной. Как гипотезу, имеющую подтверждение на практике, примем, что коэффициенты взаимности Lkm, постоянные по Онзагеру, могут нелинейно зависеть от обобщенных сил Lkm=Lkm(Xm). Для деформационных задач примем постулат Коши, в соответствии с которым в упругой области зависимость коэффициента Онзагера от деформации в одномерном случае будет следующей:

L(e)=E× (1–ke+me2).                                               (1)

Тогда поведение рассматриваемой системы (деформирование образца) полностью определяется потенциалом

и описывается вариационным уравнением

.

Исследования этой системы методами синергетики [3] окончательно приводят к уравнению в упругой области:

,

, .                                         (2)

Схематично отклик системы на внешнее воздействие показан на рисунке 1. Здесь к.т. – кри- тическая точка. Если внешнее напряжение выше напряжения в этой точке, происходит хрупкое разрушение. АВ – метастабильное состояние, в котором возможен фазовый переход второго рода с образованием диссипативных структур [2]; ВС – неустойчивый участок, на котором возможен перескок на другую устойчивую ветвь; CD – площадка текучести; АВС – зуб текучести.

Такой подход приводит к решению задач нелинейной упругости. Использование метода конечных элементов требует решения систем большой размерности и увеличенного времени счета из-за трудоемкости итерационного процесса. Возможны несколько путей повышения скорости и точности численного решения. В алгоритме применяются два из них: метод граничных элементов, позволяющий уменьшить размерность матрицы разрешающей системы, и максимальное использование аналитических вычислений [4–6].

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

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

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

Следует отметить, что используемый алгоритм позволяет проводить распараллеливание процесса вычислений.

Построение граничного интегрального уравнения

Рассматривается задача нелинейной упругости в области V с границей Г для однородного изотропного материала. Предполагается, что в упругой области коэффициент Онзагера L(e) определяется согласно постулату Коши (1). В многомерном случае L(eij) для изотропного тела зависит от инварианта тензора деформаций I2, и его можно представить в виде L(I2)=E(1–F(I2)), F(I2)=kI2– , где k и m – положительные постоянные, (2). Связь тензора напряжений sij и тензора деформаций примет вид

.                  (3)

Здесь G – модуль сдвига; u – коэффициент Пуассона; eij – компоненты тензора деформаций; e – объемная деформация; dij – символ Кронекера. В плоском случае .

Коротко можно записать sij=(1–F(I2))Eijkmekm, где  – тензор упругих постоянных.

На части границы Г1 заданы перемещения , а на части Г2 заданы силы . Чертой отмечены заданные функции.

Как известно, метод граничных элементов [7] сводит решение нелинейной упругой задачи для тела V с границей Г к решению интегрального уравнения:

(4)

Здесь x=(x1, x2) – координаты точки интегрирования; x=(x1, x2) – координаты точки наблюдения; ui(x) – i-я компонента перемещения точки тела; pi(x) – i-я компонента усилий;  – компоненты перемещений;  – компоненты усилий фундаментального решения. Индекс x у символов dГx и dVx показывает, что интегрирование проводится по переменным x.

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

Согласно методу граничных элементов [7], разобьем границу области Г на N граничных элементов Гn, в качестве которых в данном случае взяты отрезки прямых. Элементы и узлы нумеруются последовательно, против часовой стрелки. При этом область получается слева, а внешняя часть – справа. Каждому элементу ставятся в соответствие два узла: начало элемента, точка x1n, и его конец, точка x2n, в направлении обхода нумерации. Индекс показывает, что эти точки относятся к элементу с номером n.

На элементе перемещения аппроксимируются линейными функциями, усилия – кусочно-посто­янными:  .

Здесь  – значение i-й компоненты усилия на элементе с номером n;  – значение i-й компоненты перемещения в узле k элемента с номером n;  – базисные функции n элемента, , ; Ln – длина элемента. Функции подобраны так, что в одном из узлов она равна 1, в другом 0. Такая система функций называется нормализованной. Вне элемента n на границе функции  обращаются в ноль. Заметим, что функции  определены только на границе и здесь s не вектор, а скалярный параметр, меняющийся вдоль границы. Для вычисления интеграла от нелинейных слагаемых

разобьем область V на M треугольных конечных элементов Vq, q=1, 2, …, M. Возможны четыре варианта взаимного положения точки x и элемента Vq области V.

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

где S – площадь элемента; xk – координаты k узла элемента;  – значения компонент тензора напряжений в узле xk.

2.     Точка x совпадает с вершиной элемента Vn. Интеграл становится сингулярным со слабой особенностью. В этом случае интеграл вычисляется аналитически.

3.     Точка x лежит на границе.

4.     Точка x внутри элемента.

Случаи 3 и 4 сводятся к случаю 2 разрезанием элемента на вспомогательные.

Запишем (4) для узлов xn, k=1, 2, …, N, с учетом непрерывности перемещений , j=1, 2, где  – перемещение k узла элемента n в направлении j:

          (5)

Здесь ,  

                       (6)

В данном случае и далее первый индекс в обозначении перемещений будем опускать. Интегралы (5) по граничным элементам вычислены аналитически. Система (4) решается методом последовательных приближений. Первое приближение [7] ищем, отбросив нелинейные слагаемые:

…

Здесь ,  – значения компонент тензора деформаций и, соответственно, функции F(I2(e)), вычисленные по r-приближению решения. Верхние индексы показывают номер узла (n) и номер приближения (r). Решение линейной задачи описано в [6]. Блок-схема алгоритма приведена на рисунке 2.

Пример применения алгоритма

Для иллюстрации работы алгоритма рассматривается задача о растяжении неоднородной изотропной полосы с круговым включением (рис. 3). Длина полосы L=100 м, ширина h=10 м, радиус отверстия r=1 м.

Область включения обозначим V(2), соответственно обозначим границы областей: Г(1) – внешняя граница области V(1), Г(2) – внутренняя граница области и граница включения,  – упругие коэффициенты внешней области,  – упругие коэффициенты включения. 

Верхний индекс в скобках указывает, к какой области относится параметр или переменная.

На общей границе области V(1) и включения выполняются условия непрерывности перемещения:  i=1, 2, и условия рав­новесия усилий:  i=1, 2.

Использование изложенного выше подхода приводит к системе нелинейных уравнений:

      (7)

которая решается методом последовательных приближений.

На рисунке 4 приведены графики нормального sn и касательного st напряжений в сечении x1=0 на отрезке 0£x2£5 относительно растягивающего усилия. Модули сдвига заданы относительно предела текучести sp=70 MПa. Вычисления проводились при значениях относительного модуля сдвига G(1)=500, G(2)=250, коэффициент Пуассона u=0,3.

Для сравнения были проведены расчеты при значениях относительного модуля сдвига G(1)=500, G(2)=0, что соответствует отверстию, и при G(1)=500, G(2)=500, что соответствует однородной полосе. Полученные решения совпадают с теоретическими [8].

На рисунке 4 приведен график касательных напряжений в сечении а–а.

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

Литература

1.     Пригожин И. Введение в термодинамику необратимых процессов. М.: Иностранная литература, 1960. 127 с.

2.     Дьярмати И. Неравновесная термодинамика. М.: Мир, 1974. 301 с.

3.     Томпсон Дж.М. Неустойчивости и катастрофы в науке и технике. М.: Мир, 1985. 256 с.

4.     Федотов В.П., Нефедова О.А. Применение модифицированного метода граничных элементов для решения параболических задач // Вестн. Сам. гос. техн. ун-та: Сер. Физ.-мат. науки. 2011. Вып. 4 (25). С. 93–101.

5.     Федотов В.П., Спевак Л.Ф. Аналитическое интегрирование функций влияния для решения задач упругости и теории потенциала методом граничных элементов // Математическое моделирование. 2007. Т. 19. № 2. С. 87–104.

6.     Федотов В.П., Горшков А.В. Численно-аналитический метод решения задач упругости с особенностями // Вестн. Сам. гос. техн. ун-та: Сер. Физ.-мат. науки. 2005. Т. 38. C. 29–34.

7.     Бреббия К., Теллес Ж., Вроубел Л. Методы граничных элементов. М.: Мир, 1987. 524 с.

8.     Седов Л.И. Механика сплошных сред. В 2 т. М.: Наука, 1983.

References

1.     Prigogin I. Introdaction to termodynamics of irreversible processes. Springer Publ., Illinois, USA, 1955, 160 p. (Russ. ed.: Moscow, Inostrannaya literatura Publ., 1960, 127 p.).

2.     Dyarmati I. Neravnovesnaya termodinamika [Non-equilib­rium thermodynamics]. Moscow, Mir Publ., 1974, 301 p.

3.     Thompson J.M.T. Instabilities and catastrophes in science and ingineering. New York, Wiley Publ., 1982, 226 p. (Russ. ed.: Moscow, Mir Publ., 1985, 256 p.).

4.     Fedotov V.P., Nefedova О.А. Application of modified method of boundaty elements to solve parabolic problems. Vestn. Sam. gos. tekhn. un-ta: ser. Fiz.-mat. nauki [The bulletin of Samara State Tech. Univ. Physics and Mathematics Science Series]. 2011, vol. 25, pp. 93–101.

5.     Fedotov V.P., Spevak L.F. Analytical integration of influence functions to solve the elastic problems and potential theory using the boundary elements method. Matematicheskoe modeliro­vanie [Mathematical modeling]. 2007, vol. 19, no. 2, pp. 87–104.

6.     Fedotov V.P., Gorshkov A.V. Numerical-analytical me- thod of solving elastic problem with singularities. Vestn. Sam. gos. tekhn. un-ta: ser. Fiz.-mat. nauki [The bulletin of Samara State Tech. Univ. Physics and Mathematics Science Series]. 2005, vol. 38, pp. 29–34.

7.     Brebbiya C.A., Telles J.C.F.,  Wrobel L.C.  Boundary element techniques. Berlin, Springer – Verlaq Publ., Heidelberg New York, Tokyo 1984. (Russ. ed.: Moscow, Mir Publ., 1987, 524 p.).

8.     Sedov L.I. Mekhanika sploshnykh sred [Mechanics of continua medium]. Moscow, Nauka Publ., 1983, vol. 1–2.