Журнал вычислительной математики и математической физики, 2020, T. 60, № 12, стр. 2028-2049
Прохождение предельных особых точек методом продолжения решения по параметру в задачах неупругого деформирования
Е. Б. Кузнецов 1, *, С. С. Леонов 1, 2, **
1 Московский авиационный институт
(национальный исследовательский университет)
125993 Москва, Волоколамское ш., 4, Россия
2 Российский университет дружбы народов
117198 Москва, ул. Миклухо-Маклая, 6, Россия
* E-mail: kuznetsov@mai.ru
** E-mail: powerandglory@yandex.ru
Поступила в редакцию 06.02.2020
После доработки 18.06.2020
Принята к публикации 04.08.2020
Аннотация
В статье разрабатываются численные методы решения задачи Коши для систем обыкновенных дифференциальных уравнений с одной предельной особой точкой, расположенной на правой границе рассматриваемого интервала изменения аргумента. Сферой применения таких начальных задач являются неупругое деформирование металлических конструкций при различных температурно-силовых режимах в условиях ползучести, расчет прочностных характеристик и оценка остаточных деформаций при проектировании ядерных реакторов, в строительной и аэрокосмической отраслях, машиностроении. При этом численное решение подобных начальных задач сопряжено со значительными трудностями вследствие их плохой обусловленности. Традиционные явные методы могут быть использованы для задач данного класса, но только вне окрестности предельной особой точки, в которой происходит резкий рост погрешности численного решения. Это приводит к необходимости чрезмерного уменьшения шага интегрирования, что увеличивает время счета и делает явные методы вычислительно затратными. Предложено для численного решения использовать метод продолжения решения, заключающийся в замене исходного аргумента задачи на новый, при котором преобразованная задача имеет лучшую обусловленность. Однако наилучший аргумент не дает в данном случае желаемых вычислительных преимуществ, так как приводит к значительному усложнению вида исходной задачи. В качестве более адекватного подхода авторы предлагают использовать специализированный аргумент продолжения решения, называемый модифицированным наилучшим. Для задачи деформирования вплоть до разрушения трубчатых образцов из стали Х18Н10Т показаны преимущества метода продолжения решения по модифицированному наилучшему аргументу по сравнению с традиционными явными методами решения задачи Коши и наилучшей параметризацией. Достоверность результатов подтверждается сопоставлением с экспериментальными данными и результатами других авторов. Библ. 42. Фиг. 2. Табл. 14.
1. ВВЕДЕНИЕ
Во многих задачах прикладной математики и механики приходится сталкиваться с процессами, которые быстро изменяются во времени. Для их описания используются начальные и краевые задачи для систем дифференциальных, функционально-дифференциальных, дифференциально-алгебраических, интегральных уравнений. Для данного класса задач характерно наличие на кривых и поверхностях решений пограничных слоев, контрастных структур, предельных особых точек. В большинстве случаев определяющие уравнения являются нелинейными и получить точное аналитическое решение не удается. К тому же, за счет их плохой обусловленности и жесткости, для подобных задач традиционные явные численные методы становятся малоэффективными. Необходимость в использовании неявных схем и специализированных методов решения жестких задач приводит к увеличению вычислительной ресурсоемкости, что является следствием известных недостатков неявных и специализированных методов, связанных с решением нелинейных систем алгебраических и трансцендентных уравнений, процедурами повышения устойчивости метода и его адаптивности.
С вычислительной точки зрения все указанные недостатки позволяет смягчить или полностью устранить метод продолжения решения по наилучшему аргументу, известный также как метод наилучшей параметризации и arc length method. Он состоит в замене исходного аргумента на новый, отсчитываемый по касательной к интегральной кривой рассматриваемой задачи. При этом преобразованная задача становится обусловленной наилучшим образом, что устраняет все трудности, связанные с резким возрастанием скорости изменения интегральных кривых.
Идея метода продолжения решения по параметру появилась еще в работах А. Пуанкаре и У. Леверье. По сути, метод замены переменной под знаком интеграла и метод малого параметра являются отражением этой общей идеи.
В вычислительных целях впервые метод продолжения решения по параметру был использован в работе бельгийского математика М. Лаэя для построения решений трансцендентных уравнений и их систем [1], [2].
Позднее в работах советского математика Д.Ф. Давиденко был предложен метод сведения системы нелинейных уравнений к системе обыкновенных дифференциальных уравнений [3], [4]. Однако при построении замкнутых кривых этот подход требует смены параметра продолжения решения в окрестности предельных особых точек. Для устранения этого недостатка было предложено использовать параметры продолжения общего вида. Для решения нелинейных краевых задач теории упругости в работе И.И. Воровича и В.Ф. Зипаловой [5] было предложено использовать в качестве параметра продолжения длину кривой множества решений. В этой же работе впервые сделано предположение, что данный параметр продолжения будет в некотором смысле наилучшим для задач, кривые множества решений которых имеют предельные особые точки.
Доказательство этого факта было намечено в статье голландского инженера Э. Рикса [6] применительно к исследованию устойчивости прощелкивающихся и выпучивающихся конструкций. Начиная с конца 70-х годов прошлого века началось систематическое исследование вопросов применения метода наилучшего параметра к задачам математики и механики с предельными особыми точками и точками бифуркации. Отметим две важные монографии Э.И. Григолюка, В.И. Шалашилина [7] и E.L. Allgower, K. Georg [8], подводящие итог полученных результатов к концу 80-х годов прошлого века. В монографии В.И. Шалашилина и Е.Б. Кузнецова [9] впервые было доказано, что параметр, отсчитываемый по касательной к кривой множества решений системы нелинейных уравнений, доставляет задаче наилучшую обусловленность, т.е. является наилучшим. Более того, этот результат был обобщен на системы обыкновенных дифференциальных, дифференциально-алгебраических и функционально-дифференциальных уравнений. В дальнейших работах Е.Б. Кузнецова было рассмотрено обобщение полученных результатов на многомерный случай [10], а также совместно с С.Д. Красниковым предложен алгоритм прохождения точек бифуркации различной коразмерности с использованием редукции Ляпунова–Шмидта [11], [12]. Один из последних результатов связан с развитием нового подхода для систем обыкновенных дифференциальных уравнений – метода продолжения решения по модифицированному наилучшему аргументу, отсчитываемому в направлении, близком к касательному [13], [14].
Помимо указанных выше работ, метод продолжения решения по наилучшему параметру применялся к решению жестких и сверхжестких задач в статьях Н.Н. Калиткина и соавт. [15]–[18]; для исследования устойчивости панелей оболочек из изотропных [19] и ортотропных [20] материалов при геометрически нелинейном деформировании, а также для решения гиперболических систем с предельными особыми точками [21] и физически нелинейных задач [22]. Специальные виды параметров продолжения используют в своих работах Э.И. Григолюк, Е.А. Лопаницын при расчете тонких пологих оболочек с учетом конечных прогибов [23] и С.С. Гаврюшин и соавт. при расчете напряжений и деформаций сложных стержневых и оболочечных элементов конструкций [24].
В качестве иллюстрации применения метода продолжения решения используются задачи определения деформационно-прочностных характеристик элементов конструкций, работающих в условиях сложного напряженного состояния при различных температурно-силовых режимах. В большинстве случаев при исследовании поведения металлов и других материалов, таких как бетон, в холодном состоянии достаточно учитывать только упругую и пластическую деформации. Но при расчете конструкций, работающих при высоких температурах, необходим учет также ползучести материала.
Термином “ползучесть”, согласно известной монографии Ю.Н. Работнова [25, с. 9], “будем называть всю совокупность явлений, которые можно объяснить, допустив, что зависимость между напряжениями и деформациями содержит время, явно или через посредство некоторых операторов”. Более узкое определение в работе [26, с. 241] дает Н.Н. Малинин, оно заключается в следующем: “напряжения и деформации, возникшие при нагружении деталей, изменяются во времени, даже если нагрузки остаются постоянными. Это явление называют ползучестью материала. Одну сторону этого явления – изменение во времени деформаций – называют собственно ползучестью или последействием, а другую – изменение во времени напряжений – релаксацией”.
На фиг. 1 приведены кривые, изображающие зависимости деформации от времени при различных напряжениях. На кривой $\varepsilon (t)$, соответствующей напряжению ${{\sigma }_{2}}$, можно выделить три четко выраженных участка (стадии ползучести) [27, с. 20]:
I – неустановившаяся ползучесть, т.е. участок, на котором скорость ползучести монотонно уменьшается до своего наименьшего значения;
II – установившаяся ползучесть, т.е. участок, на котором скорость ползучести сохраняет постоянное наименьшее значение;
III – участок ускоряющейся ползучести, предшествующий разрушению.
В предположении, что время нагружения до заданного значения напряжения мало по сравнению с длительностью испытания, кривые $\varepsilon (t)$ начинаются со значения деформации ${{\varepsilon }_{0}}$, соответствующего мгновенному нагружению. Мгновенная деформация складывается из упругой ${{\varepsilon }^{e}}$ и пластической ${{\varepsilon }^{p}}$ составляющих. Разность между полной и мгновенной деформацией есть деформация ползучести ${{\varepsilon }^{c}}$ (в дальнейшем верхний индекс $c$ при обозначении деформации ползучести будет опускаться) [27, с. 19–20].
В зависимости от величины напряжения, на кривой $\varepsilon (t)$ могут отсутствовать некоторые стадии ползучести, как это показано на фиг. 1.
Подробный обзор литературы, посвященной ползучести металлов, приводится в работах С.А. Шестерикова и А.М. Локощенко [28]–[33]. Указанные обзоры охватывают временной период с середины 40-х годов XX века до первой декады XXI века и содержат материалы по экспериментальным работам, критериальному и кинетическому подходу к описанию ползучести, поведению металлов в агресcивных средах. Обзоры содержат широкий материал по зарубежным школам теории ползучести и публикациям наиболее значимых их представителей.
Процесс деформирования в условиях ползучести описывается начальными задачами для систем обыкновенных дифференциальных уравнений. На фиг. 1 видно, что при приближении к моменту разрушения скорость деформации ползучести возрастает, стремясь к бесконечности. Тоже происходит и при наличии упрочнения конструкции в окрестности начального момента времени на стадии неустановившейся ползучести. Это означает, что при описании всех трех стадий ползучести определяющие соотношения ползучести имеют две предельные особые точки. Традиционные явные методы решения начальной задачи могут перескакивать через момент разрушения за счет накапливаемой ошибки, но не могут преодолеть предельную особую точку в начальный момент времени, что делает их малоэффективными для данного класса задач. Использование же специальных методов решения жестких и плохо обусловленных задач, в том числе и неявных методов, гораздо более трудоемко.
Данная статья является продолжением работ [14], [34], в которых рассматривается применение метода продолжения решения по наилучшему аргументу [9, с. 50–73] и модифицированному наилучшему аргументу [13] к задачам ползучести при наличии всех трех стадий ползучести. Показано, что для подобных задач применение наилучшего аргумента позволяет получить значительные вычислительные преимущества, что выражается в уменьшении времени счета до 50 процентов при сохранении точности вычислений [34]. Еще до 50–60% можно уменьшить время счета, используя модифицированный наилучший аргумент, что достигается за счет упрощения вида определяющих соотношений по сравнению с наилучшим аргументом [14].
При отсутствии стадии упрочнения преобразование к наилучшему аргументу может не давать вычислительных преимуществ по сравнению с исходной задачей. Это происходит по причине усложнения вида преобразованной задачи. На примере задачи разрушения трубчатых образцов из стали Х18Н10Т при неоднородном напряженно-деформированном состоянии исследуем вычислительные особенности применения модифицированного наилучшего аргумента.
2. ОПРЕДЕЛЯЮЩИЕ СООТНОШЕНИЯ ПОЛЗУЧЕСТИ
Впервые ползучесть наблюдалась еще в начале XIX века в экспериментах К.Л.М.А. Навье (1826), Г.Г. Кориолиса (1830) и Л.Ж. Вика (1834). Но систематические исследования ползучести металлов начались только в 40-х годах прошлого века. В конце 50-х годов XX века Л.М. Качанов [35] и Ю.Н. Работнов [36] пришли к выводу, что использующихся в то время терминов механики деформируемого твердого тела (тензоры напряжений и деформаций, вектор перемещений) недостаточно для описания процесса длительного разрушения материалов и элементов конструкций в условиях ползучести. Ими был предложен новый подход к исследованию длительной прочности. Этот подход основан на использовании введенного параметра поврежденности $\omega $. Исходному состоянию материала (при $t = 0$) соответствует значение $\omega = 0$, при разрушении параметр $\omega $ становится равным единице. Этот подход в дальнейшем был развит в монографии Ю.Н. Работнова [25] и получил название кинетической теории ползучести.
При использовании уравнений кинетической теории к расчету ползучести и длительной прочности конструкции задачи определения напряженно-деформированного состояния и длительной прочности совмещаются. Одним из наиболее распространенных подходов к моделированию процессов ползучести и разрушения металлических конструкций является использование уравнений теории структурных параметров, см. Ю.Н. Работнова [25, с. 363–364], которые в одномерном случае с одним скалярным параметром поврежденности можно записать в виде системы, состоящей из определяющего уравнения вида
связывающего деформацию ползучести $\varepsilon $ с действующим напряжением $\sigma $, и эволюционного уравнения вида описывающего развитие поврежденности в конструкции. Здесь $t$ – время, $T$ – температура, $\sigma $ – в общем случае переменное действующее напряжение, функциональные зависимости, входящие в правые части уравнений (1) и (2), определяются по результатам эксперимента.В случае постоянной температуры $T = {\text{const}}$ запишем систему (1), (2) в виде
(3)
$\frac{{d\varepsilon }}{{dt}} = \frac{{{{f}_{1}}(\sigma )}}{{\Psi (\omega )}},\quad \frac{{d\omega }}{{dt}} = \frac{{{{f}_{2}}(\sigma )}}{{\Psi (\omega )}}.$В качестве начальных условий для системы ОДУ (3) берутся однородные
Искомыми функциями в задаче (3), (4) являются $\varepsilon (t)$ и $\omega (t)$, $\sigma (t,\varepsilon )$ – известная непрерывная функция времени и деформации ползучести, в частности постоянная величина. Решение задачи ищется в области $V = \left\{ {\left( {\varepsilon ,\omega ,t} \right)\,|\,0 \leqslant \varepsilon \leqslant \varepsilon *,\;0 \leqslant \omega \leqslant 1,\;0 \leqslant t \leqslant t{\text{*}}} \right\}$, где $\varepsilon {\text{*}}$ – значение деформации ползучести в момент разрушения, $t{\text{*}}$ – длительная прочность конструкции.
Отметим некоторые особенности системы (3). Рассматриваются процессы деформирования, для которых $\dot {\varepsilon } \geqslant 0$ и $\dot {\omega } \geqslant 0$, т.е. процессы деформирования и накопления повреждений в материале предполагаются монотонными. Функции ${{f}_{1}}(\sigma )$ и ${{f}_{2}}(\sigma )$ будем считать непрерывными, монотонными и положительными. Дополнительно будем полагать, что ${{f}_{1}}(\sigma ) = \varphi (t) \cdot {{\psi }_{1}}(\varepsilon )$ и ${{f}_{2}}(\sigma ) = \varphi (t) \cdot {{\psi }_{2}}(\varepsilon )$, где $\varphi (t)$, ${{\psi }_{1}}(\varepsilon )$ и ${{\psi }_{2}}(\varepsilon )$ являются положительными и непрерывными функциями своих аргументов. В качестве $\Psi (\omega )$ будем рассматривать неотрицательные функции, причем для неупрочняющихся материалов, т.е. таких материалов, для которых стадия неустановившейся ползучести отсутствует, будем применять монотонно убывающие функции, такие что $\Psi (1) = 0$, а для упрочняющихся материалов немонотонные унимодальные функции, для которых $\Psi (0) = 0$ и $\Psi (1) = 0$.
При сделанных предположениях в статье [37] получено решение задачи (3), (4) в интегральном виде
(5)
$\omega = H(\varepsilon ),\quad \varepsilon = {{G}^{{ - 1}}}[\Phi (t)],\quad t* = {{\Phi }^{{ - 1}}}\{ G[{{H}^{{ - 1}}}(1)]\} ;$Решение начальной задачи (3), (4) дается соотношениями (5). Но для получения аналитического решения в явном виде нужно вычислить интегралы, входящие в них. Это удается сделать не всегда, в особенности при громоздкой правой части системы (3). Проиллюстрируем процесс решения подобных задач на тестовом примере. Отметим трудности, возникающие при этом.
3. ЗАДАЧА РАСТЯЖЕНИЯ ОБРАЗЦОВ ИЗ СТАЛИ Х18Н10Т
В качестве тестовой задачи рассмотрим деформирование трубчатых образцов из нержавеющей стали Х18Н10Т под действием постоянной одноосной растягивающей нагрузки при постоянной температуре в условиях ползучести вплоть до разрушения.
3.1. Постановка задачи
Для описания данной задачи будем использовать уравнения теории структурных параметров Ю.Н. Работнова вида (3), которые, в случае отсутствия упрочнения, запишем в виде (см. [38])
(6)
$\frac{{d\varepsilon }}{{dt}} = \frac{{A \cdot {{\sigma }^{n}}}}{{{{{(1 - {{\omega }^{r}})}}^{n}}}},\quad \frac{{d\omega }}{{dt}} = \frac{{B \cdot {{\sigma }^{n}}}}{{{{{(1 - {{\omega }^{r}})}}^{n}}}}$При деформировании конструкции в условиях ползучести под действием постоянных нагрузок введем гипотезу о равномерном распределении деформации по длине образца. Тогда, в случае малых деформаций, выражение для напряжения примет вид
где ${{\sigma }_{0}}$ – постоянное начальное напряжение.Учитывая соотношение (8), перейдем от (6) к системе (см. [38])
3.2. Аналитическое решение
Поделив первое уравнение системы (6) или (9) на второе, можно получить связывающее деформацию $\varepsilon $ и параметр $\omega $ соотношение
Но даже используя соотношение (10), аналитическое решение задачи (6), (4) или (9), (4) удается получить не для всех значений параметров $n$ и $r$. В статье [37] доказано, что уже при постоянном напряжении для задачи (6), (4) справедлива
Теорема 1. Задача (6), (4) при постоянном напряжении $\sigma = {{\sigma }_{0}} = {\text{const}}$ интегрируется аналитически тогда и только тогда, когда параметры модели $n$ и $r$ удовлетворяют одному из следующих условий:
1) $n$ - натуральное число;
2) $1{\text{/}}r$ – натуральное число;
3) $1{\text{/}}r + n$ – натуральное число.
Для задачи (9), (4) приходится накладывать более строгие ограничения на значения параметров $n$ и $r$. В статье [37] доказана
Теорема 2. Задача (9), (4) интегрируется аналитически тогда и только тогда, когда параметры модели $n$ и $r$ удовлетворяют одному из следующих условий:
1) $n$ – натуральное число;
2) $r$ – натуральное число;
3) $1{\text{/}}r$ – натуральное число.
При этих условиях существует аналитическое решение задачи (9), (4), даваемое выражением (10) и равномерно сходящимися функциональными рядами
– для натуральных $r$ имеем
– для дробных $r$ имеем
3.3. Экспериментальные данные и параметры модели
Приведем результаты, полученные в работе [38]. В ней описаны испытания на ползучесть при чистом растяжении 21 образца одной плавки нержавеющей стали Х18Н10Т. Испытывались трубки с внешним диаметром 12 мм, толщиной стенки 0.5 мм и рабочей длиной 70–100 мм при постоянных растягивающей нагрузке (с начальными напряжениями ${{\sigma }_{0}}$ = 39.2, 49, 58.8 и 78.5 МПа) и температуре 850°С. Здесь и далее все величины переводятся из системы единиц МКГСС в СИ.
По экспериментальным данным в работе [38] определены характеристики ползучести для задачи (9), (4)
(11)
$A = 6.716 \cdot {{10}^{{ - 9}}}\;{{{\text{ч}}}^{{ - 1}}}\;{\text{МП}}{{{\text{а}}}^{{ - 3.2}}},\quad B = 6.381 \cdot {{10}^{{ - 8}}}\;{{{\text{ч}}}^{{ - 1}}}\;{\text{МП}}{{{\text{а}}}^{{ - 3.2}}},\quad n = 3.2,\quad r = 2.1.$В табл. 1 приведены теоретико-экспериментальные данные о процессе деформирования трубок из стали Х18Н10Т [38], где $t_{e}^{ * }$ и $\varepsilon _{e}^{ * }$ – средние для каждого напряжения ${{\sigma }_{0}}$ экспериментальные значения длительной прочности и соответствующей ему деформации ползучести, $t_{a}^{ * }$ – значение длительной прочности, рассчитанное по формуле
3.4. Явные методы с постоянным шагом
Значения материальных констант (11) для задачи (9), (4) показывают, что согласно теореме 2 аналитическое решение данной задачи получено быть не может. Но даже при существовании аналитического решения не всегда удается его найти. Эта ситуация характерна для большинства прикладных задач ползучести. Поэтому основным инструментом исследования данного класса задач являются численные методы.
Стоит отметить, что задача (9), (4) имеет особенность. В момент разрушения, т.е. при значении $\omega = 1$, правые части системы (9) теряют смысл. Таким образом, момент разрушения является предельной особой точкой, что делает рассматриваемую задачу плохо обусловленной. Это может привести к затруднениям в процессе численного решения.
Рассмотрим применение для решения задачи (9), (4) традиционных явных методов численного интегрирования задачи Коши [39, с. 277–285]:
1) метод Эйлера;
2) метод Рунге–Кутты четвертого порядка точности.
Вначале используем для численного решения постоянный шаг интегрирования.
На фиг. 2 изображены кривые ползучести для задачи (9), (4) для значений начального напряжения ${{\sigma }_{0}} = $ 39.2, 49, 58.8 и 78.5 МПа (отмеченных около соответствующих кривых), полученные методом Рунге–Кутты четвертого порядка точности с переменным шагом интегрирования при заданной точности $\theta = {{10}^{{ - 4}}}$ (результаты решения с переменным шагом интегрирования приведены в п. 3.7.1). Кривые ползучести, полученные при постоянном шаге, имеют аналогичный вид.
В табл. 2 и 3 приведены результаты численного решения задачи для четырех значений начального напряжения ${{\sigma }_{0}}$ = 39.2, 49, 58.8 и 78.5 МПа с шагом интегрирования $h = {{10}^{{ - 1}}},\;{{10}^{{ - 2}}},\;{{10}^{{ - 3}}},\;{{10}^{{ - 4}}}$ ч, где $\varepsilon _{n}^{ * }$ и $\omega _{n}^{ * }$ – расчетные значения деформации ползучести и параметра поврежденности в момент разрушения соответственно, $t_{n}^{ * }$ – расчетное значение длительной прочности конструкции, $N$ – количество шагов по независимой переменной, ${{t}_{c}}$ – время счета. Здесь ${{h}_{t}}$ – шаг интегрирования по аргументу $t$.
Таблица 2.
Решение задачи (9), (4), метод Эйлера с постоянным шагом
| ${{\sigma }_{0}},$ МПа | ${{h}_{t}}$, ч | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$, ч | ${{t}_{c}},$ c | $N$ |
|---|---|---|---|---|---|---|
| 39.2 | 10–1 | 0.0944 | 0.8972 | 52.8 | 0.017 | 528 |
| 10–2 | 0.1 | 0.95 | 52.61 | 0.063 | 5261 | |
| 10–3 | 0.1013 | 0.962 | 52.581 | 0.457 | 52 581 | |
| 10–4 | 0.1049 | 0.9961 | 52.5778 | 4.79 | 525 778 | |
| 49 | 10–1 | 0.1004 | 0.9539 | 26 | 0.013 | 260 |
| 10–2 | 0.1042 | 0.9895 | 25.78 | 0.038 | 2578 | |
| 10–3 | 0.1047 | 0.9948 | 25.749 | 0.246 | 25 749 | |
| 10–4 | 0.1045 | 0.9932 | 25.7451 | 2.196 | 257 451 | |
| 58.8 | 10–1 | 0.0999 | 0.9494 | 14.6 | 0.014 | 146 |
| 10–2 | 0.0962 | 0.9141 | 14.39 | 0.035 | 1439 | |
| 10–3 | 0.1029 | 0.9778 | 14.369 | 0.163 | 14 369 | |
| 10–4 | 0.1042 | 0.9901 | 14.3655 | 1.304 | 143 655 | |
| 78.5 | 10–1 | 0.0928 | 0.8818 | 5.9 | 0.012 | 59 |
| 10–2 | 0.1014 | 0.9635 | 5.75 | 0.016 | 575 | |
| 10–3 | 0.102 | 0.9692 | 5.725 | 0.061 | 5725 | |
| 10–4 | 0.102 | 0.9688 | 5.7218 | 0.539 | 57 218 |
Таблица 3.
Решение задачи (9), (4), метод Рунге–Кутты с постоянным шагом
| ${{\sigma }_{0}},$ МПа | ${{h}_{t}}$, ч | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$, ч | ${{t}_{c}},$ c | $N$ |
|---|---|---|---|---|---|---|
| 39.2 | 10–1 | 0.0889 | 0.8445 | 52.5 | 0.038 | 525 |
| 10–2 | 0.0961 | 0.9133 | 52.57 | 0.207 | 5257 | |
| 10–3 | 0.1012 | 0.9617 | 52.577 | 1.71 | 52 577 | |
| 10–4 | 0.1026 | 0.9743 | 52.5772 | 17.969 | 525 772 | |
| 49 | 10–1 | 0.0882 | 0.838 | 25.7 | 0.027 | 257 |
| 10–2 | 0.0955 | 0.9077 | 25.74 | 0.112 | 2574 | |
| 10–3 | 0.0994 | 0.9442 | 25.744 | 0.877 | 25 744 | |
| 10–4 | 0.1016 | 0.9654 | 25.7445 | 8.194 | 257 445 | |
| 58.8 | 10–1 | 0.0835 | 0.7936 | 14.3 | 0.027 | 143 |
| 10–2 | 0.0938 | 0.8911 | 14.36 | 0.079 | 1436 | |
| 10–3 | 0.1019 | 0.9681 | 14.365 | 0.508 | 14 365 | |
| 10–4 | 0.1029 | 0.9776 | 14.365 | 2.168 | 143 650 | |
| 78.5 | 10–1 | 0.0845 | 0.8031 | 5.7 | 0.026 | 57 |
| 10–2 | 0.0947 | 0.8996 | 5.72 | 0.036 | 572 | |
| 10–3 | 0.0976 | 0.9268 | 5.721 | 0.221 | 5721 | |
| 10–4 | 0.1039 | 0.9874 | 5.7214 | 1.97 | 57 214 |
Из полученных результатов можно сделать следующие выводы.
1. Традиционные методы численного решения задачи Коши позволяют получить решение задачи (9), (4) при использовании постоянного шага интегрирования.
2. Подойти близко к моменту разрушения не удается ни методом Эйлера, ни методом Рунге–Кутты четвертого порядка точности. При этом чем больше начальное напряжение, тем труднее подойти к моменту разрушения. Для увеличения точности решения приходится значительно уменьшать шаг интегрирования.
3. Количество шагов по аргументу $t$ для обоих рассматриваемых методов приблизительно одинаково. Однако метод Рунге–Кутты более трудоемкий в использовании, что влечет увеличение времени счета по сравнению с методом Эйлера.
4. Визуально метод Рунге–Кутты 4-го порядка точности останавливается дальше от момента разрушения, по сравнению с методом Эйлера. Однако по полученным ранее результатам [13], [14] метод Рунге–Кутты дает меньшую среднюю погрешность численного решения.
Замечание 1. Отсутствие возможности нахождения аналитического решения не позволяет получить абсолютную или относительную погрешность полученных приближенных решений. Но то что параметр поврежденности $\omega $ равен единице в момент разрушения, позволяет определять, насколько близко мы подходим к нему, что является показателем точности решения. Это и понимается под словом “точность” в сделанных выводах.
3.5. Применение наилучшего аргумента
Для устранения вычислительных трудностей, связанных с плохой обусловленностью задачи (9), (4), преобразуем ее к наилучшему аргументу, дифференциал которого удовлетворяет соотношению
(12)
$\mathop {\left( {d\lambda } \right)}\nolimits^2 = \mathop {\left( {d\varepsilon } \right)}\nolimits^2 + \mathop {\left( {d\omega } \right)}\nolimits^2 + \mathop {\left( {dt} \right)}\nolimits^2 .$(13)
$\frac{{d\varepsilon }}{{d\lambda }} = \frac{{A\sigma _{0}^{n}{{{(1 + \varepsilon )}}^{n}}}}{Q},\quad \frac{{d\omega }}{{d\lambda }} = \frac{{B\sigma _{0}^{n}{{{(1 + \varepsilon )}}^{n}}}}{Q},\quad \frac{{dt}}{{d\lambda }} = \frac{{{{{(1 - {{\omega }^{r}})}}^{n}}}}{Q},$Начальные условия (4) для преобразованной системы (13) примут вид
Задача (13), (14) обладает рядом вычислительных преимуществ по сравнению с исходной задачей (9), (4) (см. [9]).
1. Преобразованная к наилучшему аргументу (12) задача (13), (14) является наилучшим образом обусловленной.
2. Квадратичная норма правой части системы (13) равна единице. Это дает возможность решать задачу (13), (14) любыми численными методами интегрирования задачи Коши, в том числе и явными.
3. Показатель жесткости системы (13) меньше, чем у исходной системы (9).
Преобразованная задача (13), (14) решалась с использованием тех же методов решения, что и исходная задача (9), (4). Результаты приведены в табл. 4 и 5. Обозначения аналогичны используемым в табл. 2 и 3. Используется постоянный шаг интегрирования ${{h}_{\lambda }}$ по аргументу $\lambda $, равный шагу ${{h}_{t}}$.
Таблица 4.
Решение задачи (13), (14), метод Эйлера с постоянным шагом
| ${{\sigma }_{0}}$, МПа | ${{h}_{\lambda }}$ | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$, ч | ${{t}_{c}}$, c | $N$ |
|---|---|---|---|---|---|---|
| 39.2 | 10–1 | 0.1035 | 0.9835 | 52.6472 | 0.069 | 1004 |
| 10–2 | 0.1046 | 0.9938 | 52.5843 | 0.174 | 10 024 | |
| 10–3 | 0.105 | 0.9976 | 52.5779 | 1.355 | 100 220 | |
| 10–4 | 0.1052 | 0.9997 | 52.5773 | 13.564 | 1 002 178 | |
| 49 | 10–1 | 0.1017 | 0.9659 | 25.8141 | 0.019 | 492 |
| 10–2 | 0.1041 | 0.9889 | 25.7516 | 0.055 | 4909 | |
| 10–3 | 0.1052 | 0.9995 | 25.7453 | 0.677 | 49 075 | |
| 10–4 | 0.1052 | 0.9996 | 25.7446 | 6.381 | 490 721 | |
| 58.8 | 10–1 | 0.1007 | 0.9563 | 14.4341 | 0.030 | 275 |
| 10–2 | 0.104 | 0.9876 | 14.372 | 0.049 | 2740 | |
| 10–3 | 0.105 | 0.9978 | 14.3657 | 0.217 | 27 384 | |
| 10–4 | 0.1052 | 0.9994 | 14.3651 | 3.555 | 273 815 | |
| 78.5 | 10–1 | 0.0979 | 0.9303 | 5.7891 | 0.021 | 110 |
| 10–2 | 0.1028 | 0.9764 | 5.7284 | 0.031 | 1092 | |
| 10–3 | 0.1047 | 0.9951 | 5.7221 | 0.178 | 10 908 | |
| 10–4 | 0.1052 | 0.9992 | 5.7215 | 1.479 | 109 059 |
Таблица 5.
Решение задачи (13), (14), метод Рунге–Кутты с постоянным шагом
| ${{\sigma }_{0}}$, МПа | ${{h}_{\lambda }}$ | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$, ч | ${{t}_{c}}$, c | $N$ |
|---|---|---|---|---|---|---|
| 39.2 | 10–1 | 0.1037 | 0.9854 | 52.5772 | 0.103 | 1002 |
| 10–2 | 0.1043 | 0.9905 | 52.5772 | 0.565 | 10 021 | |
| 10–3 | 0.105 | 0.9978 | 52.5772 | 5.191 | 100 217 | |
| 10–4 | 0.1052 | 0.9995 | 52.5772 | 50.046 | 1 002 174 | |
| 49 | 10–1 | 0.1008 | 0.9577 | 25.7444 | 0.065 | 490 |
| 10–2 | 0.1046 | 0.9935 | 25.7446 | 0.305 | 4907 | |
| 10–3 | 0.1048 | 0.9958 | 25.7446 | 2.676 | 49 071 | |
| 10–4 | 0.1052 | 0.9991 | 25.7446 | 25.137 | 490 717 | |
| 58.8 | 10–1 | 0.0989 | 0.9398 | 14.3646 | 0.048 | 273 |
| 10–2 | 0.1045 | 0.9929 | 14.365 | 0.2014 | 2738 | |
| 10–3 | 0.105 | 0.9975 | 14.365 | 1.433 | 27 381 | |
| 10–4 | 0.1051 | 0.9985 | 14.365 | 14.246 | 273 811 | |
| 78.5 | 10–1 | 0.1026 | 0.9746 | 5.7214 | 0.060 | 109 |
| 10–2 | 0.1026 | 0.975 | 5.7214 | 0.097 | 1090 | |
| 10–3 | 0.1044 | 0.992 | 5.7214 | 0.577 | 10 905 | |
| 10–4 | 0.105 | 0.9972 | 5.7214 | 5.687 | 109 055 |
Для метода продолжения решения по наилучшему аргументу кривые ползучести имеют аналогичный с представленными на фиг. 2 вид.
Результаты решения преобразованной задачи (13), (14) показывают следующее.
1. Полученные решения непреобразованной и преобразованной задач отличаются незначительно.
2. Для преобразованной задачи удается ближе подойти к моменту разрушения, чем для исходной.
3. Количество шагов по аргументу $\lambda $ возрастает вдвое, что увеличивает время счета втрое при ${{h}_{t}} = {{h}_{\lambda }}$. При этом стоит указать, что уже при шаге ${{h}_{\lambda }} = {{10}^{{ - 3}}}$ удается получить результаты, аналогичные по точности результатам решения исходной задачи при ${{h}_{t}} = {{10}^{{ - 4}}}$. То есть переход к аргументу $\lambda $ дает возможность использовать больший шаг интегрирования, уменьшая время счета до трех раз.
4. Как для исходной задачи, так и для преобразованной задачи (13), (14), при увеличении начального напряжения наблюдается увеличение погрешности.
3.6. Применение модифицированного наилучшего аргумента
При всех достоинствах задача (13), (14) имеет и ряд недостатков:
1) система уравнений (13) имеет значительно более сложный вид по сравнению с исходной (9). Это может перекрывать все вычислительные преимущества, указанные выше;
2) во многих случаях длина интервала изменения аргумента $\lambda $ значительно больше, чем для аргумента $t$. Это приводит к увеличению количества шагов по независимой переменной для преобразованной задачи (13), (14) при использовании постоянного шага интегрирования.
Усложнение вида преобразованной задачи (13), (14) по сравнению с исходной задачей (9), (4) приводит к значительному увеличению времени счета, в особенности для метода Рунге–Кутты. В статьях [13] и [14] для упрощения вида преобразованной задачи вводится модифицированный наилучший аргумент и показаны его преимущества при решении задач ползучести для упрочняющихся конструкций. Используем модифицированный аргумент и для задачи (9), (4).
По аналогии с результатами статьи [14], для задачи (9), (4) используем модифицированный наилучший аргумент вида
(15)
$\mathop {\left( {d\kappa } \right)}\nolimits^2 = \mathop {\left( {d\varepsilon } \right)}\nolimits^2 + \mathop {\left( {d\omega } \right)}\nolimits^2 + \mathop {\left( {\frac{{dt}}{{{{{(1 - {{\omega }^{r}})}}^{n}}}}} \right)}\nolimits^2 .$В результате получим преобразованную систему
(16)
$\frac{{d\varepsilon }}{{d\kappa }} = \frac{{A\sigma _{0}^{n}{{{(1 + \varepsilon )}}^{n}}}}{{{{Q}_{1}}}},\quad \frac{{d\omega }}{{d\kappa }} = \frac{{B\sigma _{0}^{n}{{{(1 + \varepsilon )}}^{n}}}}{{{{Q}_{1}}}},\quad \frac{{dt}}{{d\kappa }} = \frac{{{{{(1 - {{\omega }^{r}})}}^{n}}}}{{{{Q}_{1}}}},$Видно, что система (16), за счет упрощения знаменателей, имеет более простой вид по сравнению с системой (13). Но если для задачи (13), (14) доказана ее наилучшая обусловленность, то для задачи (16), (14) аналогичных утверждений нет. В общем случае неизвестно, насколько хуже она обусловлена.
В статьях [13] и [14] для оценки обусловленности задач, преобразованных к модифицированному наилучшему аргументу, предложена и апробирована мера обусловленности, которая для задачи (16), (14) записывается в виде
(17)
$\Xi = \mathop {max}\limits_{\kappa \in [0,\kappa *]} \left| {\int\limits_0^\kappa {\frac{{({{{(1 - {{\omega }^{r}}(\xi ))}}^{n}} - 1)}}{{\sqrt {1 + ({{A}^{2}} + {{B}^{2}})\sigma _{0}^{{2n}}{{{[1 + (A{\text{/}}B) \cdot \omega (\xi )]}}^{{2n}}}} }}d\xi } } \right|$(18)
$\overline \Xi = \mathop {max}\limits_{\kappa \in [0,\kappa {\text{*}}]} \left| {\frac{1}{{\kappa {\text{*}}}}\int\limits_0^\kappa {\frac{{({{{(1 - {{\omega }^{r}}(\xi ))}}^{n}} - 1)}}{{\sqrt {1 + ({{A}^{2}} + {{B}^{2}})\sigma _{0}^{{2n}}{{{[1 + (A{\text{/}}B) \cdot \omega (\xi )]}}^{{2n}}}} }}d\xi } } \right|.$(19)
$\kappa {\kern 1pt} * = \int\limits_0^1 {\frac{{\sqrt {1 + ({{A}^{2}} + {{B}^{2}})\sigma _{0}^{{2n}}{{{[1 + (A{\text{/}}B) \cdot \omega ]}}^{{2n}}}} }}{{B\sigma _{0}^{n}{{{[1 + (A{\text{/}}B) \cdot \omega ]}}^{n}}}}d\omega } .$В общем случае мера обусловленности $\Xi $ и ее осредненная величина $\overline \Xi $ могут принимать значения в интервале от 0 до бесконечности. При наилучшей обусловленности $\Xi = \overline \Xi = 0$, в остальных случаях эти величины строго положительные.При увеличении значения $\Xi $ или $\overline \Xi $ обусловленность преобразованной задачи становится хуже. Если величина меры обусловлености или ее осредненного значения стремятся к бесконечности, то используемый аргумент продолжения решения $\kappa $ приводит к плохо обусловленной начальной задаче. Для оценки обусловленности конкретной задачи удобнее использовать осредненную величину $\overline \Xi $. Если имеются две задачи, параметризованные аргументами ${{\kappa }_{1}}$ и ${{\kappa }_{2}}$, для которых имеются осредненные значения меры обусловленности $\mathop {\overline \Xi }\nolimits_1 $ и $\mathop {\overline \Xi }\nolimits_2 $, то при $\mathop {\overline \Xi }\nolimits_1 > \mathop {\overline \Xi }\nolimits_2 $ задача, соответствующая аргументу ${{\kappa }_{2}}$, обусловлена лучше, при $\mathop {\overline \Xi }\nolimits_1 < \mathop {\overline \Xi }\nolimits_2 $ – хуже, а при $\mathop {\overline \Xi }\nolimits_1 = \mathop {\overline \Xi }\nolimits_2 $ – обусловленность обеих задач одинакова.
Но мера обусловленности задачи (16), (14) вычисляется в виде (17) или (18) только в том случае, если существуют входящие в формулы (17) и (18) интегралы. Для установления этого факта все точки рассматриваемой области, в которой ищется решение, разделялись на точки локальной эквивалентности, точки локальной неэквивалентности I и II рода, в которых отклонение между направлениями отсчета аргументов $\lambda $ и $\kappa $ при переходе от точки к точке изменяется на бесконечно малую величину, конечную величину и бесконечно большую величину соответственно [13]. В случае, когда все точки рассматриваемой области являются точками локальной эквивалентности или локальной неэквивалентности I рода, мера обусловленности может быть выражена в виде интеграла Римана вида (17) или (18).
Докажем, что справедлива следующая
Теорема 3. Для задачи (16), (14) все точки области $V = \left\{ {\left( {\varepsilon ,\omega ,t} \right)\,|\,0 \leqslant \left. {\varepsilon \leqslant \varepsilon *,\;0 \leqslant \omega \leqslant 1,\;0 \leqslant t \leqslant t{\text{*}}} \right\}} \right.$ являются точками локальной эквивалентности.
Доказательство. Для аргумента (15) функция $\Psi (\omega ) = {{(1 - {{\omega }^{r}})}^{n}}$ не зависит от $\varepsilon $ и неявно зависит от времени $t$. Функция $\Psi (\omega )$ не обращается в ноль ни в одной точке рассматриваемого временного интервала, кроме момента разрушения $\omega = 1$ при $t = t{\text{*}}$ (с учетом накладываемых на параметры $n$ и $r$ ограничений (7)). Таким образом, согласно теореме 3.1 статьи [13], все точки области $V$, для которых $t \in [0,t*)$, являются точками локальной эквивалентности.
Покажем, что и в правой граничной точке имеет место локальная эквивалентность. Используем для этого результаты теорем 3.4 и 3.5 о локальной эквивалентности точек рассматриваемой области V, доказанных в статье [13].
При $t = t{\text{*}}$ достаточно доказать локальную эквивалентность только в ее левой полуокрестности, т.е. равенство нулю предела
Применяя формулу Тейлора, распишем данный предел в виде
В левой полуокрестности точки $t = t{\text{*}}$ для предыдущего предела справедливы следующие предельные соотношения, которые следуют из вида функции $\Psi (\omega )$ и уравнений исходной системы (9)
Продолжая разложение в ряд Тейлора, приходим к следующему равенству:
(20)
$\begin{gathered} \mathop {lim}\limits_{\omega \to 1 - 0} \mathop {lim}\limits_{\Delta t \to + 0} \mathop {\left| {\frac{{{{d}_{1}}}}{{1 - {{\omega }^{r}}}} + \cdots + \frac{{{{d}_{N}}}}{{{{{(1 - {{\omega }^{r}})}}^{{N + (N - 1)n}}}}}\Delta {{t}^{{N - 1}}} + \; \cdots {\kern 1pt} } \right|}\nolimits^{ - 1} = \\ = \;\mathop {lim}\limits_{\omega \to 1 - 0} \mathop {lim}\limits_{\Delta \omega \to + 0} \mathop {\left| {\frac{{d_{1}^{'}}}{{1 - {{\omega }^{r}}}} + \cdots + \frac{{d_{N}^{'}}}{{{{{(1 - {{\omega }^{r}})}}^{N}}}}\Delta {{\omega }^{{N - 1}}} + \cdots {\kern 1pt} } \right|}\nolimits^{ - 1} , \\ \end{gathered} $Рассмотрим предел, стоящий справа в равенстве (20). Будем полагать вначале, что параметр $n$ не принимает целых значений. В этом случае в выражении для коэффициентов $d_{N}^{'}$ произведение $\prod\nolimits_{k = 1}^N {(n - k + 1)} $ при $N \to \infty $ стремится к бесконечности со скоростью большей, чем экспоненциальная. Произведение ${{r}^{N}}\prod\nolimits_{k = 1}^N {(n - k + 1)} $ также будет возрастать по абсолютной величине при любых значениях $r > 0$ при увеличении $N$. Из этого следует, что коэффициенты $d_{N}^{'}$ при $N \to \infty $ стремятся к бесконечности со скоростью, большей экспоненциальной.
Если параметр $n$ принимает только целые значения, то начиная с некоторого номера ${{N}_{ * }} = n + 1$ все коэффициенты $d_{N}^{'}$ равны нулю.
Аналогичные утверждения справедливы и для коэффициентов ${{d}_{N}}$.
Учитывая это, покажем, что вне зависимости от порядка взятия пределов правая и левая части равенства (20) стремятся к нулю.
Рассмотрим предел, стоящий справа в равенстве (20).
При целых значениях параметра $n$ подмодульное выражение рассматриваемого предела будет конечной суммой. Устремляя вначале $\Delta \omega \to + 0$, получаем предел
Рассмотрим случай, когда параметр $n$ не принимает целых значений. Устремим вначале $\Delta \omega \to + 0$. Как было показано выше, коэффициенты $d_{N}^{'}$ стремятся к бесконечности со скоростью, превышающей экспоненциальную. При фиксированном значении $\omega $ слагаемые рассматриваемого предела будут монотонно возрастающими по абсолютной величине. Отсюда следует, что предел
Если же вначале устремить $\omega \to 1 - 0$ или рассматривать одновременное стремление $\Delta \omega \to + 0$ и $\omega \; \to \;1 - 0$, то качественно ситуация не изменится – предел, стоящий справа в равенстве (20), будет равен нулю. И это будет справедливо для всех значений параметров $A$, $B$, $n$, $r$, удовлетворяющих условиям (7).
Используя аналогичные рассуждения можно доказать стремление к нулю предела, стоящего слева в равенстве (20), вне зависимости от порядка взятия пределов.
Таким образом, доказано, что все точки области $V$ являются точками локальной эквивалентности.
Выполнение условий теоремы 3 гарантирует, что интегралы, входящие в соотношения для мер обусловленности (17) и (18), существуют и могут быть вычислены. В табл. 6 приводятся значения меры обусловленности $\Xi $ и ее осредненного значения $\overline \Xi $.
Таблица 6.
Основные теоретико-экспериментальные данные для стали Х18Н10Т
| ${{\sigma }_{0}}$, МПа | $\Xi $ | $\kappa {\text{*}}$ | $\overline \Xi $ |
|---|---|---|---|
| 39.2 | 53.7964 | 106.374 | 0.5057 |
| 49 | 26.3415 | 52.0863 | |
| 58.8 | 14.6981 | 29.0632 | |
| 78.5 | 5.8541 | 11.5755 |
Для задачи (16), (14) величина $\overline \Xi = 0.5057$, т.е. близка к нулю, и слабо зависит от значения начального напряжения ${{\sigma }_{0}}$. Это говорит о том, что рассматриваемая задача обусловлена хуже, чем задача (13), (14), но разница в обусловленности незначительна.
Замечание 2. Для иллюстрации применения введенной меры обусловленности вычислим ее значение для исходной задачи (9), (4). Как сказано выше, эта задача плохо обусловлена. Это следует из того, что при $\omega = 1$ правая часть системы (9) теряет смысл. Покажем это, используя введенную в статье меру обусловленности. В статье авторов [13] дан и теоретически обоснован способ нахождения меры обусловленности, состоящий в вычислении суммарной длины отклонения между векторами, в направлении которых отсчитываются рассматриваемый аргумент и наилучший. Выберем в качестве аргумента продолжения решения $\mu $ время $t$. Направляющий вектор в каждой точке интегральной кривой для аргумента $\mu $ будет равен $d\bar {\mu } = \mathop {\left( {0\;0\;dt} \right)}\nolimits^{\text{т}} $, а для наилучшего аргумента – $d\bar {\lambda } = \mathop {\left( {d\varepsilon \;d\omega \;dt} \right)}\nolimits^{\text{т}} $. В этом случае мера обусловленности для аргумента $t$ примет вид интеграла
Преобразованная к аргументу $\kappa $ задача (16), (14) решалась с использованием тех же методов решения, что и исходная задача (9), (4). Результаты приведены в табл. 7 и 8. Обозначения аналогичны используемым в табл. 2 и 3. Используется постоянный шаг интегрирования ${{h}_{\kappa }}$ по аргументу $\kappa $, равный шагу ${{h}_{t}}$.
Таблица 7.
Решение задачи (16), (14), метод Эйлера с постоянным шагом
| ${{\sigma }_{0}}$, МПа | ${{h}_{\kappa }}$ | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$, ч | ${{t}_{c}}$, c | $N$ |
|---|---|---|---|---|---|---|
| 39.2 | 10–1 | 0.1052 | 0.999 | 52.6348 | 0.033 | 1063 |
| 10–2 | 0.1053 | 0.9999 | 52.583 | 0.179 | 10 637 | |
| 10–3 | 0.1053 | 1 | 52.5778 | 1.534 | 106 374 | |
| 10–4 | 0.1053 | 1 | 52.5773 | 14.409 | 1 063 739 | |
| 49 | 10–1 | 0.1053 | 0.9999 | 25.8021 | 0.022 | 521 |
| 10–2 | 0.1052 | 0.9998 | 25.7503 | 0.104 | 5208 | |
| 10–3 | 0.1053 | 1 | 25.7452 | 0.736 | 52 086 | |
| 10–4 | 0.1053 | 1 | 25.7446 | 6.955 | 520 862 | |
| 58.8 | 10–1 | 0.1049 | 0.9968 | 14.4225 | 0.019 | 290 |
| 10–2 | 0.1052 | 0.9998 | 14.3708 | 0.053 | 2906 | |
| 10–3 | 0.1053 | 1 | 14.3656 | 0.4 | 29 063 | |
| 10–4 | 0.1053 | 1 | 14.3651 | 3.931 | 290 632 | |
| 78.5 | 10–1 | 0.1043 | 0.9907 | 5.7789 | 0.017 | 115 |
| 10–2 | 0.1052 | 0.9993 | 5.7271 | 0.029 | 1157 | |
| 10–3 | 0.1053 | 0.9999 | 5.722 | 0.191 | 11 575 | |
| 10–4 | 0.1053 | 1 | 5.7215 | 1.636 | 115 755 |
Таблица 8.
Решение задачи (16), (14), метод Рунге–Кутты с постоянным шагом
| ${{\sigma }_{0}}$, МПа | ${{h}_{\kappa }}$ | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$, ч | ${{t}_{c}}$, c | $N$ |
|---|---|---|---|---|---|---|
| 39.2 | 10–1 | 0.1052 | 0.9992 | 52.5772 | 0.095 | 1063 |
| 10–2 | 0.1053 | 1 | 52.5772 | 0.593 | 10 637 | |
| 10–3 | 0.1053 | 1 | 52.5772 | 5.396 | 106 373 | |
| 10–4 | 0.1053 | 1 | 52.5772 | 53.209 | 1 063 739 | |
| 49 | 10–1 | 0.1051 | 0.9981 | 25.7446 | 0.063 | 520 |
| 10–2 | 0.1052 | 0.9999 | 25.7446 | 0.309 | 5208 | |
| 10–3 | 0.1053 | 1 | 25.7446 | 2.777 | 52 086 | |
| 10–4 | 0.1053 | 1 | 25.7446 | 26.187 | 520 862 | |
| 58.8 | 10–1 | 0.105 | 0.9974 | 14.365 | 0.048 | 290 |
| 10–2 | 0.1052 | 0.9999 | 14.365 | 0.196 | 2906 | |
| 10–3 | 0.1053 | 1 | 14.365 | 1.538 | 29 063 | |
| 10–4 | 0.1053 | 1 | 14.365 | 14.468 | 290 631 | |
| 78.5 | 10–1 | 0.1045 | 0.9923 | 5.7214 | 0.042 | 115 |
| 10–2 | 0.1052 | 0.9994 | 5.7214 | 0.095 | 1157 | |
| 10–3 | 0.1053 | 1 | 5.7214 | 0.622 | 11 575 | |
| 10–4 | 0.1053 | 1 | 5.7214 | 5.959 | 115 754 |
Для метода продолжения решения по модифицированному наилучшему аргументу кривые ползучести имеют аналогичный с представленными на фиг. 2 вид.
Из результатов решения, представленных в табл. 7 и 8, видно следующее.
1. Использование модифицированного аргумента продолжения решения позволяет подойти вплотную к моменту разрушения. Уже при шаге интегрирования, равном 0.1, удается получить результаты, сопоставимые с решением исходной задачи при шаге, равном 0.0001. А начиная с шага, равного 0.01, результаты сопоставимы с решением преобразованной к наилучшему аргументу задачи с шагом 0.0001 (при этом время счета сокращается дополнительно до 8 раз). Это показывает, что модифицированный наилучший аргумент позволяет использовать наибольший шаг интегрирования. Но стоит учитывать, что при малых шагах задача, преобразованная к модифицированному наилучшему аргументу, может требовать до двух и более раз большего времени вычислений.
2. В отличие от исходной задачи (9), (4) и преобразованной к наилучшему аргументу (13), (14) для модифицированного наилучшего аргумента решение менее зависимо от величины начального напряжения ${{\sigma }_{0}}$.
3. Введенная мера обусловленности $\Xi $ или $\overline \Xi $, значения которых даны в табл. 6, показывают, что обусловленность рассматриваемой задачи (16), (14) не намного меньше, чем наилучшая. Это подтверждается результатами расчета.
3.7. Применение переменного шага интегрирования
Для достижения высокой точности численного решения необходимо значительно измельчать шаг интегрирования, что делает применение постоянного шага малоэффективным. Кроме того, известно, что постоянный шаг интегрирования целесообразно использовать только для нежестких задач. В случае когда задача является плохо обусловленной или жесткой более эффективным становится применение переменного шага интегрирования. Существует множество методов построения неравномерных сеток, некоторые из них приведены в монографиях Н.Н. Калиткина и соавт. [40] и О.Б. Арушаняна, С.Ф. Залеткина [41]. Традиционно для смены шага используется правило Рунге, известное также как метод Рунге–Ромберга–Ричардсона.
Согласно правилу Рунге, локальную погрешность численного решения задачи Коши на $k$-м шаге, вычисленного с шагом $h$, можно записать в виде [41, с. 73–76]
где $\rho _{k}^{h}$ – главная часть локальной погрешности, вычисляемая по формуле(21)
$\rho _{k}^{h} = \frac{{{{{\left\| {{\mathbf{y}}_{k}^{{2h}} - {\mathbf{y}}_{k}^{h}} \right\|}}_{2}}}}{{{{2}^{p}} - 1}}.$Применяя значение (21) как оценку локальной погрешности, можно реализовать процедуру смены шага. Задавая параметр точности численного решения $\theta $ (величину предельной допустимой погрешности), при $\rho _{k}^{h} > \theta $ шаг уменьшается вдвое. Процедура уменьшения шага повторяется до момента, когда $\rho _{k}^{h} < \theta $. Если значение $\rho _{k}^{h} < \theta {\text{/}}{{2}^{p}}$, то шаг увеличивается вдвое [41, с. 82–84].
3.7.1. Исходная задача. В табл. 9 и 10 даны результаты решения задачи (9), (4) с переменным шагом интегрирования. Точность $\theta = {{10}^{{ - 4}}}$, величина начального шага ${{h}_{t}} = 1$ ч.
Таблица 9.
Решение задачи (9), (4), метод Эйлера с переменным шагом, точность $\theta = {{10}^{{ - 4}}}$
| ${{\sigma }_{0}}$, МПа | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$ | ${{t}_{c}}$, мc | $N$ |
|---|---|---|---|---|---|
| 39.2 | 0.1053 | 1 | 52.8426 | 19.589 | 454 |
| 49 | 0.1053 | 1 | 25.8743 | 18.121 | 449 |
| 58.8 | 0.1053 | 1 | 14.4379 | 17.408 | 451 |
| 78.5 | 0.1053 | 0.9999 | 5.7502 | 18.809 | 455 |
Таблица 10.
Решение задачи (9), (4), метод Рунге–Кутты с переменным шагом, точность $\theta = {{10}^{{ - 4}}}$
| ${{\sigma }_{0}}$, МПа | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$ | ${{t}_{c}}$, мc | $N$ |
|---|---|---|---|---|---|
| 39.2 | 0.1047 | 0.9949 | 52.5754 | 10.490 | 48 |
| 49 | 0.1042 | 0.9897 | 25.7446 | 9.595 | 45 |
| 58.8 | 0.1048 | 0.9961 | 14.3646 | 12.108 | 54 |
| 78.5 | 0.1044 | 0.9919 | 5.7213 | 9.921 | 38 |
Для всех методов с переменным шагом интегрирования кривые ползучести совпадают с представленными на фиг. 2.
Сравнивая полученные результаты с табл. 2 и 3, можно отметить следующее:
1) использование переменного шага интегрирования позволяет на порядки уменьшить время счета для рассматриваемой задачи, позволяя получить при этом решение заданной точности $\theta $;
2) при отсутствии аналитического решения правило Рунге является фактически единственным способом оценить погрешность решения. По результатам предыдущих работ [34], [14] можно отметить, что для задач ползучести средняя величина и медиана погрешности решения хорошо согласуются с оценкой локальной погрешности по правилу Рунге.
3.7.2. Наилучший аргумент. Для задач, преобразованных к аргументам $\lambda $ или $\kappa $, применение традиционной оценки погрешности по правилу Рунге (21) может быть сопряжено со значительным возрастанием количества вычислений вследствие существенных ограничений на изменение исходного аргумента $t$. Для предотвращения этого на практике часто используется модификация оценки (21), в которой не учитывается вносимая при вычислении значений $t$ погрешность. Эта модификация значительно упрощает вычислительный процесс, хотя и менее точна. В дальнейшем будем называть данный способ оценки локальной погрешности модифицированным правилом Рунге.
Замечание 3. Правомерность использования такой модификации можно доказать, используя упрощенную оценку локальной погрешности параметризованных задач, данную в работе [9, с. 64–65]. Например, для метода Эйлера традиционная оценка (21) и ее модификация отличаются не более чем на $\theta {\text{/}}2$.
Результаты применения модифицированного правила Рунге для задачи (13), (14), преобразованной к наилучшему аргументу $\lambda $, даны в табл. 11 и 12. Точность $\theta $ принималась равной ${{10}^{{ - 4}}}$, а начальный шаг ${{h}_{\lambda }} = 1$.
Таблица 11.
Решение задачи (13), (14), метод Эйлера с переменным шагом, точность $\theta = {{10}^{{ - 4}}}$
| ${{\sigma }_{0}}$, МПа | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$ | ${{t}_{c}}$, мc | $N$ |
|---|---|---|---|---|---|
| 39.2 | 0.1053 | 0.9999 | 52.983 | 15.831 | 266 |
| 49 | 0.1052 | 0.9998 | 25.9440 | 16.099 | 265 |
| 58.8 | 0.1052 | 0.9998 | 14.4752 | 16.084 | 265 |
| 78.5 | 0.1052 | 0.9997 | 5.7657 | 15.855 | 262 |
Таблица 12.
Решение задачи (13), (14), метод Рунге–Кутты с переменным шагом, точность $\theta = {{10}^{{ - 4}}}$
| ${{\sigma }_{0}}$, МПа | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$ | ${{t}_{c}}$, мc | $N$ |
|---|---|---|---|---|---|
| 39.2 | 0.1046 | 0.9941 | 52.5793 | 6.591 | 26 |
| 49 | 0.1044 | 0.9915 | 25.7458 | 5.959 | 24 |
| 58.8 | 0.0984 | 0.9346 | 14.3647 | 3.726 | 16 |
| 78.5 | 0.1043 | 0.9909 | 5.7215 | 6.060 | 23 |
Анализируя результаты решения задачи (13), (14) с переменным шагом интегрирования, можно отметить, что по сравнению с исходной задачей (9), (4) удалось следующее:
1) сократить время счета в среднем от 20 до 50 процентов. Но стоит отметить, что для задач ползучести имеет место и противоположный результат. В статье [42] показано, что усложнение вида преобразованной системы может приводить к увеличению времени счета по сравнению с исходной задачей;
2) уменьшить количество шагов при численном решении до 50 процентов. Для метода Рунге–Кутты это и привело к значительному сокращению времени счета. Для метода Эйлера уменьшение количества шагов на 40 процентов сопровождается уменьшением времени счета только на 20%, что объясняется усложнением вида преобразованной системы.
3.7.3. Модифицированный наилучший аргумент. Для задачи (16), (14) не удалось получить адекватного решения явным методом Эйлера с переменным шагом при точности $\theta = {{10}^{{ - 4}}}$. По этой причине значение $\theta $ было уменьшено до $3 \times {{10}^{{ - 5}}}$. Результаты решения при данном значении точности для метода Эйлера даны в табл. 13. Начальное значение шага ${{h}_{\kappa }} = 1$.
Таблица 13.
Решение задачи (16), (14), метод Эйлера с переменным шагом, точность $\theta = 3 \times {{10}^{{ - 5}}}$
| ${{\sigma }_{0}}$, МПа | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$ | ${{t}_{c}}$, мc | $N$ |
|---|---|---|---|---|---|
| 39.2 | 0.1053 | 1 | 53.0761 | 9.783 | 147 |
| 49 | 0.1049 | 0.9967 | 25.9773 | 9.548 | 151 |
| 58.8 | 0.1048 | 0.9959 | 14.5086 | 5.672 | 124 |
| 78.5 | 0.105 | 0.9976 | 5.7585 | 11.953 | 177 |
Для метода Рунге–Кутты точность $\theta $ принималась равной ${{10}^{{ - 4}}}$. Результаты даны в табл. 14. Начальный шаг ${{h}_{\kappa }}$ также принимался равным единице.
Таблица 14.
Решение задачи (16), (14), метод Рунге–Кутты с переменным шагом, точность $\theta = {{10}^{{ - 4}}}$
| ${{\sigma }_{0}}$, МПа | $\varepsilon _{n}^{ * }$ | $\omega _{n}^{ * }$ | $t_{n}^{ * }$ | ${{t}_{c}}$, мc | $N$ |
|---|---|---|---|---|---|
| 39.2 | 0.1048 | 0.9959 | 52.5772 | 5.330 | 28 |
| 49 | 0.0968 | 0.9200 | 25.7421 | 4.589 | 25 |
| 58.8 | 0.0997 | 0.9473 | 14.3646 | 3.617 | 21 |
| 78.5 | 0.1045 | 0.9923 | 5.7215 | 2.974 | 16 |
Для преобразованной задачи (16), (14) отметим следующее:
1) по сравнению с задачей (13), (14), преобразованной к наилучшему аргументу, сокращение количества шагов по аргументу $\kappa $ до 40–50 процентов и упрощение вида системы (16) приводят к уменьшению времени счета от 25 до 65 процентов (и до трех раз по сравнению с исходной задачей (9), (4)). Снижение времени счета происходит даже при повышении требуемой точности $\theta $;
2) ухудшение обусловленности по сравнению с наилучшей отражается в остановке процесса вычислений в более дальней от момента разрушения точке по сравнению с задачей (13), (14). В целом же полученные результаты хорошо согласуются с полученными ранее.
Замечание 4. Нельзя не отметить увеличение погрешности приближенных решений задач (13), (14) при ${{\sigma }_{0}} = 58.8$ МПа и (16), (14) при ${{\sigma }_{0}} = 49;\;58.8$ МПа, полученных методом Рунге–Кутты с переменным шагом, по сравнению с заданным значением $\theta $. Это можно частично отнести на счет дополнительной погрешности, вносимой в модифицированное правило Рунге. Но такое возрастание погрешности можно встретить и при использовании традиционного правила Рунге, в особенности для плохо обусловленных задач и задач с контрастными структурами. Из полученных результатов видно, что в ряде случаев при оценке локальной погрешности по правилу Рунге или его модификации можно получить завышенное или заниженное значение погрешности. Если для наилучшего аргумента есть некоторые теоретические результаты в данном направлении, то для модифицированного наилучшего аргумента таких результатов нет. Фактически применение модифицированного метода Рунге для задачи (16), (14) пока не обосновано. Получение уточненных оценок локальной погрешности для преобразованных задач является одним из приоритетных направлений исследований.
4. ВЫВОДЫ И ЗАКЛЮЧЕНИЕ
В статье разработаны новые методы численного решения задачи Коши для систем обыкновенных дифференциальных уравнений с одной предельной особой точкой на правой границе интервала изменения аргумента. Показано, что традиционные явные методы испытывают вычислительные трудности при решении подобных задач. Для явных методов это происходит за счет необходимости значительного уменьшения величины шага интегрирования вблизи предельной особой точки, что не позволяет близко подойти к правой границе при умеренных значениях шага интегрирования. Используемый для устранения недостатков применения явных схем метод продолжения решения по наилучшему аргументу из-за усложнения вида исходной системы уравнений также сталкивается со значительным увеличением вычислительной сложности, а значит и времени счета. Для устранения указанных недостатков наилучшей параметризации предложен новый вид модифицированного наилучшего аргумента, апробация которого проведена на задаче деформирования трубчатых образцов из нержавеющей стали Х18Н10Т под действием постоянной одноосной растягивающей нагрузки при постоянной температуре в условиях ползучести вплоть до разрушения. Для моделирования этой задачи используются уравнения теории структурных параметров Ю.Н. Работнова в форме задачи Коши для системы двух обыкновенных дифференциальных уравнений с одной предельной особой точкой в момент разрушения. Ранее доказано, что в общем случае аналитическое решение данной задачи получить нельзя. При этом ее численное решение сопряжено с рядом трудностей, связанных с плохой обусловленностью системы уравнений. Такая ситуация характерна для задач ползучести, что говорит о сложности их решения. В статье проводится всестороннее исследование особенностей численного решения задачи деформирования трубчатых образцов и даются рекомендации по выбору вида аргумента продолжения решения, а также оценка обусловленности полученной преобразованной задачи.
Использование традиционных явных методов Эйлера и Рунге–Кутты четвертого порядка точности с постоянным шагом для решения исходной задачи (9), (4) показало возможность получения ее приближенного решения, но только вне некоторой окрестности предельной особой точки, являющейся моментом разрушения. На практике необходимо подходить к моменту разрушения достаточно близко, что требует использования очень маленького шага интегрирования, что влечет увеличение вычислительных затрат. Это делает традиционные явные методы с постоянным шагом малоэффективными. Для устранения этого недостатка исходная задача преобразуется к наилучшему и модифицированному наилучшему аргументам.
Использование наилучшего аргумента позволяет получить наилучшим образом обусловленную систему уравнений, что дает возможность подойти к моменту разрушения значительно ближе и использовать больший шаг интегрирования по сравнению с исходной задачей. Это позволяет сократить время счета до трех раз. Модифицированный наилучший аргумент за счет упрощения вида преобразованной системы дает уменьшение времени счета еще до 9 раз, позволяя использовать наибольший шаг интегрирования из всех рассматриваемых методов.
Помимо постоянного шага интегрирования, использовался переменный, изменяемый по правилу Рунге для исходной задачи (9), (4) и его модификации для преобразованных задач (13), (14) и (16), (14). Для исходной задачи применение переменного шага позволяет уменьшить время счета на порядки и близко подойти к моменту разрушения. При отсутствии аналитического решения правило Рунге является фактически единственным способом контроля точности. Несмотря на эффективность применения переменного шага к исходной задаче, использование наилучшего и модифицированного наилучшего аргумента позволяет дополнительно сократить время счета до 50 процентов и до трех раз соответственно.
Использование метода продолжения решения по наилучшему аргументу позволяет значительно упростить вычислительный процесс, связанный с численным решением плохо обусловленных начальных задач. Но в ряде случаев сложность полученной преобразованной системы может перекрывать преимущества, связанные с устранением предельной особой точки. В этом случае хорошо работает модифицированный наилучший аргумент, упрощая вид преобразованной системы без значительного увеличения погрешности. Возможный рост погрешности позволяет оценить введенная мера обусловленности, которая отражает отклонение направления отсчета модифицированного аргумента от касательного. Для рассматриваемой задачи мера обусловленности мала, что говорит о близости обусловленности преобразованной к модифицированному аргументу задачи Коши к наилучшей.
Несмотря на многие преимущества метода продолжения решения по аргументам различного вида, также существуют его недостатки. Их устранение сопряжено с решением следующих задач.
1. В статье основное внимание уделяется конструированию модифицированных аргументов для задач ползучести. Определяющие уравнения ползучести явно содержат слагаемые или множители, содержащие обнуляющиеся знаменатели. Эта особенность записи исходной системы уравнений существенно используется при выборе модифицированного наилучшего аргумента. Однако предложенный метод конструирования модифицированных аргументов может быть без существенных изменений обобщен и на более широкий класс задач, хотя он и не является универсальным. Для примера рассмотрим следующую задачу Коши:
имеющую решение Данная задача, преобразованная к наилучшему аргументу $\lambda $, будет иметь видТаким образом, требуются дополнительные исследования, направленные на уточнение и развитие способов построения модифицированных наилучших аргументов для различных классов жестких начальных задач.
2. Открытым вопросом для преобразованных задач остается уточнение метода оценки локальной погрешности. Традиционное правило Рунге может значительно усложнять вычислительный процесс, давая завышенную точность. Этот недостаток устраняет предложенное авторами модифицированное правило Рунге, исключающее учет погрешности вычисления значений исходного аргумента $t$. Но данный подход обоснован только для наилучшего аргумента. Одними из приоритетных задач в данном направлении являются разработка метода оценивания локальной погрешности для метода продолжения решения по аргументам произвольного вида и разработка на его основе механизма смены шага интегрирования.
3. Для введенной меры обусловленности требуют уточнения интервалы хорошей и плохой обусловленности преобразованной системы уравнений. На данный момент известно, что при значениях меры обусловленности или ее осреднения близких к нулю преобразованная система имеет хорошую обусловленность, а при стремлении меры обусловленности к бесконечности – плохую. Но мера обусловленности и ее осредненное значение могут принимать любые положительные значения. В дальнейших исследованиях необходимо каждому значению меры обусловленности сопоставить степень обусловленности рассматриваемой преобразованной системы.
Авторы выражают благодарность Александре Дмитриевне Бересневой за помощь в проведении вычислений.
Список литературы
Lahaye M.E. Une metode de resolution d’une categorie d’equations transcendentes // Compter Rendus hebdomataires des seances de L’Academie des sciences. 1934. V. 198. №. 21. P. 1840–1842.
Lahaye M.E. Solution of system of transcendental equations // Acad. Roy. Belg. Bull. Cl. Sci. 1948. V. 5. P. 805–822.
Давиденко Д.Ф. Об одном новом методе численного решения систем нелинейных уравнений // Докл. АН СССР. 1953. Т. 88. № 4. С. 601–602.
Давиденко Д.Ф. О приближенном решении систем нелинейных уравнений // Украинский матем. ж. 1953. Т. 5. № 2. С. 196–206.
Ворович И.И., Зипалова В.Ф. К решению нелинейных краевых задач теории упругости методом перехода к задаче Коши // Прикл. матем. и механ. 1965. Т. 29. Вып. 5. С. 894–901. (Перевод: Vorovich I.I., Zipalova V.F. On the solution of nonlinear boundary value problems of the theory of elasticity by a method of transformation to an initial value Cauchy problem // J. of Appl. Math. and Mech. 1965. V. 29. № 5. P. 1055–1063.).
Рикс Э. Применение метода Ньютона к задаче упругой устойчивости // Прикл. механ. 1972. № 5. С. 204–210. (Исходная статья: Riks E. The application of Newton’s method to the problem of elastic stability // Trans. ASME. 1972. № 4. P. 1060–1065.).
Григолюк Э.И., Шалашилин В.И. Метод продолжения решения по параметру в нелинейных задачах механики твердого деформируемого тела. М.: Наука, 1988. (Перевод: Grigolyuk E.I., Shalashilin V.I. Problems of Nonlinear Deformation: The Continuation Method Applied to Nonlinear Problems in Solid Mechanics. Dordrecht: Kluwer Academic Publishers, 1991.).
Allgower E.L., Georg K. Introduction to numerical continuation methods. Berlin Heidelberg: Springer-Verlag, 1990.
Шалашилин В.И., Кузнецов Е.Б. Метод продолжения решения по параметру и наилучшая параметризация в прикладной математике и механике. М.: Эдиториал УРСС, 1999. (Перевод: Shalashilin V.I., Kuznetsov E.B. Parametric Continuation and Optimal Parametrization in Applied Mathematics and Mechanics. Dordrecht/Boston/London: Kluwer Academic Publishers, 2003.).
Кузнецов Е.Б. Многомерная параметризация и численное решение систем нелинейных уравнений // Ж. вычисл. матем. и матем. физ. 2010. Т. 50. № 2. С. 255–267. (Перевод: Kuznetsov E.B. Multidimensional parametrization and numerical solution of systems of nonlinear equations // Comput. Math. Math. Phys. 2010. V. 50. Iss. 2. P. 244–255.).
Красников С.Д., Кузнецов Е.Б. Численное продолжение решения в особых точках коразмерности единица // Ж. вычисл. матем. и матем. физ. 2015. Т. 55. № 11. С. 1835–1856. (Перевод: Kuznetsov E.B., Krasnikov S.D. Numerical continuation of solution at singular points of codimension one // Comput. Math. Math. Phys. 2015. V. 55. Iss. 11. P. 1802–1822.).
Красников С.Д., Кузнецов Е.Б. Численное продолжение решения в особых точках высокой коразмерности для систем нелинейных алгебраических или трансцендентных уравнений // Ж. вычисл. матем. и матем. физ. 2016. Т. 56. № 9. С. 1571–1585. (Перевод: Kuznetsov E.B., Krasnikov S.D. Numerical continuation of solution at a singular point of high codimension for systems of nonlinear algebraic or transcendental equations // Comput. Math. Math. Phys. 2016. V. 56. Iss. 9. P. 1551–1564.).
Кузнецов Е.Б., Леонов С.С. Параметризация задачи Коши для систем обыкновенных дифференциальных уравнений с предельными особыми точками // Ж. вычисл. матем. и матем. физ. 2017. Т. 57. № 6. С. 934–957. (Перевод: Kuznetsov E.B., Leonov S.S. Parametrization of the Cauchy problem for systems of ordinary differential equations with limiting singular points // Comput. Math. Math. Phys. 2017. V. 57. Iss. 6. P. 931–952.).
Кузнецов Е.Б., Леонов С.С. Примеры параметризации задачи Коши для систем обыкновенных дифференциальных уравнений с предельными особыми точками // Ж. вычисл. матем. и матем. физ. 2018. Т. 58. № 6. С. 914–933. (Перевод: Kuznetsov E.B., Leonov S.S. Examples of parametrization of the Cauchy problem for systems of ordinary differential equations with limiting singular points // Comput. Math. Math. Phys. 2018. V. 58. Iss. 6. P. 881–897.).
Калиткин Н.Н., Пошивайло И.П. Решение задачи Коши для жестких систем с гарантированной точностью методом длины дуги // Матем. моделирование. 2014. Т. 26. № 7. С. 3–18. (Перевод: Kalitkin N.N., Poshivaylo I.P. Arc length method of solving Cauchy problem with guaranteed accuracy for stiff systems // Math. Models and Comput. Simul. 2015. V. 7. Iss. 1. P. 24–35.).
Белов А.А., Калиткин Н.Н. Численные методы решения задач Коши с контрастными структурами // Моделирование и анализ информ. систем. 2016. Т. 23. № 5. С. 529–538.
Белов А.А., Калиткин Н.Н. Особенности расчета контрастных структур в задачах Коши // Матем. моделирование. 2016. Т. 28. № 10. С. 97–109. (Перевод: Belov A.A., Kalitkin N.N. Features of calculating contrast structures in the Cauchy problem // Math. Models and Comput. Simul. 2017. V. 9. Iss. 3. P. 281–291.).
Белов А.А., Калиткин Н.Н. Экономичные методы численного интегрирования задачи Коши для жестких систем ОДУ // Дифференц. ур-ния. 2019. Т. 55. № 7. С. 907–918. (Перевод: Belov A.A., Kalitkin N.N. Efficient Numerical Integration Methods for the Cauchy Problem for Stiff Systems of Ordinary Differential Equations // Differential Equations. 2019. V. 55. Iss. 7. P. 871–883.).
Москаленко Л.П. Методика исследования устойчивости пологих ребристых оболочек на основе метода продолжения решения по наилучшему параметру // Вестн. гражданских инженеров. 2011. № 4(29). С. 161–164.
Semenov A.A. Strength and stability of geometrically nonlinear orthotropic shell structures // Thin-Walled Structures. 2016. V. 106. P. 428–436.
May S., Vignollet J., de Borst R. A new arc-length control method based on the rates of the internal and the dissipated energy // Eng. Comput. 2016. V. 33. Iss. 1. P. 100–115.
Wang X., Ma T.-B., Ren H.-L., Ning J.-G. A local pseudo arc-length method for hyperbolic conservation laws // Acta Mechanica Sinica. 2015. V. 30. № 6. P. 956–965.
Григолюк Э.И., Лопаницын Е.А. Конечные прогибы, устойчивость и закритическое поведение тонких пологих оболочек. М.: Изд-во МАМИ, 2004.
Гаврюшин С.С., Барышникова О.О., Борискин О.Ф. Численный анализ элементов конструкций машин и приборов. М.: Изд-во МГТУ им. Н.Э. Баумана, 2014.
Работнов Ю.Н. Ползучесть элементов конструкций. М.: Наука, 2014. (Перевод: Rabotnov Yu.N. Creep problems in structural members. Amsterdam/London: North-Holland Publishing Company, 1969.).
Малинин Н.Н. Прикладная теория пластичности и ползучести. М.: Машиностр., 1975.
Локощенко А.М. Ползучесть и длительная прочность металлов. М.: Физматлит, 2016. (Перевод: Lokoshchenko A.M. Creep and Creep Rupture of Metals. Boca Raton–London–New York: CRC Press Taylor & Francis Group, 2018.).
Локощенко А.М., Шестериков С.А. Ползучесть // Итоги науки. Серия “Механика”. М.: ВИНИТИ, 1965. С. 177–227.
Шестериков С.А., Локощенко А.М. Ползучесть и длительная прочность металлов // Итоги науки. Серия “Механика деформируемого твердого тела”. М.: ВИНИТИ, 1980. Т. 13. С. 3–104.
Локощенко А.М. Ползучесть и длительная прочность металлов в агрессивных средах (обзор) // Физизико-химическая механика материалов. 2001. № 4. С. 27–41.
Локощенко А.М. Длительная прочность металлов при сложном напряженном состоянии (обзор) // Изв. РАН. Механика твердого тела. 2012. № 3. С. 116–136. (Перевод: Lokoshchenko A.M. Long-Term Strength of Metals in Complex Stress State (a Survey) // Mechanics of Solids. 2012. V. 47. № 3. P. 357–372.).
Локощенко А.М. Применение кинетической теории при анализе длительного высокотемпературного разрушения металлов в условиях сложного напряженного состояния (обзор) // Прикладная механика и техническая физика. 2012. Т. 53. № 4. С. 149–164. (Перевод: Lokoshchenko A.M. Application of kinetic theory to the analysis of high-temperature creep rupture of metals under complex stress (review) // Journal of Applied Mechanics and Technical Physics. 2012. V. 53. Iss. 4. P. 599–610.).
Локощенко А.М. Результаты исследований ползучести и длительной прочности металлов в Научно-исследовательском институте механики Московского государственного университета им. М.В. Ломоносова (к юбилею Ю.Н. Работнова) // Прикладная механика и техническая физика. 2014. Т. 55. № 1. С. 144–165. (Перевод: Lokoshchenko A.M. Results of studying creep and long-term strength of metals at the Institute of Mechanics at the Lomonosov Moscow State University (To Yu.N. Rabotnov’s Anniversary) // Journal of Applied Mechanics and Technical Physics. 2014. V. 55. Iss. 1. P. 118–135.).
Кузнецов Е.Б., Леонов С.С. Методика выбора функций определяющих уравнений ползучести и длительной прочности с одним скалярным параметром поврежденности // Прикладная механика и техническая физика. 2016. Т. 57. № 2. С. 202–211. (Перевод: Kuznetsov E.B., Leonov S.S. Technique for selecting the functions of the constitutive equations of creep and long-term strength with one scalar damage parameter // Journal of Applied Mechanics and Technical Physics. 2016. V. 57. Iss. 2. P. 369–377.).
Качанов Л.М. О времени разрушения в условиях ползучести // Изв. АН СССР. Отд. технических наук. 1958. № 8. С. 26–31.
Работнов Ю.Н. О механизме длительного разрушения // Вопросы прочности материалов и конструкций: сб. статей. Изд-во АН СССР. Москва, 1959. С. 5–7.
Кузнецов Е.Б., Леонов С.С. Об аналитическом решении одной задачи ползучести // Ж. Средневолжского математического общества. 2018. Т. 20. № 3. С. 282–294.
Локощенко А.М., Шестериков С.А. Методика описания ползучести и длительной прочности при чистом растяжении // Прикладная механика и техническая физика. 1980. Т. 21. № 3. С. 155–159. (Перевод: Lokoshchenko A.M., Shesterikov S.A. Method for description of creep and long-term strength with pure elongation // Journal of Applied Mechanics and Technical Physics. 1980. V. 21. Iss. 3. P. 414–417.).
Калиткин Н.Н. Численные методы. СПб.: БХВ-Петербург, 2011.
Калиткин Н.Н., Альшин А.Б., Альшина Е.А., Рогов Б.В. Вычисления на квазиравномерных сетках. М.: Физматлит, 2005.
Арушанян О.Б., Залеткин С.Ф. Численное решение обыкновенных дифференциальных уравнений на Фортране. М.: Изд-во МГУ, 1990.
Кузнецов Е.Б., Леонов С.С. Чистый изгиб балки из разномодульного материала в условиях ползучести // Вестник Южно-Уральского государственного университета. Серия “Математическое моделирование и программирование”. 2013. Т. 6. № 4. С. 26–38.
Дополнительные материалы отсутствуют.
Инструменты
Журнал вычислительной математики и математической физики
Пн.-пт.: 10.00-19.00