Способы подбора материальных констант для эффекта стабилизации мартенсита в микроструктурной модели сплавов с памятью формы

Сплавы с памятью формы (СПФ) являются функциональными материалами, испытывающими термоупругие мартенситные превращения (МП) при изменении температуры и/или напряжения. В ходе прямого МП высокосимметричная аустенитная фаза переходит в низкосимметричную мартенситную, при обратном превращении имеет место возврат к исходной аустенитной структуре. Если обратное МП служит результатом нагрева, то реализуется эффект памяти формы — восстановление предварительной деформации или полученной в ходе предшествующего прямого МП при охлаждении под нагрузкой, или же заданной в виде активной деформации образца в мартенситном состоянии. Если обратное превращение приходится на разгрузку, то материал ведёт себя как псевдоупругий — восстанавливается деформация, вызванная предшествующим прямым МП в процессе нагружения [1, 2] .

Рис. 1. Схематичное изображение зависимости деформации от температуры, иллюстрирующее ЭСМ после предварительной активной деформации в мартенситном состоянии; A _s ¹ и A ¹ f , A _s ² и A f ² – температура начала и окончания обратного МП в 1-м и 2-м циклах охлаждения–нагрева, M _s ² и M f ² – температура начала и окончания прямого МП во 2-м цикле

Важнейшими характеристиками СПФ являются температура начала и окончания прямого МП — M _s , M f , и обратного МП — A _s , A f . Эти параметры задают температурные интервалы фазовых переходов и играют ключевую роль при моделировании термомеханического отклика материала. СПФ находят широкое применение в промышленности и медицине. При этом их эксплуатационные свойства зависят от множества факторов, один из которых — эффект стабилизации мартенсита (ЭСМ) — обнаруживает себя как увеличение температуры обратного МП, вызванного предварительной деформацией. В повторном цикле охлаждение–нагрев в свободном состоянии данный эффект не наблюдается (Рис. 1) .

ЭСМ присущ разным СПФ: медным сплавам Cu–Zn–Al [3 –6] , Cu–Al–Be [7 –9] , Cu–Ni–Al [10 –12] , сплавам на основе железа [13] , сплавам Гейслера [14, 15] , никелиду титана [13, 16 –31] . При этом выделяют химическую (гомогенную) и механическую (гетерогенную) составляющие стабилизации [3, 5] .

Статья опубликована в открытом доступе по лицензии CC BY 4.0

Формирование представлений о химической стабилизации мартенсита связывается с исследованиями β -сплавов на основе меди, а в их систематизацию существенный вклад внесли работы С. Кустова и соавторов [3 –9] . В [5] на основе анализа калориметрических данных предложено разделять стабилизацию на две принципиально разные составляющие — химическую и механическую. Химическая стабилизация трактуется как увеличение температуры обратного превращения вследствие атомного переупорядочения в мартенсите, приводящего к снижению его химической свободной энергии G _M . Из-за этого повышается равновесная температура превращения T 0 , происходит «чистый сдвиг» температуры начала и конца обратного МП без существенного расширения температурного интервала превращения; химическая стабилизация возможна при изотермической выдержке образцов в мартенситном состоянии [32 –34] .

Причине возникновения ЭСМ в никелиде титана посвящено множество публикаций [16 –31] . В одной из первых работ [17] эффект связывается с накоплением дислокаций и вакансий, затрудняющих движение мартенситных границ. В работах [19, 21, 22] показано, что различные механизмы деформации по-разному влияют на ЭСМ, и это не позволяет дать однозначного толкования. В работе [23] ЭСМ объяснён внутренними напряжениями, возникающими при пластическом скольжении вследствие несовместности деформаций зёрен и их границ. Исходя из этого предложена реологическая модель, описывающая основные проявления эффекта. Однако последующее наблюдение ЭСМ в монокристаллах [24] поставило под сомнение решающую роль границ зёрен и пластического течения.

ЭСМ также исследовался А.А. Мовчаном и соавторами [13, 35] . В работе [13] анализируются экспериментальные данные для СПФ на основе железа и никелида титана. Отмечается, что температура начала обратного превращения в никелиде титана существенно возрастает с накоплением фазово-структурной деформации даже при нулевых напряжениях. Авторы связывают ЭСМ именно с фазово-структурной, а не с пластической деформацией, и предполагают, что при приближении предварительной деформации к предельному значению скорость повышения температуры обратного превращения должна стремиться к нулю. Для СПФ на основе железа в этой работе выявлена противоположная тенденция: температура уменьшалась при росте предварительной деформации. В продолжении этой работы, в [35] , для описания эффекта предложено в рассматриваемую термодинамическую модель к потенциалу Гиббса добавить слагаемое, связанное с фазово-структурной деформацией. Этот приём позволил описать изменение температуры начала обратного превращения при нулевых напряжениях.

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

Предлагаемая вниманию читателя работа направлена на описание ЭСМ в сплавах TiNi, поэтому механическая составляющая эффекта представляет наибольший интерес. В ней обсуждаются результаты продолжения исследований, начатых в [36 –38] ; здесь, как и в [28] , принимается предположение, что основной причиной ЭСМ является повреждение межмартенситных границ, которое затрудняет обратное превращение и исчезает после полного перехода материала в аустенит, так что при повторном цикле охлаждение–нагрев сдвиг ни температуры A _s , ни температуры A f не наблюдается. В работах [36 –38] микроструктурная модель СПФ дополнена уравнениями, позволяющими описывать ЭСМ для различных сплавов. Для их интегрирования необходимо знать не только базовые характеристики исследуемого материала, но и специфические параметры, влияющие на сдвиг температуры обратного превращения. Подбор этих параметров в [36 –38] осуществлён исключительно вручную, с большими затратами по времени. В данной работе рассматривается несколько способов формализации и ускорения подбора констант (выбор вручную, метод Нелдера–Мида, генетический алгоритм, нейросеть) и выполняется их сравнение.

• 2.    Модель

Модель СПФ, для которой производился подбор констант, подробно описана в работах [36 –39] . Основная суть, важная для понимания результатов данной статьи, заключается в том, что прямое и обратное превращения в используемой модели происходят при достижении обобщённой термодинамической силой некоего предела F ^fr , условно называемого «силой трения». Согласно гипотезе из работы [28] причиной ЭСМ является повреждение границ мартенсита, поэтому в модель вводится переменная ζ — мера повреждённости границ, изменение которой рассчитывается по формулам:

dC-b ^(Фм - ^фсг^^н(Фм - ^ф сг* ) drH(dr) m ^dZ = ^k 1 (1 - Ф crit )(1 + k з Z 2 ) ( ), ⁽¹⁾

ζ -ζ dZ = dФм H(dФм), (2)

Ф м (1 + k 3 Z 2 )

dZ = 0, (3)

Z 1 = k 2

( Ф м -Ф сг^ У НФ м - Ф сгй)~

r

I ^-Ф crit

В выражениях (1) - (4) Ф _м означает объёмную долю мартенсита, H ( • ) — функцию Хевисайда, а dФ _M и dr —

N

| Φ q - Φ M | ) характеризует q =1

2(N - 1)

приращения соответствующих величин; переменная r = г/Фм (где r = степень ориентированности мартенсита: г = 0—текстура хаотическая, г = 1 — полностью монодоменизированный мартенсит; Φq — доля q-го мартенситного варианта, N — число мартенситных вариантов. Уравнение (1) описывает изменение переменной ζ при переориентации мартенсита, уравнение (2) — при прямом МП, а уравнение (3) — при обратном превращении, формула (4) задаёт «устойчивое» значение повреждённости ζ1 при заданных ΦM и г. В формулы входят материальные константы модели k1, k2, k3, Фcrit, при этом k1 определяет интенсивность накопления повреждённости при переориентации мартенсита, а k2 — при прямом превращении; k3 характеризует скорость выхода параметра ζ на уровень насыщения, Φcrit — критическую объёмную долю мартенсита, при достижении которой в процессе прямого МП начинается рост повреждённости. На основе экспериментов из [28] в данной работе выбрано её значение Фcrit = 0.3; при помощи автоматизированных алгоритмов величина Фcrit не подбиралась.

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

fMse = F fr(1 + Z k), где F fr — базовая сила сопротивления движению неповреждённых границ, k — константа материала, которая в данной работе не подбирается ввиду отсутствия ясной физической интерпретации, взято k = 2. Далее величина FMfrSE входит в критерий обратного превращения как дополнительное ограничение для термодинамической движущей силы. Поэтому увеличение ζ повышает порог начала обратного превращения и в температурном представлении приводит к сдвигу температуры начала обратного МП As в сторону более высоких значений. Таким образом, алгоритм сдвига As реализуется не в виде отдельного эмпирического уравнения для температуры, а непосредственно через изменение пороговой силы FMfrSE, зависящей от ζ.

• 3.    Результаты подбора констант

Константы модели k 1 , k 2 , k 3 находятся для сплавов Ti 50 Ni 50 и Ti 49 Ni 51 на основе данных (Табл. 1) , полученных в работе [28] при реальных испытаниях трёх видов:

1-й вид — активная деформация в мартенситном состоянии: деформирование → разгрузка → нагрев;

2-й вид — охлаждение под нагрузкой: нагрузка в аустенитном состоянии → охлаждение до мартенсита → разгрузка → нагрев;

3-й вид — наведение мартенсита напряжениями: нагрузка в аустените при температуре T , где M s < T < A s , → разгрузка → нагрев.

В данной работе для идентификации параметров модели применяются разные способы: ручной выбор, вычисление детерминированным методом Нелдера–Мида и посредством реализации стохастического генетического алгоритма, определение нейросетевым методом на основе архитектуры Deep Sets. После

Таблица 1. Остаточная деформация после разгрузки при разном сдвиге температуры начала обратного МП в реальных экспериментах [28] для сплава Ti 50 Ni 50

Номер испытания j Вид эксперимента Остаточная деформация eu, % Сдвиг температуры ∆As,К 1 1-й 2.55 27 2 4.13 39 3 5.87 48 4 7.96 58 5 13.09 63 6 2-й 2.14 1 7 3.17 2 8 6.06 10 9 13.66 12 10 16.79 10 11 3-й 2.34 -8 12 4.04 1 13 5.87 10 14 8.15 17 15 12.46 30 описания всех способов приведены графики, построенные по установленным числовым значениям параметров, позволяющие выполнить сравнение подходов и оценить эффективность каждого из них на примерах сплавов Ti50Ni50 и Ti49Ni51.

• 3.1.    Ручной подбор

• 3.2.    Метод Нелдера–Мида

При подборе вручную учитывается следующее: так как параметр k 1 отвечает преимущественно за рост ζ при переориентации, то он значительно влияет на сдвиг температуры в эксперименте 1-го вида и совсем незначительно во 2-м; k 2 , наоборот, больше сказывается на результатах эксперимента 2-го вида и меньше — на результатах 1-го. Однако обе константы существенно воздействуют на эксперимент 3-го вида. Таким образом, основная сложность определения параметров модели ручным способом заключается в создании согласованного набора констант, который хорошо описывает любые смешанные виды экспериментов. Получены следующие значения констант для СПФ Ti ₅₀ Ni ₅₀ : k 1 = 5.6 ; k ₂ =4.0 ; k ₃ = 1.0 , и для СПФ Ti ₄₉ Ni ₅₁ : k 1 = 2.2 ; k ₂ = 4.0 ; k ₃ = 0.3 .

Метод Нелдера–Мида — оптимизационный метод, предназначенный для минимизации скалярной функции и не требующий информации о градиентах и производных [40] . Это делает его подходящим для задач, где функция имеет вид «чёрного ящика», например, при поиске материальных параметров, когда каждое значение функции требует вычислений в сторонней модели. Метод уже применялся для подбора других констант этой же модели в работе [41] .

В рамках данной работы метод Нелдера–Мида используется для минимизации невязки — f itness , между теоретическими (расчётными) и реальными (из натурных экспериментов) параметрами стабилизации мартенсита. Невязки рассчитываются по результатам n = 15 реальных и вычислительных экспериментов:

fitness =

n

1У (^AS,exp n j=1

^{∆A j} s, th

где j — номер эксперимента, n =15 — общее число экспериментов, AA S exp — сдвиг температуры A _s в j -м реальном эксперименте, а ∆A ^j _{s, th} — соответствующий расчётный сдвиг. Данные о реальных сдвигах температуры ∆A ^j _{s, exp} берутся из работы [28] . Для построения зависимостей сдвига температуры начала обратного МП ∆A _s от остаточной деформации после разгрузки e u используется по пять точек для каждого вида испытания (Табл. 1) . Расчёты производятся численно по модели с такими же e u , посчитаны теоретические сдвиги ∆A ^j _{s, th} при допуске по деформации 10%.

В общем случае алгоритм Нелдера-Мида работает с симплексом, у которого m +1 вершин в m -мерном пространстве. В рассматриваемой задаче оптимизируются три параметра: k 1 , k 2 , k 3 , поэтому симплексом является тетраэдр с вершинами x 1 , x ₂ , x ₃ , x ₄ e R ³ . Одна итерация метода состоит из шагов:

• 1)    В вершинах тетраэдра вычисляются значения функции f itness и по её возрастанию вершины упорядочиваются:

fitness ( x 1 ) ^ fitness ( x 2 ) ^ fitness ( x3) ^ fitness^ ^ )

Здесь x 1 соответствует минимальному, а x 4 — максимальному значению целевой функции в текущем симплексе.

• 2)    По значениям x 1 , x ₂ , x ₃ находится центр тяжести x c = (x 1 +x ₂ +x ₃ )/3 и строится отражённая точка x r = x c + (x c — X 4 ) .

• 3)    Если fitness ( x r ) < fitness ( x 1 ) , выполняется расширение x _e = x c +2(x r — x c ) . В следующий симплекс вместо x ₄ включается x e при условии fitness ( x e ) < fitness ( x r ) , иначе вместо x ₄ в работу вступает точка x r .

• 4)    Если fitness ( x 1 ) <  fitness ( x r ) <  fitness(x ₃ ) , отражённая точка x r в следующем симплексе заменяет точку x 4 без дополнительных преобразований.

• 5)    Если отражение не принято, осуществляется сжатие. При fitness ( x r ) >  fitness ( x 4 ) используется внутреннее сжатие x in = x c — (x c — x ₄ )/2 ; при fitness ( x 3 ) <  fitness ( x r ) < fitness ( x 4) — внешнее сжатие x out = x c + (x r — x c ) . Полученная точка вступает в работу в следующем симплексе, если в ней значение функции меньше сравниваемого значения: fitness ( x in ) < fitness ( x 4 ) или fitness ( x out ) < fitness ( x r ) .

• 6)    Если сжатие не принято, все вершины, кроме x 1 , смещаются к x 1 :

x } ^ew = x 1 + (x i — x 1 )/2,      i = 2,3,4.

Блок-схема процесса изображена на рисунке 2. На рисунке 3 описанные шаги показаны наглядно на примере симплекса с числом вершин m+1 = 3 . Процесс повторяется до выполнения одного из критериев сходимости: max^ | fitness(x i ) — fitness(x 1 ) | <г , или fitness ( x 1 ) <е, или до достижения предела по числу итераций. Здесь ε — заданный допуск сходимости.

Для k 1 , k ₂ , k ₃ опытным путём были выбраны естественные ограничения: k 1 e [0,10] , k ₂ e [0,20] , k ₃ e [0,5] . Начальный симплекс выглядел так: { k 1 ,k ₂ ,k ₃ } = { 5.0,10.0,2.5 } , { 6.5,10.0,2.5 } , { 5.0,12.5,2.5 } , { 5.0,10.0,3.2 } .

Методом Нелдера–Мида найдены значения констант: для сплава Ti 50 Ni 50 — k 1 = 5.98 , k 2 = 4.63 , k 3 = 1.21 при невязке f itness = 5.14 K; для сплаваTi 49 Ni 51 — k 1 = 1.46 , k 2 = 3.10 , k 3 = 0.00 при невязке f itness = 9.10 K.

Начало

Отражение ( x_r )

Нахождение центра тяжести симплекса без точки x m+1

Вычисление fitness , нахождение точки x_m+1

fitness(xr)

fitness(x)?

fitness(xr ) ≥ fitness(xm+1)?

Внутреннее сжатие ( x_in )

in fitness(xin) < fitness(xm+1)?

fitness(xr) < fitness(x1)?

Внешнее сжатие ( x ) out

Сокращение нет

fitness(x_out) fitness(x_r) ?

Конец

Рис. 2. Блок-схема, описывающая принцип работы метода Нелдера–Мида

Расширение ( x_e )

fitness(x_e) fitness(x_r) ?

да

Сформировать симплекс следующей итерации да

Соблюдён критерий остановки?

да

Рис. 3. Визуализация операций метода Нелдера–Мида в двумерном случае

• 3.3.    Генетический алгоритм

Генетический алгоритм относится к классу эволюционных стохастических методов оптимизации, имитирующих механизмы естественного отбора: наследование, скрещивание и мутации. В отличие от классических детерминированных методов, он обладает высокой устойчивостью к локальным минимумам и хорошо работает в условиях сложных, разрывных или шумных ландшафтов целевой функции, что особенно важно для задач, в которых каждая оценка функции связана с запуском крупного вычислительного модуля. Аналогичные генетические алгоритмы используются для подбора констант в различных моделях: ферментативного выпота, качества воды и конститутивных моделях механики [42 –44] .

Генетический алгоритм оперирует популяцией кандидатов (особей) — параметрических векторов {k 1 ,k ₂ ,k ₃ } , компоненты которых называются генами. Они итеративно улучшаются путём применения трёх основных генетических операторов:

• 1)    Отбор — выбор лучших решений.

• 2)    Скрещивание — рекомбинация параметров двух кандидатов.

• 3)    Мутация — локальная случайная модификация параметров.

На 1-м шаге алгоритма создаётся случайная популяция заданного размера. Каждый кандидат сразу оценивается по значению функции приспособленности — той же функции f itness , которая описана выше. Функция f itness вызывает 15 запусков микроструктурной модели для вычисления ∆A ^t _s ^h .

Далее проводится турнирный отбор: из трёх случайно выбранных особей побеждает та, у которой значение f itness минимально (лучшая приближённость к эксперименту). Победители турнира становятся «родителями» для скрещивания и получения особей следующей популяции. Скрещивание срабатывает с определённой вероятностью, при этом гены «детей» {k 1 ,k ₂ ,k ₃ } рассчитываются по формуле:

Ci = а • Pi + (1 -аур?, где Ci = ki — ген «ребёнка», Pi = ki — ген первого «родителя», Pi2 = ki — ген второго «родителя», а — случайная величина от 0 до 1. Если скрещивания не происходит, «ребёнок» полностью наследует гены первого «родителя»: Ci = P1.

С некоторой вероятностью гены «детей» испытывают гауссовскую мутацию — их значения изменяются в определённых пределах:

knew = ki +Ri •s •e, где kinew — мутировавший ген, ki — исходный ген, Ri — диапазон вариации гена (ограничения на ki), s — диапазон мутаций, ξi — стандартная случайная гауссовская величина со средним значением 0 и стандартным отклонением 1. Благодаря такому механизму обеспечивается локальное исследование окрестности решений.

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

В данной работе выбраны следующие параметры генетического алгоритма: размер популяции — 20; число поколений — 25; вероятность скрещивания — 0.9; вероятность мутации — 0.2; диапазон мутаций — 0.1.

На основе генетического алгоритма подобраны следующие константы: для сплава Ti ₅₀ Ni ₅₀ — k 1 = 6.68 , k ₂ =4.92 , k ₃ = 1.43 при значении fitness = 5.14 K; для сплава Ti ₄₉ Ni ₅₁ — k 1 = 1.49 , k ₂ = 3.11 , k ₃ = 0.00 при значении fitness = 9.12 K. Зависимости сдвига температуры AA s от деформации после разгрузки показаны на рисунках после описания всех способов подбора параметров.

• 3.4.    Нейросеть

• 3.4.1.    Генерация датасета. Обучающий датасет, представляемый в виде таблицы, генерируется при помощи многократного запуска модели СПФ и состоит из блоков данных, полученных в 15 вычислительных экспериментах с одним и тем же для каждого блока набором значений констант k 1 , k 2 , k 3 . Таблица 2 — один из блоков со

• 3.4.2.    Описание модели. Нейросеть решает задачу вида:

Нейросети применяются для подбора констант в различных моделях. Так, в работе Майсснера с соавторами [45] выполнен подбор параметров упругопластической модели MAT024 в программном комплексе LS-DYNA на основе полносвязной нейросети прямого распространения — Feed-Forward Artificial Neural Network (FFANN), обученной решать обратную задачу подбора параметров по диаграмме нагружения. Показано, что нейросетевое определение параметров по кривой текучести обеспечивает воспроизведение экспериментальных диаграмм сила– перемещение с точностью, сопоставимой, а в ряде случаев превосходящей точность классического итерационного метода оптимизации (LS-OPT), при существенно меньших вычислительных затратах. Каупом и соавторами [46] термодинамические параметры макроскопической одномерной модели СПФ подбираются с помощью полносвязной нейросети прямого распространения, обученной в постановке обратной задачи по циклическим диаграммам напряжение–деформация. Показано, что параметры теплоёмкости, теплообмена и скрытой теплоты, определённые нейросетью, обеспечивают высокую точность воспроизведения экспериментальных откликов, в том числе эффектов скоростной и амплитудной зависимости. В работе Ленцена и соавторов [47] выполнена идентификация параметров термомеханической модели СПФ с использованием многослойной FFANN, восстанавливающей до десяти механических и термодинамических параметров по экспериментальным циклическим кривым. Графики, построенные по данным численных расчётов со значениями параметров, полученными нейросетью, корректно воспроизводят экспериментальные диаграммы нагружения и разгрузки, при этом точность ограничена главным образом допущениями, принятыми в самой модели. В [48] разработан гибридный метод подбора параметров континуальной модели бародезии (модели конститутивной теории) для сыпучих материалов, основанный на нейросети с понижением размерности, обученной на данных численного моделирования бессеточным методом MESHFREE, предназначенным для решения задач гидро- и аэродинамики, а также механики сплошных сред. Показано, что нейросеть восстанавливает параметры модели со средней относительной ошибкой порядка 10 ^-3 и обеспечивает практически полное совпадение расчётных и эталонных напряжений при верификационных расчётах. Таким образом, нейросетевой подход хорошо показывает себя при подборе констант различных моделей и может рассматриваться как альтернатива классическим методам оптимизации.

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

Таблица 2. Остаточная деформация после разгрузки, сдвиг температуры начала и конца обратного МП в модельных экспериментах

Номер испытания j Вид эксперимента k1 k2 k3 Остаточная деформация eu, % Сдвиг температуры ∆As,К Сдвиг температуры ∆Af, К 1 1-й 0.0 8.7 0.9 2.57 -5.0 1.6 2 4.12 -5.4 2.5 3 5.87 -5.9 4.5 4 7.94 -6.6 6.9 5 13.07 -7.6 7.3 6 2-й 2.13 13.0 17.5 7 3.18 18.3 31.1 8 6.08 32.0 56.3 9 13.67 52.0 72.5 10 16.82 52.0 72.5 11 3-й 2.34 -1.3 11.5 12 4.05 0.8 54.8 13 5.86 12.0 63.8 14 8.15 48.3 67.7 15 12.46 48.3 67.7 значениями констант k1 = 0.0, k2 = 8.7, k3 =0.9. Первые пять строк относятся к экспериментам 1-го вида, вторые пять строк — к экспериментам 2-го вида, последующие пять строк — к экспериментам 3-го вида. Для каждого испытания записаны сдвиги температуры AAs, AAf и остаточная деформация eu. Всего создаётся D = 500 таких блоков с наборами констант, полученными методом Коробова [49], который выбран для более равномерного распределения постоянных в заданных диапазонах значений.

X = ( (AA 1 ,AA f ,e U ),...,(AA S ⁵ ,AA f ⁵ ,e ^ ) ) ^ (к 1 ,к 2 ,к з ).

Здесь X — входной блок из 15 численных экспериментов. Обучению подвергается нейросеть с архитектурой Deep Sets [50] . Её использование обусловлено инвариантностью к перестановкам элементов входного набора: результат работы сети не зависит от порядка, в котором перечислены эксперименты в датасете. Для рассматриваемой задачи это важно, поскольку порядок перечисления экспериментов в обучающей выборке не должен сам по себе влиять на значения идентифицируемых параметров. Входным набором в рассматриваемой задаче являются ∆A s и ∆A f — сдвиги температуры начала и окончания обратного МП, и e u — остаточная деформация после разгрузки. В отличие от полносвязных нейросетей, в которых данные всех экспериментов должны обрабатываться как единый вектор из 3 • 15 = 45 элементов с различными весами для разных экспериментов, Deep Sets обрабатывает данные каждого эксперимента ( AA^l s ,AA f ,e U ) посредством одной и той же функции ф с общими весами, которая превращает исходные данные в 64-мерный вектор признаков — h _i . Полученные наборы усредняются и подаются во внешнюю нейросеть ρ , благодаря чему каждый эксперимент вносит вклад в общее описание материала, а итоговое представление получается из совокупности данных всех экспериментов. Такой подход исключает неосмысленную зависимость от индекса эксперимента и улучшает переносимость модели на новые наборы испытаний.

Таким образом, Deep Sets сначала анализирует каждую строку отдельно, затем объединяет результаты. Это позволяет описать модель в виде:

/ 1 15      \

F(х)= р( 15£ф(21)) ’ где F — отображение, реализуемое полной нейросетевой моделью, zi — вектор из трёх значений (AAs,AAf ,eu) для i-го эксперимента в блоке, ф(-) — внутренняя нейросеть, обрабатывающая каждую строку в отдельности, р(9 — внешняя нейросеть, преобразующая усреднённые выходные признаки ф(-) в параметры k1, k2, k3.

Каждая строка состоит из трёх чисел:

Z i = (AA S ,AA f ,e U) .

Нейросеть преобразует её в вектор:

hi = ф(Zi), где φ — это последовательность линейных и нелинейных преобразований: ф(z) = a[W2(a(Wiz+bi))+b2^.

Здесь W 1 , W ₂ — матрицы весов; b 1 , b ₂ — смещения, подбирающиеся алгоритмом обучения; ст ( • ) = x/(1 + e ^x) — функция активации SiLU (Sigmoid Linear Unit) [51] .

После обработки всех строк сеть получает 15 векторов: h 1 ,h 2 ,...,h 15 , которые усредняются так:

H =15 ^hi , H e R 64 .

Финальный блок нейросети ρ принимает на вход вектор H (Рис. 4) и предсказывает параметры:

(ki ,k2 k )= p(H),    p(H ) = W5 a(W4 a (W3 H+b3)+64) + b5, где Wi — матрицы весов, bi — смещения, σ — функция активации SiLU.

Полная модель выглядит следующим образом:

1 ¹⁵

^{( k} 1 ^,k 2 ^,k 3 )= pl 15 £v⁽ z i )

i =1

Рис. 4. Структура нейросети с архитектурой Deep Sets

• 3.4.3.    Обучение. Для каждого сгенерированного блока экспериментальных характеристик:

Xbiock =((4A1 ,AAf,eU ,..., Jls5,AAf ,e^)), нейросеть делает предсказание:

^{( k} 1 ^,k 2 ^,k 3 ) = F ^{( X} block ).

Поскольку истинные параметры ( k 1 ,k ₂ ,k ₃ ) известны при генерации датасета, обучение сводится к минимизации расхождения между предсказанными и истинными значениями. В качестве целевой функции используется среднеквадратичное отклонение предсказанных значений параметров от истинных после их предварительной нормализации по ранее заданным пределам:

L =

D 3

— XX nrm

3Dl^l^\kl,d d=1ℓ=1

-

nrm kℓ,d

k ℓ,d - k ℓ ^min max min .

k _ℓ    - k _ℓ

Здесь d — номер блока датасета; I e { 1,2,3 } — номер подбираемой константы; верхний индекс nrm означает нормированное значение; величины вида k _ℓ ⁿ _, ^r _d ^m соответствуют предсказанным нейросетью значениям; k _ℓ ^min и k _ℓ ^max — нижняя и верхняя границы диапазона значений константы k ℓ ; функция L служит критерием качества приближения создаваемого полной нейросетевой моделью обратного отображения F : X _block ^ ( k 1 ,k ₂ ,k ₃ ) . Для повышения устойчивости обучения исходный датасет разделяется на обучающую и валидационную выборки. Обучающая часть требуется непосредственно для обновления параметров модели, тогда как валидационная применяется для контроля обобщающей способности сети и выбора момента остановки обучения. Нейросеть последовательно проходит все обучающие блоки, вычисляет функцию ошибки и корректирует параметры. Один полный проход по всему обучающему набору данных называется эпохой. Обучение включает сотни эпох, что позволяет модели постепенно улучшать аппроксимацию.

В архитектуре нейросети предусмотрен механизм отключения нейронов (dropout), который в процессе обучения случайным образом перекрывает выход для части нейронов, снижает склонность нейронной сети к переобучению, повышает устойчивость к шуму в данных, обеспечивает более надежное восстановление параметров (k 1 ,k ₂ ,k 3 ) .

Минимизация L осуществляется методом стохастического градиентного спуска с адаптивными моментами AdamW. Алгоритм метода хранит две величины, называемые первым и вторым моментами:

m i = ^e imi 1 ^{+ (1 —} Mg i ,    v t = ^e 2 v i 1 + ^{(1 —} AXg t ) .

Здесь t — номер шага обучения; i — номер линейного слоя (1-й и 2-й — слои сети φ , 3-й–5-й — слои сети ρ ); в 1 = 0.9 и в 2 = 0.999 — коэффициенты сглаживания; m ⁰ = v i = 0 ; gi — градиент функции L по матрицам весов или векторам смещений слоя i , вычисляемый методом обратного распространения ошибки. Формально данный процесс эквивалентен последовательному применению правила цепочки к композиции отображений, входящих в архитектуру модели, что может быть записано в терминах якобианов слоёв [52] .

С использованием моментов параметры модели обновляются по правилу:

t

W t+1 = W it - п -- nAW^,                              (5)

v vi + EAdam bt+1=bi - n     mi.                                         (6)

V v i ^{+ E} Adam

В формулах (5) и (6) приняты обозначения: W _i ^t и b _i ^t — матрица весов и вектор смещений слоя i на шаге обучения t ; i G { 1,2,3,4,5 } ; n — скорость обучения; A — коэффициент весового штрафа; m t = m t /(1 - e t ) ; v it = v t /(1 - в 2 ) ; E Adam = 10 ^-8 . Начальная скорость обучения равняется 3 • 10 ^-4 и по мере обучения модели может автоматически уменьшаться до предела в 10 ^-6 при отсутствии улучшения значения функции потерь на валидационной выборке в течение 10 последовательных эпох. Таким образом, веса нейросети последовательно адаптируются так, чтобы ошибка вычисления L уменьшалась на обучающем и валидационном поднаборах.

Для предотвращения переобучения используется ранняя остановка. К этому приходится прибегать, если значение валидационной ошибки не улучшается в течение заданного числа эпох. Процесс обучения прекращается, и сохраняется модель с наилучшим качеством аппроксимации. Финальные параметры {W i ,b i } 5=1 определяют нелинейную аппроксимацию обратного перехода от наблюдаемых характеристик блока к искомым константам модели.

• 3.4.4.    Инференс. При подборе констант по данным реального эксперимента формируется блок из 15 строк: z 1 ,...,z ₁₅ , вычисляется ф(z _i ) для всех строк, строится усреднённый вектор H , запускается финальный блок нейросети р ( Н ) и получаются параметры ( k 1 ,k ₂ ,k 3 ) — этот процесс называется инференсом. Инференс имеет дело с ансамблем моделей, обученных с различными начальными весами. Для каждого результата рассчитывается функция f itness и определяется набор констант с её наименьшим значением.

Нейросетевым методом найдены следующие константы: для сплава Ti ₅ Q Ni ₅ Q — k1 = 9.18 ; k ₂ = 4.96 ; k ₃ = 2.29 при значении fitness = 6.00 K; для сплава Ti ₄₉ Ni ₅₁ — k 1 = 1.48 ; k ₂ = 0.00 ; k ₃ = 0.00 при значении fitness = 13.60 K.

Зависимости сдвига температуры ∆A _s от деформации после разгрузки для всех методов подбора констант (ручного, метода Нелдера–Мида, генетического алгоритма и нейросети Deep Sets) показаны на рисунках 5 и 6 для двух разных СПФ.

e„, %

Рис. 5. Зависимости сдвига температуры начала обратного МП при трёх видах экспериментов для сплава Ti 50 Ni 50 : М – рассчитанные по модели с константами, подобранными одним из способов (вручную ( а ), методом Нелдера–Мида ( б ), при помощи генетического алгоритма ( в ), с использованием нейросети с архитектурой Deep Sets ( г )); Э – установленные в [28] путём натурных испытаний

Рис. 5. Продолжение

( в )

^ %

( г )

0-

( а )

-□-Вид 1 (Э) -О-Вид 2 (Э) -Л- Вид 3 (Э) -■- Вид 1 (М) -•-Вид2(М) -А-ВидЗ(М)

* 30

0-

6    3    6    9    12 15 18

<3 %

-□-Вид 1 (Э) -О- Вид 2 (9) - △ - Вид 3 (Э) -■- Вид 1 (М) -•-Вид2(М) -А- Вид 3 (М)

0    3    6    9    12 15 18

0-

( в )

^, %

—□— Вид 1 (3)

—О— Вид 2 (Э)

-i- Вид 3 (3)

-■- Вид 1 (М)

• Вид 2(М)

-А- Вид 3 (М)

0-

( г )

0    3    6    9   12 15 18

-□- Вид 1 О)

-О-Вид 2 (Э)

-Л-Вил 3 0)

-■-Вид 1 (М)

-•- Вид2(М)

-А- Вид 3 (М)

О 3    6    9    12 15 18

е„,% е„,%

Рис. 6. Зависимости сдвига температуры начала обратного МП при трёх видах экспериментов для сплава Ti 49 Ni 51 : М – рассчитанные по модели с константами, подобранными одним из способов (вручную ( а ), методом Нелдера–Мида ( б ), при помощи генетического алгоритма ( в ), с использованием нейросети с архитектурой Deep Sets ( г )) и Э – установленные в [28] путём натурных испытаний

• 4.    Сравнительный анализ

Для сплава Ti 50 Ni 50 метод Нелдера–Мида и генетический алгоритм обеспечивают как качественное, так и количественное соответствие расчётных и экспериментальных кривых (Рис. 5, Табл. 3) . Оба подхода дают значения функции f itness меньшие, чем ручной подбор, и достигают близких минимумов за сравнимое время. Нейросетевой подход также позволяет получить удовлетворительное описание эксперимента, однако уступает двум оптимизационным методам по значению функции невязки и оказывается чувствительным к начальным весам: в ансамбле из 100 обученных моделей только 60 дают реалистичный результат, тогда как остальные предсказывают константы на границах заданных диапазонов.

Таблица 3. Найденные значения параметров модели при идентификации по данным реальных экспериментов для сплава с памятью формы Ti 50 Ni 50

Способ подбора параметров модели	Параметры			Невязка
	k 1	k 2	k 3	f itness , К
Вручную	5.60	4.00	1.00	5.61
Метод Нелдера–Мида	5.98	4.63	1.21	5.14
Генетический алгоритм	6.68	4.92	1.43	5.14
Нейросеть класса Deep Sets	9.18	4.96	2.29	6.00

Отдельного пояснения требует последняя экспериментальная точка, соответствующая деформации 12.5% при 3-м виде предварительной деформации на рисунке 5. Для неё все четыре набора констант дают расчётный сдвиг температуры, заметно отличающийся от экспериментального значения. Это связано с ограничением в само´ й используемой модели: эволюционные уравнения (1) и (2) устроены так, что при росте предварительной деформации повреждённость границ мартенсита, а вместе с ней и сдвиг температуры обратного превращения выходят на насыщение. За скорость выхода к насыщению отвечает константа k ₃ , в знаменателе этих уравнений. Согласно экспериментальным данным, насыщение сдвига температуры наблюдается при 1-м и 2-м видах предварительной деформации, однако для 3-го вида последняя точка такой тенденции не показывает. Поэтому все варианты подбора, будучи в рамках одной и той же модели, воспроизводят насыщение и для 3-го вида деформирования. Дополнительной причиной расхождения может служить различие в способах определения температуры превращения: в исходной работе [28] она определяется методом касательных по графикам зависимости деформации от температуры, тогда как при автоматизированном подборе используется допуск 10% от восстанавливаемой деформации. Такой критерий берется вследствие большого количества обрабатываемых расчётных кривых и возможных локальных артефактов модели. В большинстве случаев он даёт значения, близкие к найденным методом касательных, однако не устраняет полностью влияние локальных артефактов модельного отклика; этим, в частности, объясняется небольшая немонотонность расчётных кривых для экспериментов 1-го вида на рисунке 5.

Для сплава Ti 49 Ni 51 метод Нелдера–Мида и генетический алгоритм также дают меньшие значения f itness по сравнению с ручным подбором (Табл. 4) , однако количественное соответствие эксперименту оказывается заметно хуже, чем для Ti 50 Ni 50 (Рис. 6) . Это также частично объясняется различием способов определения характеристических температур. Кроме того, как показано в работе [37] , более точное описание экспериментальных зависимостей для этого сплава возможно при использовании расширенного набора подбираемых констант k , k 1 , k 2 , k 3 . Нейросетевой подход для Ti 49 Ni 51 оказался наименее устойчивым: удовлетворительного количественного описания экспериментальных данных для 2-го и 3-го видов экспериментов получить не удалось.

Таблица 4. Найденные значения параметров модели при идентификации по данным модельных экспериментов для сплава с памятью формы Ti 49 Ni 51

Способ подбора параметров модели	Параметры			Невязка
	k 1	k 2	k 3	f itness , К
Вручную	2.20	4.00	0.30	11.09
Метод Нелдера–Мида	1.46	3.10	0.00	9.10
Генетический алгоритм	1.49	3.11	0.00	9.12
Нейросеть класса Deep Sets	1.48	0.00	0.00	13.60

С точки зрения устойчивости полученного решения важно отметить, что метод Нелдера–Мида является локальным методом оптимизации, поэтому выбор начального симплекса может влиять как на скорость сходимости, так и на найденный набор констант. При другом начальном симплексе алгоритм может прийти к другому локальному минимуму функции f itness . Нахождение глобального минимума предпочтительно, однако для рассматриваемой задачи принципиально важна не единственность набора k 1 , k 2 , k 3 , а достижение требуемой точности описания экспериментальных зависимостей. Сравнение числовых данных этого метода с результатами генетического алгоритма и нейросетевого подхода показывает, что разные наборы материальных констант могут давать близкие значения f itness и приемлемое соответствие экспериментальным данным. Поэтому локальный минимум может быть достаточным для практических задач, если обеспечивается требуемая точность расчёта. Одновременно такая неединственность решения указывает на ограниченную идентифицируемость параметров в текущей постановке.

Сравнение методов подбора является неполным без учёта вычислительных затрат. Один вызов целевой функции f itness требует выполнения серии из 15 вычислительных экспериментов. Для генетического алгоритма с используемыми параметрами (20 особей и 25 поколений) общее число вызовов целевой функции составляет 500, что соответствует 7500 запускам микроструктурной модели. При этом значение функции f itness , совпадающее с итоговым до второго знака после запятой, достигается уже после 2–3 тысяч запусков модели. Метод Нелдера– Мида при выбранном начальном симплексе демонстрирует сопоставимую скорость сходимости: стабилизация значения целевой функции наблюдается после порядка 2–3 тысяч запусков модели. Суммарные временные затраты автоматизированного подбора составляют 8–12 часов на доступных вычислительных мощностях (процессор Intel Core i7-11370H) при последовательном выполнении расчётов.

Нейросетевой подход требует предварительной генерации обучающего датасета. В данной работе она включает создание 500 блоков данных, то есть 7500 запусков микроструктурной модели, и занимает около 24 часов. При этом обучение и инференс при использовании графического процессора NVIDIA RTX 3050 Mobile выполняются существенно быстрее генерации данных и не нуждаются в дополнительных запусках микроструктурной модели. Ручной подбор, как правило, требует меньшего числа запусков модели, однако его продолжительность зависит от опыта исследователя, выбранной последовательности пробных расчётов и критериев остановки. Поэтому при ручном подборе корректнее говорить не о фиксированном числе запусков, а о сопоставимых временны´ х затратах при меньшей воспроизводимости процедуры.

Таким образом, основные вычислительные затраты в рассматриваемой задаче определяется временем работы микроструктурной модели, а не временем самого алгоритма оптимизации. Значительное ускорение может быть достигнуто за счёт распараллеливания вычислений, в частности, одновременного выполнения 15 расчётов, необходимых для одной оценки f itness . Альтернативный подход заключается в замене микроструктурной модели нейросетевой аппроксимацией, обученной на её данных; в этом случае процедура последующего подбора параметров может быть ускорена на порядки.

• 5.    Заключение

В работе рассмотрены три автоматизированных подхода к подбору материальных констант модели эффекта стабилизации мартенсита: метод Нелдера–Мида, генетический алгоритм и нейросеть с архитектурой Deep Sets. Для сплава Ti 50 Ni 50 метод Нелдера–Мида и генетический алгоритм дают наилучшее количественное соответствие экспериментальным данным и превосходят ручную калибровку по значению функции невязки. Для сплава Ti 49 Ni 51 эти методы обеспечивают в основном качественное соответствие и показывают, что участвующий в расчётах набор подобранных параметров недостаточен для столь же точного количественного описания.

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

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

Работа выполнена при финансовой поддержке Санкт-Петербургского государственного университета, проект 116636233.