Известия РАН. Механика твердого тела, 2022, № 2, стр. 72-81

1. Введение. Многие детали машин при эксплуатации подвергаются немонотонной нерегулярной нагрузке. Подобные воздействия испытывают детали самолетов, вертолетов, экскаваторов, ветрогенераторов, наземного транспорта. Традиционным подходом к оценке усталостной прочности в подобных случаях, используемым многими конструкторскими бюро, является метод подсчета циклов, который подразумевает разбиение истории нагружения на отдельные полуциклы и определение вклада каждого из них в поврежденность без учета их очередности [1–3]. Однако очередность циклов нагрузки может значительно влиять на время до разрушения [4, 5]. Существует эволюционный подход к описанию накопления усталостных микроповреждений, который позволяет естественным образом учесть это влияние. Он заключается в записи эволюционного уравнения для параметра поврежденности в форме $\dot {D}(t) = F({\mathbf{\sigma }}(\tau ),t)$, где ${\mathbf{\sigma }}(t)$ есть история изменения тензора напряжений, F – функционал процесса. Такой подход естественным образом применим и в случае многоосного непропорционального нагружения. Данное эволюционное уравнение интегрируется по времени, пока $D$ не достигнет единицы, что означает разрушение. Некоторые эволюционные модели предложены в [6, 7]. Модель Оттосена [7] предназначена для описания многоцикловой усталости, не связана с уравнениями движения и допускает точное решение на регулярном циклическом нагружении, что существенно упрощает процедуру ее идентификации. В настоящей работе предпринимается усовершенствование этой модели для возможности учета специфических эффектов, наблюдаемых в различных программах нерегулярного циклического нагружения [5, 8, 9].

2. Модель. Предлагаются следующие уравнения:

В основе данной модели лежит понятие поверхности выносливости в пространстве напряжений, задаваемой уравнением $\beta = 0$. Пока напряженное состояние лежит внутри поверхности ($\beta \leqslant 0$, рис. 1,a), прироста параметра поврежденности $D$ не происходит, а когда напряженное состояние выходит за ее пределы ($\beta > 0$, рис. 1,b), параметр $D$ начинает расти в зависимости от расстояния до поверхности и направления изменения напряжений. Кроме того, поверхность выносливости может менять размер и перемещаться, приспосабливаясь к нагружению. Местоположение поверхности задается тензорнозначным параметром α.

Выражение для $\beta $, приведенное в (2.1), имеет следующую структуру. Первое слагаемое определяет вклад девиаторной части напряженного состояния и является расстоянием между центром поверхности выносливости α и девиатором напряжений s = σ – – $\frac{1}{3}{{I}_{1}}{\mathbf{I}}$, где ${{I}_{1}} = {{\sigma }_{{ii}}}$ есть первый инвариант тензора напряжений, I – единичный тензор. Норма симметричного тензора определяется как ${\text{||A||}} = \sqrt {{\mathbf{A}}:{\mathbf{A}}} $. Второе слагаемое учитывает вклад гидростатической части и имеет вид материальной функции $g$ от первого инварианта тензора напряжений ${{I}_{1}}$. Третье слагаемое задает радиус поверхности выносливости, который зависит от ее положения α в пространстве напряжений. Последнее слагаемое описывает сокращение радиуса поверхности выносливости как функцию максимального размаха напряжений в течение истории нагружения до текущего момента времени $\Delta {{\sigma }_{{{\text{max}}}}}(t) = \mathop {\max }\limits_{{{t}_{1}} < {{t}_{2}} \leqslant t} \left( {{\text{||}}{\mathbf{\sigma }}({{t}_{1}}) - {\mathbf{\sigma }}({{t}_{2}}){\text{||}}} \right)$. Обезразмеривающий коэффициент σ_–1 есть предел выносливости при регулярном симметричном нагружении.

Закон движения поверхности выносливости в (2.1) задается выражением для скорости ${\mathbf{\dot {\alpha }}}$ движения ее центра. Здесь $h(\beta ) \geqslant 0$ есть материальная функция, задающая скорость движения поверхности выносливости в зависимости от значения $\beta $. Функция Хевисайда $H(\beta )$ исполняет роль переключателя: пока напряженное состояние находится внутри поверхности выносливости, она не движется. Материальная функция f определяет зависимость ${\mathbf{\dot {\alpha }}}$ от угла между направлением ${{\dot {s}}}$ движения точки, изображающей напряженное состояние, и направлением ${\mathbf{s}} - {\mathbf{\alpha }}$ напряженного состояния относительно центра поверхности. Данный угол определяется как $\cos (\angle ({\mathbf{\dot {s}}},{\mathbf{s}} - {\mathbf{\alpha }}))$ = = $({\mathbf{\dot {s}}}:({\mathbf{s}} - {\mathbf{\alpha }})){\text{/}}({\text{||}}{\mathbf{\dot {s}}}{\text{||||}}{\mathbf{s}} - {\mathbf{\alpha }}{\text{||}})$; $0 \geqslant f \geqslant 1$, причем если эти направления совпадают ($\cos (\angle )$ = = 1), то f = 1 и если они противоположно направлены ($\cos (\angle ) = - 1$), то $f = 0$. Поверхность выносливости движется в направлении напряженного состояния, поскольку ${\mathbf{\dot {\alpha }}}$ пропорционально ${\mathbf{s}} - {\mathbf{\alpha }}$.

Последнее уравнение в (2.1) есть эволюционное уравнение для параметра поврежденности. Здесь $p(\beta ) \geqslant 0$ есть материальная функция, отвечающая за скорость накопления поврежденности. По структуре это уравнение отличается от уравнения для ${\mathbf{\dot {\alpha }}}$ только отсутствием тензорнозначного множителя и присутствием множителя $1{\text{/}}{{\left( {1 - F({\mathbf{\sigma }})} \right)}^{{{{L}_{4}}}}}$, где $F({\mathbf{\sigma }}) = 1$ задает поверхность прочности, а L₄ – материальный параметр. Таким образом, при приближении к поверхности прочности скорость прироста поврежденности неограниченно возрастает. Скорость прироста поврежденности зависит от расстояния до поверхности выносливости (регулируется функцией $p(\beta )$) и направления изменения напряженного состояния (регулируется функцией f, той же, что и в выражении для ${\mathbf{\dot {\alpha }}}$).

Скорость прироста поврежденности принята независимой от текущей поврежденности, так как выражение с зависимостью $\dot {D} = G(D)F({\mathbf{\sigma }}(t))$ сводится к $\dot {D} = F({\mathbf{\sigma }}(t))$ заменой переменной D. Свобода такой перенормировки допускается несвязанностью модели (2.1) с уравнениями движения среды.

Модель (2.1) склерономна, т.е. независима от скорости изменения нагрузки ${\mathbf{\sigma }}(t)$, что достигается наличием множителя ${\text{||}}{\mathbf{\dot {s}}}{\text{||}}$ в правой части выражений для ${\mathbf{\dot {\alpha }}}$ и $\dot {D}$. Для уравнений (2.1) важна только форма пути в пространстве напряжений, а в случае одноосного нагружения – только значения максимумов и минимумов нагрузки.

По сравнению с моделью Оттосена

где $C,K$ есть материальные константы, модель (2.1) имеет следующие отличия:

1. Выражение для $\beta $ усовершенствовано, чтобы описывать одноосную и сдвиговую диаграммы Хейга произвольной формы. Модель (2.2) описывает одноосную и сдвиговую диаграммы Хейга для одноосного нагружения, а в случае сдвигового нагружения она предсказывает одинаковую амплитуду выносливости ${{\tau }_{{\text{a}}}}$ для любого среднего напряжения цикла ${{\tau }_{{\text{m}}}}$, в том числе и за пределами прочности материала.

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

3. В правой части эволюционных уравнений модели (2.2) $\dot {\beta }$ заменено на ${\text{||}}{\mathbf{\dot {s}}}{\text{||}}$. На траекториях напряжений $\beta = {\text{const}} > 0$ модель (2.2) не накапливает поврежденность, несмотря на то, что напряженное состояние все время находится за пределами поверхности выносливости.

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

5. В правой части эволюционных уравнений (2.2) присутствует простой переключатель $H(\dot {\beta })$, который отключает накопление поврежденности и движение поверхности выносливости, если напряженное состояние приближается к поверхности выносливости снаружи. Он был заменен на материальную функцию $f(\angle ({\mathbf{\dot {s}}},{\mathbf{s}} - {\mathbf{\alpha }}))$, чтобы увеличить возможность модели описывать результаты экспериментов на непропорциональное нагружение.

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

3. Идентификация модели. Материальные функции $g({{I}_{1}})$ и $r({\text{||}}{\mathbf{\alpha }}{\text{||}})$ в выражении для $\beta $ можно идентифицировать по двум диаграммам Хейга зависимостей пределов выносливости при одноосном и сдвиговом нагружениях от среднего напряжения в цикле.

На первом этапе рассматривается регулярное сдвиговое нагружение со средним значением ${{\tau }_{m}}$ и амплитудой ${{\tau }_{a}}$. Ненулевыми компонентами тензора напряжений являются ${{\sigma }_{{12}}} = {{\sigma }_{{21}}}$; те же ненулевые компоненты ${{\alpha }_{{12}}} = {{\alpha }_{{21}}}$ имеет и тензор ${\mathbf{\alpha }}$, так как ${\mathbf{\dot {\alpha }}}$ пропорционально ${\text{s}} - {\mathbf{\alpha }}$ (2.1). Такое нагружение лежит на диаграмме Хейга, если максимальное ${{\tau }_{{{\text{max}}}}} = {{\tau }_{{\text{m}}}} + {{\tau }_{{\text{a}}}}$ и минимальное ${{\tau }_{{{\text{min}}}}} = {{\tau }_{{\text{m}}}} - {{\tau }_{{\text{a}}}}$ напряжения цикла лежат на поверхности выносливости $\beta = 0$. Тогда весь цикл лежит внутри поверхности выносливости и касается ее границ крайними значениями, подставляя которые в выражение β = 0 (2.1), получим

где ${{\alpha }_{0}}$ есть компонента ${{\alpha }_{{12}}}$ тензора α. В ходе такого циклического нагружения ${{\alpha }_{0}}$ не изменяется, так как цикл лежит внутри поверхности выносливости. Исключив ${{\alpha }_{0}}$ из (3.1), получим где функция ${{\tau }_{a}}({{\tau }_{{\text{m}}}})$ задает сдвиговую диаграмму Хейга. Первое выражение в (3.2) позволяет построить предсказываемую моделью диаграмму Хейга по известной функции $r(x)$, а второе, наоборот, – функцию $r(x)$ по заданной, например, экспериментальной, сдвиговой диаграмме Хейга.

На втором этапе рассматривается регулярное одноосное нагружение со средним значением ${{\sigma }_{{\text{m}}}}$ и амплитудой ${{\sigma }_{{\text{a}}}}$. Тензор напряжений имеет единственную ненулевую компоненту ${{\sigma }_{{11}}}$, а значит тензор α имеет ненулевые диагональные элементы ${{\alpha }_{{11}}}$ и α₂₂ = = ${{\alpha }_{{33}}} = - \frac{1}{2}{{\alpha }_{{11}}}$, что следует из пропорциональности ${\mathbf{\dot {\alpha }}}$ и ${\text{s}} - {\mathbf{\alpha }}$ (2.1). Как и ранее, нагружение лежит на диаграмме Хейга если максимальное ${{\sigma }_{{{\text{max}}}}} = {{\sigma }_{{\text{m}}}} + {{\sigma }_{{\text{a}}}}$ и минимальное ${{\sigma }_{{{\text{min}}}}} = {{\sigma }_{{\text{m}}}} - {{\sigma }_{{\text{a}}}}$ напряжения цикла лежат на поверхности выносливости $\beta = 0$, подставляя которые в выражение для $\beta $ (2.1). и обозначив ${{\alpha }_{{11}}}$ как ${{\alpha }_{0}}$, получим

Исключив из (3.3) ${{\alpha }_{0}}$, получим

Если известна функция ${{\sigma }_{{\text{a}}}}({{\sigma }_{{\text{m}}}})$, задающая диаграмму Хейга для одноосного нагружения, то выражение (3.4) можно рассматривать как уравнение относительно функции $g(x)$. Это функциональное уравнение достаточно сложное, поэтому удобнее пойти другим путем. Если заданной считать функцию $g(x)$, то выражение (3.4) будет определять неявную зависимость между ${{\sigma }_{{\text{a}}}}$ и ${{\sigma }_{{\text{m}}}}$. Пусть экспериментальная диаграмма Хейга задана в формате набора из $n$ пар $\{ \sigma _{{\text{m}}}^{i},\sigma _{{\text{a}}}^{i}\} ,\;i = 1..n$. Будем искать $g(x)$ как функцию, параметризованную набором параметров ${{k}_{i}},\;i = 1..m$, достаточно гибко регулирующих ее форму. Функция g(x) должна быть монотонно возрастающей, так как повышение ${{I}_{1}}$ негативно влияет на усталостную прочность. При некотором конкретном наборе параметров ${{k}_{i}}$ выражение (3.4) сопоставляет значениям $\sigma _{{\text{m}}}^{i}$ значения $\hat {\sigma }_{{\text{a}}}^{i}$, отличающиеся от экспериментальных $\sigma _{{\text{a}}}^{i}$. Параметры ${{k}_{i}}$ подбираются так, что невязка между прогнозируемыми $\hat {\sigma }_{{\text{a}}}^{i}$ и экспериментальными $\sigma _{{\text{a}}}^{i}$ минимальна.

Модель (2.1) далее идентифицируется для алюминиевого сплава AA2024-T4. Экспериментальная диаграмма Хейга [10] (рис. 2,a) удовлетворительно описывается функциями $g({{I}_{1}}),$ $r({\text{||}}{\mathbf{\alpha }}{\text{||}})$ изображенными на рис. 3. Ввиду отсутствия сдвиговой диаграммы Хейга для данного материала была использована кривая, похожая по форме на данные по другим алюминиевым сплавам, а сдвиговой предел выносливости при симметричном нагружении был принят ${{\tau }_{{ - 1}}} = 0.55{{\sigma }_{{ - 1}}}$ [11]. Функция g искалась в виде аппроксимации, 4 параметра которой были найдены минимизацией невязки с экспериментальными данными симплекс-методом.

В качестве критерия прочности был выбран критерий Шлейхера [12]

который приводится к виду $F({\mathbf{\sigma }}) = 1$ при

Параметры a_i выражаются через пределы прочности на растяжение ${{\sigma }_{{{\text{tu}}}}}$, сжатие ${{\sigma }_{{{\text{cu}}}}} > 0$ и сдвиг ${{\tau }_{{\text{u}}}}$:

Для AA2024-T4 пределы прочности были приняты ${{\sigma }_{{{\text{tu}}}}} = 474$ МПа ${{\tau }_{{\text{u}}}} = 290$ МПа [13] и ${{\sigma }_{{{\text{cu}}}}} = 550$ МПа из оценочных соображений.

Следующую группу материальных функций предлагается искать в виде:

где ${{k}_{1}},\;{{b}_{1}},\;{{k}_{2}},\;{{L}_{2}},\;{{b}_{2}},\;{{k}_{3}},\;{{L}_{3}},\;{{b}_{3}}$ – материальные константы. Функция $w(x)$ представляет собой сглаженную кусочно линейную функцию, принимающую значение 0 при $x \leqslant {{b}_{1}}$ и имеющую наклон k₁ при $x > {{b}_{1}}$. Восемь констант (3.7) и константа ${{L}_{4}}$ (2.1), связанная с поверхностью прочности, определяются по экспериментальным кривым Велера и некоторым экспериментам на нерегулярное одноосное нагружение. Рекомендуется использовать две и более кривых Велера для различных коэффициентов асимметрии нагружения $R = {{\sigma }_{{{\text{min}}}}}{\text{/}}{{\sigma }_{{{\text{max}}}}}$. Записывается невязка между предсказаниями модели и экспериментальными данными, а затем находится набор параметров, соответствующий ее минимуму. Сделать это можно, например, симплекс-методом.

Для AA2024-T4 были использованы две экспериментальные диаграммы Велера с симметричным и отнулевым одноосными нагружениями [13] (рис. 4,a), а также данные по нерегулярному одноосному нагружению [9], представляющие собой эксперименты на усталостное разрушение при наложении периодических перегрузок большой амплитуды на регулярное нагружение малой амплитуды (рис. 4,b). Сплав AA2024-T351, исследованный в [9], отличается от AA2024-T4 лишь небольшим различием в температурной обработке и имеет очень близкие механические свойства [14]. Невязка записывалась как

Первая часть невязки относится к диаграммам Велера; здесь $N_{i}^{{exp}}$ есть экспериментальное значение числа циклов до разрушения на i-м нагружении, ${{N}_{i}}$ – предсказание модели на том же нагружении, $n$ – общее число экспериментальных точек, использованных в невязке (обозначены полыми значками на рис. 4,a). Вторая часть невязки относится к нерегулярному нагружению [9], где $q_{{10/50}}^{{exp}} = 2.14$ есть отношение экспериментальных значений циклов до разрушения при ${{n}_{{{\text{small}}}}} = 10$ и ${{n}_{{{\text{small}}}}}$ = 50 для нагружения, изображенного на рис. 4,b при $\sigma _{{\text{a}}}^{{{\text{small}}}} = 110$ МПа, $\sigma _{{{\text{min}}}}^{{{\text{big}}}} = - 414$ МПа, $\sigma _{{{\text{max}}}}^{{{\text{big}}}}$ = = 359 МПа; ${{q}_{{10/50}}}$ – отношение, предсказываемое моделью. Путем минимизации невязки симплекс-методом были найдены значения констант для AA2024-T4:

Этому набору констант соответствуют кривые Велера, изображенные на рис. 4,a, и ${{q}_{{10/50}}} = 2.07$.

Последняя материальная функция $f(\angle ({\mathbf{\dot {s}}},{\mathbf{s}} - {\mathbf{\alpha }}))$, отвечающая за накопление поврежденности в зависимости от направления изменения напряжений, разыскивается в виде

Материальная константа L₅ определяется из экспериментов на непропорциональное циклическое нагружение.

Для AA2024-T4 использованы экспериментальные данные [8] для нагружения, состоящего из комбинации циклического одноосного растяжения–сжатия и циклического сдвига, имеющих разность фаз $\theta $:

Наименьшее отклонение предсказаний модели от экспериментальных данных достигается при ${{L}_{5}} = 1.1$. Сравнение экспериментальных данных [8] и результатов расчетов по модели приведено на рис. 5,a; идентифицированная материальная функция $f(\angle )$ изображена на рис. 5,b. Плохое соответствие в случае ${{\sigma }_{{\text{a}}}} = 142$ МПа, ${{\tau }_{{\text{a}}}} = 82$ МПа, $\theta = \pi {\text{/}}2$ может быть следствием завышенного предела выносливости, который обычно не является ярко выраженным для алюминиевых сплавов.

4. Пример расчета. Модель использована для прогноза накопления поврежденности при нерегулярном одноосном нагружении элемента подвески авиадвигателя в течение полетного цикла. На рис. 6 изображена динамика накопления поврежденности на этой истории нагружения, а также движение поверхности выносливости. Расчет по модели (2.1) с набором констант (3.9) показывает, что деталь выдержит ${{N}_{{\text{f}}}} = 4998$ таких полетов, тогда как метод подсчета циклов в комбинации с формулой Одинга [15] эквивалентного отнулевого цикла $\sigma _{a}^{{{\text{eq}}}} = \left( {\begin{array}{*{20}{c}} {(0.2{{\sigma }_{{\text{m}}}} + {{\sigma }_{{\text{a}}}}){\text{/}}\sqrt 2 ,}&{{{\sigma }_{{\text{m}}}} < 0} \\ {\sqrt {{{\sigma }_{{\text{a}}}}({{\sigma }_{{\text{a}}}} + {{\sigma }_{{\text{m}}}})} {\text{/}}\sqrt 2 ,}&{{{\sigma }_{{\text{m}}}} \geqslant 0} \end{array}} \right.$ и основанный на диаграмме Велера (рис. 3,a), дает результат ${{N}_{{\text{f}}}} = 8226$. Метод подсчета циклов практически не учитывает наличие циклов небольшой амплитуды в истории, даже если они имеют высокое среднее напряжение (в истории на рис. 6, большой цикл с ${{\sigma }_{{{\text{min}}}}}$ = = –234 МПа, ${{\sigma }_{{{\text{max}}}}} = 421$ МПа дает 99.7%-й вклад в накопление поврежденности).

Одного эксперимента на нерегулярное нагружение, типа изображенного на рис. 3,b, оказывается недостаточно для полной идентификации модели. Так, набор констант

предсказывает почти идентичные кривые Велера и величину ${{q}_{{10/50}}}$, что и ранее найденный набор (3.9), но дает отличающееся время жизни на нерегулярном нагружении рис. 6: ${{N}_{{\text{f}}}} = 5607$. Недостающим экспериментом может служить нагружение, состоящее из полуциклов большой амплитуды $\sigma _{{\text{a}}}^{{{\text{big}}}}$ и циклов малой амплитуды $\sigma _{{\text{a}}}^{{{\text{small}}}}$ между ними (рис. 7). При ${{n}_{{{\text{small}}}}} = 8,$ $\sigma _{{\text{a}}}^{{{\text{small}}}} = 50$ МПа, $\sigma _{{\text{a}}}^{{{\text{big}}}} = 180$ МПа, наборы констант (3.9) и (4.1) предсказывают соответственно ${{N}_{{\text{f}}}} = 33725$ и ${{N}_{{\text{f}}}} = 50756$. Имея результаты подобного эксперимента, можно было бы выбрать наиболее подходящий набор констант или же сразу добавить эти результаты в невязку (3.8).

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