CC BY
8
https://doi.org/10.15350/17270529.2022.1.2
УДК 519.688:538.91
Изотермическая фазово-полевая модель твердотельных превращений с упругими напряжениями
12 2 1
В. Г. Лебедев ' , В. А. Копытов , В. И. Ладьянов
1 Удмуртский федеральный исследовательский центр УрО РАН, Россия, 426067, Ижевск, ул. Т. Барамзиной, 34
2 Удмуртский государственный университет, Россия, 426034, Ижевск, ул. Университетская, 1
Аннотация. В статье предложена термодинамически согласованная модель фазовых переходов в твердых средах для случая чистых веществ. В качестве источника упругих напряжений рассмотрено изменение плотности сплошной среды при фазовом превращении, зависимость плотности фаз от температуры получена из аппроксимации экспериментальных данных. В таком подходе плотность вещества каждой фазы описывается своим уравнением, предполагающими, что при движении межфазной границы вещество убывающей фазы захватывается растущей фазой. Влияние изменения ориентации кристаллической структуры в данном подходе не рассматривается. Предлагаемая модель основана на методе фазового поля и описывает взаимодействие фазовых превращений с упругими напряжениями, вызванными изменением плотности вещества в фазах. Термодинамическое описание отдельных фаз основано на экспериментальных потенциалах Гиббса. Предельные случаи полученной модели соответствуют изотермической фазово-полевой модели в чистых веществах и модели термоупругих напряжений. Модель использована для численного анализа роста сферического зародыша феррита из аустенита и может являться основой для учета напряжений в процессах формирования внутренней структуры, в частности дендритов, в металлических растворах. Представлены графики профилей фазового поля и смещений деформации, полученные в результате численного расчета.
Ключевые слова: фазовое поле, микроструктура, упругие напряжения, неравновесная термодинамика, быстрые процессы затвердевания.
И Владимир Лебедев, e-mail: lvg@udsu.ru
Isothermal Phase-Field Model of Solid-State Transformations with Elastic Stresses
1'2 2 1 Vladimir G. Lebedev ' , Vitaly A. Kopytov , Vladimir I. Ladyanov
1 Udmurt Federal Research Center UB RAS (34, T. Baramzina St., Izhevsk, 426067, Russian Federation)
2 Udmurt State University (1, Universitetskaya St., Izhevsk, 426034, Russian Federation)
Summary. The article proposes the thermodynamically consistent model of phase transitions in solid media for the case of pure substances. As a source of elastic stresses, the change in the density of a continuous medium during the phase transformation is considered; the dependence of the phase densities on temperature is obtained from the approximation of the experimental data. In this approach, the density of the substance of each phase is described by its own equation assuming that when the interphase moves, the substance of the decreasing phase is captured by the growing phase. The effect of the change of the crystal structure orientation is not considered in this approach. The proposed model is based on the phase field method and describes the interaction of phase transformations with elastic stresses caused by changes in the density of the matter in phases. The thermodynamic description of individual phases is based on the experimental Gibbs potentials. The limiting cases of the obtained model correspond to the isothermal phase-field model in pure substances and the model of thermo-elastic stresses. The model is used for the numerical analysis of the growth of a spherical ferrite embryo from austenite and can be the basis for taking into account stresses in the processes of the internal structure formation in metal solutions, in particular, dendrites. The diagrams of the phase field profiles and strain displacements obtained as a result of the numerical calculation are presented.
Keywords: phase field, microstructure, elastic stresses, nonequilibrium thermodynamics, fast solidification processes.
И Vladimir Lebedev, e-mail: lvg@udsu.ru
ВВЕДЕНИЕ
Микроструктура металлов определяется формой, размерами, относительным количеством и взаимным расположением кристаллов (дендритов) отдельных фаз или их совокупности в виде зерен. От микроструктуры во многом зависят физические и механические свойства металлов. Получение нужной микроструктуры является одной из основных задач современной металлургии. Однако факторы, контролирующие формирование микроструктуры, весьма разнообразны и лежат в огромном диапазоне от атомных масштабов до макроскопических условий на производстве. Для описания сложных взаимодействий между различными факторами на мезоскопических масштабах может быть использован фазово-полевой подход [1]. Применимость метода фазового поля к описанию фазовых превращений в твердых телах была показана в работе [2]. Существенным фактором, влияющим на формирование микроструктуры чистых веществ, являются внутренние напряжения [3], возникающие при изменении плотности вещества в разных фазах и изменении ориентации кристаллической решетки. Первый фактор обычно считается более существенным, чем второй [4], поэтому изменением ориентации кристаллической решетки в первом приближении можно пренебречь. Более того, в целях простоты будем полагать деформации чисто упругими.
Целью настоящей работы является вывод самосогласованных фазово-полевых уравнений, в отличие от феноменологических подходов [3, 5], исследующих процесс фазового перехода на фоне заданных напряжений. Полученные уравнения включают в себя как описание фазовых переходов в чистых веществах на основе работы [6], так и учет изменения плотности вещества фаз на основе разделения законов сохранения вещества по фазам, аналогично работе [7]. Полученная модель проверяется на задаче роста сферического зародыша, описывающей рост а-фазы железа (феррита) из у-фазы (аустенита) при произвольном переохлаждении с помощью численного моделирования сферически-симметричной задачи.
ТЕРМОДИНАМИКА, ПЛОТНОСТЬ И УПРУГИЕ СВОЙСТВА ФАЗ
Для термодинамического описания фаз использованы потенциалы Гиббса [8] для чистого железа. Термодинамически а- и у-фазы железа описываются одним и тем же мольным потенциалом Гиббса BCC-фазы (ОЦК) [8], которая в некотором температурном интервале становится метастабильной. Устойчивой в этом интервале температур становится у-фаза (FCC - ГЦК). Разность потенциалов Гиббса фаз BCC и FCC показана на рис. 1.
500
-3000 -
600 800 1000 1200 1400 1600 1800 2000
Temperature, K
Рис. 1. Разница энергий Гиббса для ОЦК и ГЦК фаз железа без магнитного вклада [8] Fig. 1. Gibbs energy difference for BCC and FCC phases of Fe without magnetic contribution [8]
Поэтому данное исследование применимо для полного температурного диапазона, при котором железо находится в одной из отмеченных твердых фаз. Более того, в работе не будут разделяться и другие физические характеристики а- и у-фаз железа, полагая, что они относятся к одной и той же фазе ВСС. Учёт магнитного вклада в энергию Гиббса [8] приводит к дополнительному провалу на кривой рис. 1 вблизи температуры Кюри (1043 К).
Для зависимости плотности вещества фаз от температуры были использованы следующие соотношения [9]:
ра(Т) = 7952.85 -0.266• Т- 5.62• 10-5 • Т2,
рг (Т) = 8237.93 - 0.506 • Т, (1)
где размерности [Т] = К, [р] = кг/м . Графики зависимостей (1) представлены на рис. 2.
Зависимость упругих модулей фаз от температуры приведена на рис. 3. Поскольку данные по упругим модулям в основном рассматриваются в литературе без учета существующих фаз, и в то же время достаточно противоречивы [10, 11], приведенные данные (рис. 3) во многом являются результатом экстраполяции температурных данных по упругим модулям из работ [10, 11], с учетом их связи между собой и модулем Юнга Е :
Е =-935-,
' 33, + о,
Коэффициенты Пуассона в каждой из фаз выражаются через упругие модули как:
1 33 - 20
V, =--■-L.
г 2 3Д. + о
8100
^ _I_I_I_I_I_I_I_I_I_I_1_Ll_
с 500 750 1000 1250 1500 1750
Temperature, K
Рис. 2. Экстраполяция плотностей а- и у-фаз железа по данным работы [9]
Fig. 2. Extrapolation of the densities of a- and y-phases of iron according to the work [9]
900 1200 Temperature, K
Рис. 3. Экстраполяция модулей сжатия B и сдвига G для а- и у-фаз [10, 11]
Fig. 3. Extrapolation of compression modules B and shear G for a- and y-phases [10, 11]
ФОРМАЛИЗМ ФАЗОВОГО ПОЛЯ
Для описания состояний сплошной среды затвердевающей системы вводится параметр ф (фазовое поле), принимающий значения ( = 0 в фазе у и ( = 1 в фазе а. Однозначность локального фазового состояния нарушается вблизи границы раздела фаз, обладающей малой, но конечной толщиной. Объемные свойства фаз внутри диффузной границы интерполируются функциями вида:
р(() = (2(3 - 2(), g(ф) = (2(1 -()2, (2)
выбор которых обусловлен требованием устойчивости фаз [12].
Потенциалом, определяющим релаксацию неравновесной изотермической системы, является потенциал Гиббса всей системы:
а, (*, т) = \ (ое + оп )ау = \ оуау, (3)
где Ое - средняя по фазам объемная плотность равновесного потенциала Гиббса, величина Оп является неравновесной плотностью энергии Гиббса, связанной с процессом фазового перехода и упругостью среды, Ок - объемная плотность потенциала Гиббса.
Предполагая изменение плотности вещества в процессе фазового перехода, запишем среднюю по фазам плотность равновесного потенциала Гиббса в точке пространства г в момент времени ^ как интерполяцию по мольным плотностям потенциалов Гиббса О(г) (Т) каждой из фаз:
Ое = - X & ()Р (г, *)О (Т) + Wg(О), (4)
№ '=а,у
где Ж - высота потенциального барьера между состояниями а и у, р(г, ¿) - плотность вещества 1-й фазы в момент времени I в пространственной точке г, л - молярная масса (М = г/моль). Интерполяционные функции &(() определены как [12]
йа (() = Р((), & (() = 1 - Р(() = Р(1 -(). (5)
Объемная неравновесная плотность энергии:
Оп = 1 ((2 +&(У()2 )+ее1, (6)
учитывает вклад неравновесных эффектов [13], связанных с релаксацией скорости изменения
фазового поля ( = д(/ &, наличия межфазной границы <г(У()2 (V - оператор градиента)
[14], и плотности упругой энергии ее1. В выражении (6) у и о - кинетические коэффициенты,
которые будут определены далее.
Понимая давление в упругой среде как результат всестороннего сжатия [4]:
Р = -ДУ-и, (7)
где В0 - объемный модуль сжатия, запишем плотность упругой энергии с учетом работы давления как
о __1
ее, = В1 (V»? + Оо X Ц2 - - Хрл^ВУ'»
9, , охи; —хрлрвУ-», (8)
2 №
где О0 - упругий модуль сдвига, Ог - мольные объемы фаз при температуре Т [15]. Вектор и в выражениях (7 - 8) обозначает вектор смещений (сдвига) упругой среды, (V- и) - дивергенция вектора смещений, а тензор
и,, = -' 2
( ди. ди; ^ —'- + —1
Кдх; дх' у
-и
3 1
1
является бесследовой (ии = 0) частью (иначе - девиатором) тензора деформации, - символ Кронекера. Коэффициенты В0 и О0 являются, соответственно, интерполяцией упругих модулей объемного сжатия ( В ) и сдвига ( О ) отдельных фаз:
30 = £6 ф) в01 , О0 = £6 (фО .
■=а,у ■=а,у
УРАВНЕНИЯ МАССОПЕРЕНОСА
Рассмотрим законы сохранения в сплошной среде со смещением и(г, Х). Объемные плотности физических величин можно разделить на два типа: пропорциональные количеству вещества, которые увлекаются потоком, далее будем называть их "адвективными" (увлекаемыми), и не зависящие от количества вещества (т. е. "неадвективные"). В этом случае объемную плотность энергии Гиббса О (3) разделим на адвективный вклад
а(а) = - , (9)
м
где
о = О,-П,В0, (V-и), (10)
Неадвективный вклад равен:
= 1 [уф-+^ф)2)+Wg+В (V• и)2+О0 (11)
Полная производная по времени от адвективной величины Аа (г, Х) в системе отсчета Лагранжа, связанной со скоростью переноса V, будет определяться соотношением
йА(а) = А (а) + (^)А(а).
йХ
где как и ранее А(а) = дА(а) / дХ. Полная производная по времени от неадвективной величины А(и) будет определяться только ее частной производной.
Под скоростью переноса будем понимать величину V = и, пренебрегая, в силу малости смещений деформации, конвективным вкладом. Тогда дифференциальный закон сохранения для объемной плотности р в системе отсчета Лагранжа запишется в виде:
р = -р^. (12)
После дифференцирования средней по фазам плотности вещества
находим:
р = £6Р
2
где 66 = д6 (ф)/ дф. Подставляя (14) в закон сохранения вещества (12), получаем:
=-Тбб'рФ.
(13)
(14)
(15)
1 2 1 Последнее показывает, что изменение плотности вещества определяется как движением границы раздела фаз ( ~ ф ), так и потоками ( ~ 6iРV ), возникающими благодаря
деформации упругой среды.
Поскольку уравнение (12) является законом сохранения для средних величин, оно выполняется локально для каждой точки упругой среды. Но уравнение (15) не определяет по отдельности изменение плотности каждой из фаз. С его помощью можно учесть только полное изменение плотности во всех фазах. Следовательно, необходимо переопределить
уравнение (15), записав подобные уравнения в каждой из фаз так, чтобы в сумме они воспроизводили уравнение (15) в каждой фазе отдельно [7].
Для случая сосуществования двух фаз, с учетом (4), имеем:
X £'гРгФ = (Ру - Ра )р'(ф)ф (16)
I
В системе отсчета равномерно движущейся вдоль оси х со скоростью V межфазной границы, профиль фазового поля имеет вид ф(¿) = ф(х - vt). Тогда последнее выражение сводится к выражению
X & рф = -у>(Рг - Ра )р'(ф)ф'(г). (17)
I
Правая часть выражения (17) соответствует потоку вещества через фронт благодаря разнице плотностей фаз. Теперь заметим, что такой поток вещества, при движении фронта в направлении фазы у, может изменить только плотность ра. Это происходит потому, что
фронт захватывает вещество с плотностью ру, которое после прохождения фронта должно
стать частью фазы, описываемой полем фа. Наоборот, если фронт движется по направлению
к фазе а, то будет меняться плотность ру. Соответственно, перераспределение вещества на
фронте будет определяться движением фронта, что можно учесть, записав законы сохранения вещества для каждой из фаз в виде:
= -Арр'(фф - р&У -V, (18)
где Ар = р-р, ва =в(ф), ву =0(-ф), 0(х)является функцией Хэвисайда:
Г- х > 0
в( x) = \ (19)
[о, x < 0 v 7
Благодаря тому факту, что 0(ф) + 9(-ф) = 1,
ф = (в(ф) + в(-ф))ф, (20) сумма двух уравнений (18) воспроизводит закон сохранения вещества (12).
УРАВНЕНИЯ РЕЛАКСАЦИИ
В качестве управляющего функционала изотермической системы выберем энергию Гиббса всей системы (3), объемная плотность которой Ог разделена на сумму адвективных
( О(а)) (9) и неадвективных ( О (п) ) (11) вкладов по плотности вещества.
Отметим, что дифференцирование адвективных слагаемых под знаком интеграла должно учитывать изменение меры интегрирования [ 16]
^ = \{&О(а> + О(а\ - V + -О(п> \у. (21)
& ^ & ж у
Тогда, дифференцируя выражение (3) по времени в системе отсчета Лагранжа, связанной с потоком вещества, и требуя монотонного возрастания энтропии, получаем
г! ^^ "~(а)л_ гг(п) +
= п- хр£о~а+о
& ^Ц, „ ч
^ ' . (22) жу < о
+- Хкр, + р1 + р£У- V О
ц , у
После использования закона сохранения (18) и интегрирования полученного выражения по частям (последнее позволяет вынести скорости из-под пространственных производных), производная энергии Гиббса системы по времени (22) может быть представлена в виде:
dGs r(.SGs . SG V ^
—s = II p—^ + u—^ dV < 0,
dt J ^ Sp Su J
где
SG
' = -aV2p + ( + AFp'(p) + Wg '(p), (23)
Sp
^ = V,(BoV-u)-2V, (GUj )+V,
SU:
1 ZpABo; «j M j
(24)
Где обозначению AF в выражении (23) соответствует комбинация
\С
AF =
paFa-prFr-Ap^j^FI (25)
с компонентами Г в каждой фазе, равными:
Г = О' - В Д и) (26)
Пользуясь произволом выбора, существующим в неравновесной термодинамике [17], относительно потоков (для сохраняющихся параметров порядка) и скоростей (для несохраняющихся параметров порядка), определим скорости как
• .. О - „ О
5р И
Выбор мобильностей Мр > 0,Мм > 0 гарантирует выполнение условия монотонного
убывания энергии Гиббса (22) и определяет совместно с соотношениями (18) релаксационную динамику фазовых превращений системы в чистых веществах с учетом внутренних упругих напряжений, возникающих за счет различной плотности фаз.
Определяя времена релаксации как т = М у, запишем уравнение для фазового поля в
виде:
ТрР + р = Мр [а¥2р - АГр ' (р) - Wg' (р)\ (27)
Поскольку скорость звука в конденсированной среде существенно превышает характерные скорости фазовых превращений, считаем, что не только Ми >>М , но и
Ми . Поэтому, в силу и /Ми ^ 0, для деформации упругой среды вместо динамического уравнения имеем стационарное уравнение вида:
V (в0 V • u)+ 2(V • G0 V)u + (rot(G0 rot u)) = V
1 Zp&b я j p j
(28)
где
~ 2 Bo = Bo Go-
Физически, условие U /Mu ^ o, означает, что на характерных масштабах времени формирования упругих напряжений t' = Mut, скорости деформации очень малы. Поделив уравнение (18) почленно на множитель Mu и переходя к новой переменной времени t', можно пренебречь вкладом скорости, что дает:
Ztpt =-App' (р(. (29)
Из последнего уравнения следует, что ру = const при фазовом переходе r ^ о, а ра меняется по мере фазового перехода, согласно
LP о. = -ра - Pr У(Р)р. (30)
ПРЕДЕЛЬНЫЕ СЛУЧАИ
Полная локально-неравновесная модель твердотельных фазовых переходов в однокомпонентных системах описывается системой уравнений (27), (28), (29). В качестве проверки модели рассмотрим предельные ситуации, известные для фазовых переходов в однокомпонентных системах без напряжений и для термоупругих напряжений в одной фазе.
Пренебрегая изменением объема при фазовом переходе (Др = 0 и V- u = 0), находим, что движущие силы (26), (27) сводятся к разности потенциалов Гиббса, что соответствует локально-неравновесной изотермической модели [6].
В случае неоднородно нагретой однофазной системы, считая температуру внешним заданным параметром, медленно меняющегося по пространству вблизи некоторого среднего значения T0, правая часть в уравнении (28), с учетом теплового расширения, сведется к градиенту выражения:
Ро-(1 + c(T - To ))■ ад>/р, где B0, р, fi0 - (объемный модуль упругости, плотность вещества, мольный объем) при температуре T0 в данной фазе, р0О0 / ц = 1. Параметр а - коэффициент теплового расширения. В результате источником напряжений в уравнении (28) становится величина B(Tq )cVT , что после сокращения коэффициентов приводит к известным выражениям [18].
РАВНОВЕСИЕ СФЕРИЧЕСКОГО ЗАРОДЫША а-ФАЗЫ В МАТРИЦЕ у-ФАЗЫ
Рассмотрим предельный случай равновесия двух фаз, считая, что зародыш а-фазы сферической формы радиуса R находится в равновесии с окружающей его фазой у и температура во всем пространстве одинакова и равна T . Плотности вещества при данной
температуре обозначим как р° и р°. Соответствующие уравнения для плотностей вырождаются в тождества и оба множителя да = 0. В качестве независимых функций остается профиль фазового поля q0(r) и u(r) - радиальное смещение, которое является единственной ненулевой компонентой поля смещений при условии сферической симметрии.
Расписывая слагаемое (V ■ G0 V)u как сумму производной по направлению и лапласиана:
(V-GoV)u =(VGo-V)u + GoV2u,
c учетом тождества
V2u = V(V- u)- rot (rot u),
и равенства нулю ротора смещений rot u = 0 из-за их радиальности поля, находим, что уравнения (27) и (28) сведутся к уравнениям:
Ч dr2 (- W^g = 0,
r dr ^ dr) dq dq ^^
d (B V-u - B )-4 UdG0 = 0, dr r dr
где учтено, что p0 Q° / ß = 1 и
B = ЦB0l + 4G0l).
i=aj V 3 у
Первое из уравнений (31) определяет термодинамическое равновесие в системе, а второе - механическое. С учетом ограниченности смещений при r ^ 0
lim (r 2u )= 0,
и сферической симметрии запишем граничное условие в виде u(0) = 0.
Отсутствие внешних сил при r ^ да приводит к условию:
lim
L \du 2u 11 — V )--hV -
V Y> л 7
7 dr 7 r
= 0,
где v7 - коэффициент Пуассона у-фазы, приводит к условию u'(r)|г_ко = 0. С учетом граничных условий на фазовое поле:
lim f dp\ = 0 lim f dp\ = 0 (32)
^ dr j dr j
получаем краевую задачу для первого из уравнений (31).
Расписывая условия термодинамического равновесия, находим:
о0 + W^g-AG0 & +(V-и)АЛ0± = 0, (33)
r dp dp dp
где p' = dp/dr ,
AG0 = - (paGaT )-prGr(T0)), AD0 = 1 (00^0)-B^T ))
,0 =
Р
Умножая уравнение (33) на р\ с учетом полных производных, запишем
1 р'2 + ^р'2 -Wdg-АО0& + и)АО0^ = 0. (34)
2 ёт т ёт ёт ёт
Чтобы выделить характерный пространственный масштаб явления, сделаем замену
т
определяя масштаб Ь таким образом, чтобы коэффициент при производной dg / стал равен 4 (выбор множителя связан с наиболее простым видом стационарного решения). Отсюда находим:
8 = . 1^. (35)
\W
Соответственно, вместо уравнения (33) имеем:
ёр'2 + 4р'2 - 4^ - АО0 & + (V • и)А00 ^ = 0, (36)
где теперь р' = ёр/, АО0 = АО0 / W, АО0 = АО0 / W .
После интегрирования в бесконечных пределах, из-за граничных условий (32) и начального состояния среды, подстановка дает нулевой вклад:
(р2 - 4g(р))( = 0. (37)
В результате приходим к соотношению:
Ф = |[р^ - АО0 + АО0(V •«)]^ = 0, (38)
из которого видно, что термодинамическое равновесие (у - а) фаз при наличии зародышей может существовать при любой температуре, ниже температуры равновесного фазового перехода, за счет компенсации поверхностной энергии, с учетом кривизны и упругой энергии. Из того, что производные ёр /и р'2 отличны от нуля только вблизи границы фаз, легко видеть, что условие (38) определяет положение этой границы.
Чтобы сохранить решение в виде гиперболического тангенса, аналогичное одномерной задаче затвердевания в работе [12], потребуем, чтобы в уравнении (36) независимо обращались в ноль выражение
d(q'2 -4g(q))= 0, (39)
dq
и оставшаяся часть уравнения (36) выполнялась не глобально, как (38), а локально:
q + (до0 (V ■ u)-AG0 = 0. (40)
q dp
Тогда имеем
q'2 - 4g(q) = const, (41)
из всех решений которого граничным условиям удовлетворяет только функция
q0 = 1 (1 - tanh(q-4 )). (42)
Используя данное решение, находим, что вклад поверхностной энергии
ад ад
<r|(Vq)2 r 2dr = х\r 2S(r - r0 )dr, (43)
0 0
где 5(r - r0) - дельта-функция Дирака, x - коэффициент поверхностного натяжения. После замены r ^Sq находим значение коэффициента о как
( = 3Sx. (44)
Для границы у-а, используя % = 1 J/m2, S = 109m [19], находим из (35) и (44)
( = 3дх= 3 -10-9J/m, W = 2( = 6 -109 J/m3 (45)
S
Вычисляя производные, находим
dq p(q0)= - 4cosh43(q-4)="3(q0 )2, поэтому от уравнения (40) остается соотношение
AO0 (V- u)-AG0 = —, (46)
3q
которое можно рассматривать как дополнительное условие ко второму из уравнений (31). Расписывая дивергенцию (V - u) в сферической системе координат, находим
\ du 2
dq q
d ¡^ \ d2u 2du 2 -(V- u) =-7 +---7u.
dqK ' dq2 qdq q2
Тогда второе из уравнений (31) можно переписать как:
2 i \ d u du 3 „ AG dp 1 +2- F
■ +
dq1 dqyq в dq
где
=fu (q), (47)
0 , /I £-1 Л \П f <2Л/~0 , Л е-Л
fu (q)=■
3AG0 + 4q~l f 1 AG dp
3AO0
3AGu + q dp_
q ' ~ в dqy 0 3An J~B~dq
+ 2-
+
AB0 - ABj
(48)
РОСТ СФЕРИЧЕСКОГО ЗАРОДЫША а-ФАЗЫ В МАТРИЦЕ у-ФАЗЫ
Рассмотрим рост в железе сферически симметричного зародыша а-фазы в окружающей матрице у-фазы. Из-за наличия сферической симметрии можно считать фазовое поле ф и деформацию и функциями радиуса и времени, а уравнения (27) и (28) записать в сферических координатах.
Уравнения для фазового поля примет вид:
zjp + ф = Мо
1 Я Я * - Я r2 Яр- AFp' (р) - Wg '(p) r or or
(49)
с начальными условиями, соответствующими стационарному решению, вида
р (r,0) = 0:
p(r,0) = 1 (1 - tanh((r - r0)/^)),
где r - начальное положение диффузной границы, ё - ширина этой границы. Граничные условия (32) переписаны с учетом конечности интервала для численного расчета в виде: (Яр/Яг)г=0 = 0 и (Яр/Яг)r=z = 0, где z0 - характерный размер зерна.
Численное исследование уравнения (29) рассматривалось методом установления, поэтому оставлена производная по итерационному времени, но учтено, что благодаря радикальности, rot u = 0. В результате уравнение для деформаций примет вид:
Ur = Mu
V
4
B0 +- G0
1 Я 2 Я ~ 2 Я
Т—r — ur -B0— иг--
r Яг Яг r Яг
1 ^P^B j Q j p j
(50)
с начальным условием и(г,0) = 0 и граничными и(0, г) = 0, и'(^, г) = —и / ^. Выбор второго граничного условия обусловлен необходимостью убывания напряжения на бесконечности.
ЧИСЛЕННОЕ РЕШЕНИЕ
Поскольку уравнение (50) имеет вид параболического уравнения, для его решения можно выбрать неявную противопоточную схему [20]:
u"+1 - и"
At
V
4
B0 +- G0
^un++\ - 2и
Ar
n+1 . „ n+1 „ n+1 „ n+1 Л
n +Un-1 Un - uj-1 2 + 2
rAr
2
- Bo- ujr -
r r
1 V ' » о Яр
Pj j j Яг
(51)
где п - сетка по времени, ] - сетка по координате, ДХ - шаг по временной сетке, Дг - шаг по сетке координат. Эта схема является, безусловно, монотонной и устойчивой по Фридрихсу, поэтому ограничения на шаги по сетке накладывает только уравнение (49). Для его решения может быть использована градиентно-устойчивая схема, гарантирующая невозрастание полной энергии Гиббса [21]:
Аг2 / 2
р( ^ +
At + 2т
(2V>-(W + 6AF )р)
(n+1)
= р(й)
At2
At + 2т
(v^ - 3(W + 2AFр + 2р3 )и)
+
2 At 2т
At + 2т
P
(n)
(52)
где
2
p(n+i) =-р(n) +2 (р(п+1) -р(п))
Данная схема не гарантирует устойчивость на малых отклонениях, но, как было показано в [21], для устойчивого решения можно принять Аг « Аг2. Проверка устойчивости выполнялась эмпирически.
Для численного решения необходимо перейти к безразмерным переменным не только для времени и радиуса (^ ^ Т', г ^ г0£ ), но и для энергии Гиббса, модулей сжатия и
сдвига, значений плотности. Для простоты поделим их на значения этих величин для а-фазы при 500 К.
Применяя к схемам (51) и (52) метод прогонки при температуре 500 К со значениями констант, приведёнными в таблице, были получены графики для фазового поля (рис. 4) и деформации (рис. 5). Синим цветом (рис. 4) обозначено положение диффузной границы в начальный момент времени и в стационарном положении, когда движущая сила АР = 0 (25). На рис. 5 приведен профиль деформации в стационарном положении. Из графика видно, что зародыш а-фазы растянут, а матрица у-фазы сжата, что соответствует ожидаемому результату, а данные графики качественно соответствуют результатам, полученным в [22].
л со (б в
Радиус, м / Radius, m Рис. 4. Профиль фазового поля
Fig. 4. Phase field profile
с
си
m
ш и
a
01
и
I
01
01
Радиус, м / Radius, m Рис. 5. Профиль смещений при деформации
Fig. 5. Displacement profile during deformation
На рис. 6 приведён график зависимости радиуса зерна железа от температуры. Провал на графике так же соответствует температуре Кюри. Минимальный размер при данной температуре можно объяснить тем, что деформация и за меньшее время компенсирует энергию Гиббса в движущей силе АР (25).
d
a R
су
и д
а Р
Temperature, K
Рис. 6. Зависимость размера зерна железа от температуры
Fig. 6. Dependence of iron grain size on temperature
Таблица - Значения величин в уравнениях (51), (52) при численном решении
Table - Values of quantities in equations (51), (52) for numerical solution
ЗАКЛЮЧЕНИЕ
В работе разработана теоретическая модель, записанная в виде системы согласованных уравнений в частных производных, описывающей фазовый переход двухфазной системы чистого железа. Полученная модель проверена на задаче роста сферического зародыша а-фазы в матрице у-фазы железа, полученные результату качественно соответствуют ранее известным результатам [22]. Данная модель может являться основой для учета напряжений в процессах формирования внутренней структуры, в частности, дендритов в металлических растворах.
Величина Quantities Название Name Значение Value
W Высота потенциального барьера The height of the potential barrier 0.8
z 0 Конечный размер зародыша The final size of the embryo 104 m
T Температура Temperature 10"3, c
5 Полуширина диффузной границы The half-width of the diffuse interface 109, m
T Время релаксации Relaxation time 0.001, c
с Коэффициент "диффузии" для фазового поля The "diffusion" coefficient for the phase field 0.3
rj z0 Отношение начального и конечного размера зародыша The ratio of the initial and final size of the embryo 0.01
СПИСОК ЛИТЕРАТУРЫ
1. Provatas N., Elder K. Phase-Field Methods in Materials Science and Engineering. Weinheim, Germany: Wiley-VCH Verlag GmbH & Co. KGaA, 2010. 298 p. http://dx.doi.org/10.1002/9783527631520
2. Chen L.-Q. Phase-field models for microstructure evolution // Annual Review of Materials Research, 2002, vol. 32, pp. 113-114.
https://doi.org/10.1146/annurev.matsci.32.112001.132041
3. Slutsker J., Thornton K., Roytburd A. L., Warren J. A., McFadden G. B. Phase field modeling of solidification under stress // Physical Review B, 2006, vol. 74, 014103. https://doi.org/10.1103/PhysRevB.74.014103
4. Kosevich A. M. The Crystal Lattice: Phonons, Solitons, Dislocations. Berlin: Wiley-VCH, 1999. 326 p.
5. Steinbach I., Apel M. Multi phase field model for solid state transformation with elastic strain // Physica D: Nonlinear Phenomena, 2006, vol. 217, pp. 153-160. https://doi.org/10.1016/j.physd.2006.04.001
6. Лебедев В. Г., Лебедева А. А., Галенко П. К. О
мезоскопическом описании локально-неравновесных
процессов затвердевания чистых веществ // Письма в
Журнал экспериментальной и теоретической физики. 2015. Т. 101, № 2. С. 143-147.
REFERENCES
1. Provatas N., Elder K. Phase-Field Methods in Materials Science and Engineering. Weinheim, Germany: Wiley-VCH Verlag GmbH & Co. KGaA, 2010. 298 p. http://dx.doi.org/10.1002/9783527631520
2. Chen L.-Q. Phase-field models for microstructure evolution. Annual Review of Materials Research, 2002, vol. 32, pp. 113-114.
https://doi.org/10.1146/annurev.matsci. 32.112001.132041
3. Slutsker J., Thornton K., Roytburd A. L., Warren J. A., McFadden G. B. Phase field modeling of solidification under stress. Physical Review B, 2006, vol. 74, 014103. https://doi.org/10.1103/PhysRevB.74.014103
4. Kosevich A. M. The Crystal Lattice. Berlin: Wiley-VCH, 1999. 326 p.
5. Steinbach I., Apel M. Multi phase field model for solid state transformation with elastic strain. Physica D: Nonlinear Phenomena, 2006, vol. 217, pp. 153-160. https://doi.org/10.1016/j.physd.2006.04.001
6. Lebedev V. G., Lebedeva A. A., Galenko P. K. On the Mesoscopic Description of Locally Nonequilibrium Solidification of Pure Substances. Journal of Experimental and Theoretical Physics Letters, 2015, vol. 101, iss. 2, pp. 136-140. https://doi.org/10.1134/S0021364015020101
7. Лебедев В. Г., Галенко П. К. Высокоскоростное затвердевание и плавление концентрированных растворов и конструкция параллельных Хиллерта // Расплавы. 2016. № 5. С. 422-433.
8. Dinsdale A. T. SGTE data for pure elements // Calphad, 1991, vol. 15, pp. 317-425. https://doi.org/10.1016/0364-5916(91)90030-N
9. Miettinen J. Calculation of Solidification-Related Thermophysical Properties for Steels // Metallurgical and Materials Transactions B, 1997, vol. 28B, pp. 281-297. https://doi.org/10.1007/s11663-997-0095-2
10. Ledbetter H. M., Reed R. P. Elastic Properties of Metals and Alloys, I. Iron, Nickel, and Iron-Nickel Alloys // Journal of Physical and Chemical Reference Data, 1973, vol. 2, pp. 531-617. https://doi.org/10.1063/1.3253127
11. Su D., He Y.-L, Liu J.-Q., Lu X.-G. Establishment of the Elastic Property Database of Fe-base Alloys // Proceedings of the First International Conference on Information Sciences, Machinery, Materials and Energy (ICISMME 2015), 2015, pp. 1840-1850. https://dx.doi.org/10.2991/icismme-15.2015.377
12. Kessler D. Sharp interface limits of a thermodynamically consistent solutal phase field model // Journal of Crystal Crowth, 2001, vol. 224, pp. 175-186. https://doi.org/10.1016/s0022-0248(01)00814-4
13. Galenko P., Jou D. Diffuse-interface model for rapid phase transformations in nonequilibrium systems // Physical Review E, 2005, vol. 71, 046125. https://doi.org/10.1103/PhysRevE.71.046125
14. Cahn J. W., Hilliard J. E. Free Energy of a Nonuniform System. I. Interfacial Free Energy // The Journal of Chemical Physics, 1958, vol. 28, pp. 258-267. https://doi.org/10.1063/L1744102
15. Мирзоев Ф. Х., Панченко В. Я., Шелепин Л. А. Лазерное управление процессами в твердом теле // Успехи физических наук, 1996. № 166. С. 3-32. https://doi.org/10.3367/UFNr.0166.199601a.0003
16. Дымников В. П. Избранные главы геофизической гидродинамики. М.: ИВМ РАН, 1998. 97 с.
17. Жоу Д., Касас-Баскес Х., Лебон Дж. Расширенная необратимая термодинамика. Пер. с англ. 3 издания под ред. П.К. Галенко. М.-Ижевск: НИЦ "Регулярная и хаотическая динамика"; Институт компьютерных исследований, 2006. 528 с.
18. Ландау Л. Д., Лифшиц Е. М. Теоретическая физика. В 10-ти томах. Т. VII. Теория упругости: Учебное пособие. 4-е изд. испр. и доп. М.: Наука, 1987. 248 с.
19. Loginova I., Odqvist J., Amberg G., Agren J. The phase-field approach and solute drag modeling of the transition to massive у ^ a transformation in binary Fe-C alloys //Acta Materialia, 2003, vol. 51, pp. 1327-1339. https://doi.org/10.1016/S1359-6454(02)00527-X
7. Lebedev V. G., Galenko P. K. High-rate solidification and melting of concentrated solutions and the Hillert parallel construction. Russian Metallurgy (Metally), 2016, vol. 8, pp. 785-782.
https://doi.org/10.1134/S0036029516080097
8. Dinsdale A. T. SGTE data for pure elements. Calphad, 1991, vol. 15, pp. 317-425. https://doi.org/10.1016/0364-5916(91)90030-N
9. Miettinen J. Calculation of Solidification-Related Thermophysical Properties for Steels. Metallurgical and Materials Transactions B, 1997, vol. 28B, pp. 281-297. https://doi.org/10.1007/s11663-997-0095-2
10. Ledbetter H. M., Reed R. P. Elastic Properties of Metals and Alloys, I. Iron, Nickel, and Iron-Nickel Alloys.
Journal of Physical and Chemical Reference Data, 1973, vol. 2, pp. 531-617. https://doi.org/10.1063/L3253127
11. Su D., He Y.-L, Liu J.-Q., Lu X.-G. Establishment of the Elastic Property Database of Fe-base Alloys.
Proceedings of the First International Conference on Information Sciences, Machinery, Materials and Energy (ICISMME 2015), 2015, pp. 1840-1850. https://dx.doi.org/10.2991/icismme-15.2015.377
12. Kessler D. Sharp interface limits of a thermodynamically consistent solutal phase field model.
Journal of Crystal Crowth, 2001, vol. 224, pp. 175-186. https://doi.org/10.1016/s0022-0248(01)00814-4
13. Galenko P., Jou D. Diffuse-interface model for rapid phase transformations in nonequilibrium systems.
Physical Review E, 2005, vol. 71, 046125. https://doi.org/10.1103/PhysRevE.71.046125
14. Cahn J. W., Hilliard J. E. Free Energy of a Nonuniform System. I. Interfacial Free Energy. The
Journal of Chemical Physics, 1958, vol. 28, pp. 258-267. https://doi.org/10.1063/1. 1744102
15. Mirzoev F. Kh., Panchenko V. Ya., Shelepin L. A. Laser control processes in solids. Physics-Uspekhi, 1996, vol. 39, no. 1, pp. 1-29.
https://doi.org/10.1070/PU1996v039n01ABEH000125
16. Dymnikov V. P. Izbrannye glavy geofizicheskoy gidrodinamiki [Selected Chapters of Geophysical Fluid Dynamics]. Moscow: IVM RAN Publ., 1998, 97 p.
17. Jou D., Casas-Vázquez J., Lebon G. Extended irreversible thermodynamics. Third, revised and enlarged ed. Springer. 2001; 4th ed. N. Y. Dordrecht-Heidelberg-London: Springer, 2010. 483 p. https://doi.org/10.1007/978-90-481-3074-0
18. Landau L. D., Lifshits E. M. Teoreticheskaya fizika. V 10-ti tomakh. T. VII. Teoriya uprugosti [Theoretical physics. In 10 volumes. Vol. VII. Theory of elasticity]. Ucheb. posobie. 4-e izd. ispr. i dop. Moscow: Nauka Publ., 1987. 248 p.
19. Loginova I., Odqvist J., Amberg G., Agren J. The phase-field approach and solute drag modeling of the transition to massive y ^ a transformation in binary Fe-C alloys. Acta Materialia, 2003, vol. 51, pp. 1327-1339. https://doi.org/10.1016/S1359-6454(02)00527-X
20. Самарский А. А., Николаев Е. С. Методы решения сеточных уравнений. М.: Наука, 1978. 592 c.
21. Лебедев В. Г., Обухов А. А., Бовин В. П., Ладьянов В. И. Градиентно -устойчивый алгоритм для динамики фазового поля при моделировании направленного затвердевания в Si-As // Химическая физика и мезоскопия. 2021. Т. 23, № 3, pp. 312-324. https://doi.org/10.15350/17270529.202L3.28
22. Steinbach I., Apel M. Multi phase field model for solid state transformation with elastic strain // Physica D, Nonlinear Phenomena, 2006, vol. 217, pp. 153-160. https://doi.org/10.1016/j.physd.2006.04.001
20. Samarskiy A. A., Nikolaev E. S. Metody resheniya setochnykh uravneniy [Methods for solving grid equations]. Moscow: Nauka Publ., 1978. 592 p.
21. Lebedev V. G., Obukhov A. A., Bovin V. P., Lad'yanov V. I. Gradientno-ustoychivyy algoritm dlya dinamiki fazovogo polya pri modelirovanii napravlennogo zatverdevaniya v Si-As [Gradient-stable algorithm for phase field dynamics in the simulation of directional solidification in Si-As]. Khimicheskayafizika i mezoskopiya [Chemical Physics and Mesoscopy], 2021, vol. 23, no. 3, pp. 312-324. (In Russian). https://doi.org/10.15350/17270529.202L3.28
22. Steinbach I., Apel M. Multi phase field model for solid state transformation with elastic strain. Physica D,
Nonlinear Phenomena, 2006, vol. 217, pp. 153-160. https://doi.org/10.1016/j.physd.2006.04.001
Поступила 28.01.2022; принята к опубликованию 24.02.2022 Received 28 January 2022; accepted 24 February 2022
Информация об авторах
Лебедев Владимир Геннадьевич, кандидат физико-математических наук, ведущий научный сотрудник НЦМФМ УдмФИЦ УрО РАН; заведующий кафедрой теоретической физики ИМИТИФ УдГУ, Ижевск, Российская Федерация, e-mail: lvg@udsu.ru
Копытов Виталий Александрович, студент кафедры теоретической физики ИМИТИФ УдГУ, Ижевск, Российская Федерация
Ладьянов Владимир Иванович, доктор физико-математических наук, руководитель Научного центра МФМ УдмФИЦ УрО РАН, Ижевск, Российская Федерация
Information about the authors
Vladimir G. Lebedev, Cand. Sci. (Phys.-Math.), IMITIF Leading Researcher, Scientific Center of the MPM, Udmurt Federal Research Center UB RAS; Head of the Department of Theoretical Physics Udmurt State University, Izhevsk, Russian Federation, e-mail: lvg@udsu.ru
Vitaly A. Kopytov, Student of the Department of Theoretical Physics, IMITIF Udmurt State University, Izhevsk, Russian Federation
Vladimir I. Ladyanov, Dr. Sci. (Phys.-Math.), Head of the Scientific Center of the MPM Udmurt Federal Research Center UB RAS, Izhevsk, Russian Federation
