Письма в ЖЭТФ, том 115, вып. 4, с. 256 - 261

Динамика перераспределения примеси на границах

фаз растворов: фазово-полевой подход

В.Г.Лебедев¹⁾

Удмуртский федеральный исследовательский центр Уральского отделения РАН, 426067 Ижевск, Россия

Институт математики, информационных технологий и физики,

Удмуртский государственный университет, 426034 Ижевск, Россия

Поступила в редакцию 20 декабря 2021 г.

После переработки 28 декабря 2021 г.

Принята к публикации 28 декабря 2021 г.

Предложена термодинамически согласованная фазово-полевая модель локально-неравновесного за-

твердевания двухфазных концентрированных бинарных растворов. В основе модели лежит разделение

концентрационного поля примеси на независимые поля концентраций, формально определенные во всем

доступном пространстве, для каждой из фаз. Уравнения динамики примеси следуют из разделения пол-

ного закона сохранения примеси на отдельные фазы. Уравнение фазового поля и выражения для потоков

получены из общих принципов неравновесной термодинамики. Особенностью полученных уравнений яв-

ляется захват примеси растущей фазой при движении границы раздела. Диффузионный поток разделен

по процессам и по фазам. Полученная модель аппробирована с помощью численного моделирования од-

номерной задачи направленного затвердевания раствора Si-As.

DOI: 10.31857/S1234567822040085

Исследование фазовых превращений в материа-

границе и определения движущих сил фазовых пре-

лах остается актуальным направлением до настоя-

вращений в растворах [13-15]. Однако все они обла-

щего времени [1, 2]. Теоретическое описание превра-

дают теми или иными недостатками, существующи-

щений в растворах, находящихся в неравновесных

ми до сегодняшнего дня. В данной работе развива-

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

ется идея о введении диссипативных источников, но,

примеси вблизи границ фаз [3]. Основы современ-

в отличие от работы [14], согласованных с процесса-

ного понимания превращений в концентрированных

ми тепло-массо-переноса на границе фаз. Динамика

растворах с резкой границей между фазами предло-

примеси в рамках предложенной модели аппробиру-

жены Хиллертом [4]. Для систем с диффузной гра-

ется на численном моделировании фазового перехода

ницей такая задача далека от завершения. Она яв-

в расплаве Si-As в задаче направленного затвердева-

ляется существенной как для правильного описания

ния.

процессов структурообразования [5] в концентриро-

Рассмотрим неравновесное затвердевание бинар-

ванных растворах [6], так и, в конечном счете, для

ного раствора, состоящего из растворителя и некото-

разработки компьютерных систем прогнозирования

рого количества примеси с мольной концентрацией x

микроструктуры [7]. Именно на понятии диффузной

при постоянной температуре T и постоянном давле-

границы (точнее метода фазового поля [8]) основаны

нии. Будем полагать, что в системе сосуществуют две

все последние достижения в численном моделирова-

различные термодинамически равновесные фазы α

нии фазовых переходов [9, 10]. По принципу соот-

(α = L, S), ассоциированные с жидкой (L) и твердой

ветствия, движущие силы в картине резкой границы

(S) фазами. Мольные концентрации примеси в фа-

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

зах обозначим как x^α, полагая, что каждая из функ-

фазово-полевой картины.

ций определена во всем пространстве. Пренебрегая

Фактически последней работой, традиционно не

изменением плотности вещества при фазовых пере-

заострявшей вопрос о перераспределении примеси на

ходах (ρ/µ = const), и, соответственно, внутренними

границе фаз для фазового поля [11], была EFKP-

напряжениями в системе (считая, что давление внут-

модель [12]. На сегодня, в рамках фазового поля, из-

ри фаз постоянно), для описания термодинамически

вестно несколько вариантов взаимодействия фаз на

равновесных состояний раствора используем супер-

позицию объемных плотностей потенциалов Гиббса

¹⁾e-mail: lvg@udsu.ru

фаз G^α(x^α, T):

256

Письма в ЖЭТФ том 115 вып. 3 - 4

2022

Динамика перераспределения примеси на границах фаз растворов. . .

257

∑

G_eq = p^α(ϕ)Gα(x^α, T) + Wg(ϕ),

(1)

где использованы следующие обозначения: p^α′

≡

≡ ∂p^α/∂ϕ, g^′ ≡ ∂g/∂ϕ, ϕ ≡ ∂ϕ/∂t, xα ≡ ∂xα/∂t.

Далее для определенности примем, что жидкой фа-

где

G^α = ρG^α/µ,

G^α - объемная плотность, а G^α -

зе (L) соответствует ϕ = 0, а твердой (S) - ϕ = 1.

мольная плотность потенциала Гиббса фазы α. Мо-

Разделим закон сохранения примеси на отдель-

дельные интерполяционные функции

ные фазы как:

g(ϕ) = ϕ²(1 - ϕ)²

∑

p^α x^α = -ϕ˙θ^α

x^βp^β′ - p^α∇ · J - q^α.

(4)

p^S(ϕ) = 1 - p^L(ϕ) = p(ϕ) ≡ ϕ²(3 - 2ϕ),

позволяют описывать термодинамические характе-

Если отбросить источники и диффузию, то слага-

ристики многофазных систем с помощью переменной

емые с θ^α в соотношениях (4) описывают захват

фазового поля ϕ [16], являющейся параметром по-

примеси движущимся фронтом. Действительно, при

рядка фазового превращения. Функция W g(ϕ) опре-

движении фазового поля, вся примесь из убывающей

деляет потенциальный барьер между равновесными

фазы попадает в растущую фазу. Поэтому далее для

состояниями, выбор функций p^α(ϕ) связан с требо-

краткости будем называть полученную модель как

ванием монотонности на интервале [0, 1] и граничны-

модель с захватом примеси, или сокращенно - ICM

ми условиями равенства нулю производных p^α по ϕ

(impurity capture model). Диссипативные источники

на концах единичного интервала, необходимыми для

q^α связаны с произволом, возникающим из-за разде-

устойчивости термодинамических фаз [17].

ления полного закона сохранения на отдельные фазы

Критическая динамика термодинамических пере-

в ICM. При выполнении условий

менных и фазового поля определяется функциона-

лом энергии Гиббса локально-равновесной системы

θ^S + θ^L = 1,

q^L + q^S = 0,

вида (1), дополненных, в силу аддитивности, нерав-

новесными вкладами: поверхностной энергией грани-

сумма уравнений (4) по фазам приводит к общему

цы раздела между фазами, описываемой градиент-

закону сохранения примеси во всем пространстве

ным слагаемым¹² ε²(∇ϕ)² [16] и локальной неравно-

∂ξ

весностью, представленной квадратичными выраже-

+ ∇ · J = 0,

(5)

∂t

ниями [18] для скорости изменения параметра поряд-

ка для несохраняющихся величин (ϕ˙ для фазового

где J - диффузионный поток примеси, ξ - средняя

поля) и потоков для сохраняющихся величин (поток

по фазам концентрация примеси.

J для концентрации примеси). Полная энергия Гибб-

∑

са двухфазного раствора G^sol будет иметь вид:

ξ= p^αx^α.

(6)

∫[

]

G^sol(t)= G_eq +

ε²(∇ϕ)² +

βJ²+

γϕ2 dV,

(2)

Функции θ^α ≡ θ(ϕ˙^α) в уравнении (4) определяют-

ся функцией Хэвисайда

где dV - элемент объема в пространстве. Коэффици-



ент ε² характеризует поверхностную энергию грани-



x > 0,

цы раздела фаз, в то время как кинетические коэф-

θ(x) =

x < 0,

(7)



фициенты β и γ определяют влияние потоков и ско-

x ≡ 0,

ростей изменения несохраняющегося параметра по-

рядка на релаксацию неравновесной системы.

поэтому условие ϕα(θS + θL) = ϕα выполняется по

Вывод уравнений динамики основан на требова-

определению.

нии монотонного уменьшения энергии Гиббса со вре-

Подставляя уравнения (4) в производную (3),

менем. Предполагая, для простоты, бесконечность

из требования монотонного убывания функционала

системы в пространстве, находим, что дифференци-

энергии Гиббса при релаксации к равновесию опре-

рование энергии Гиббса (2) по времени приводит к

деляем диффузионные потоки и уравнение для фа-

выражению:

зового поля:

∫ [(∑

)

(∑

G^sys(t) =

G^αp^α′ + W g^′

ϕ+

τ_D

J+ J = -M_D∇

p^αµα ,

(8)

]

(3)

∑

µ^αp^α x^α -ϕ˙ε²∇²ϕ + βJ

J+ γϕ˙

ϕ dV,

(

)

τ_ϕ

ϕ+ϕ˙ =M_ϕ

ε²∇²ϕ - Wg^′ - ΔΩp^′

(9)

Письма в ЖЭТФ том 115 вып. 3 - 4

2022

258

В.Г.Лебедев

где движушие силы миграции границы равны:

Для численного моделирования полученной ICM-

∑

модели, в целях простоты, пренебрежем локально-

ΔΩ = G^S (x^S ) - G^L(x^L) - Δx θ^αµ^α(x^α),

(10)

неравновесным вкладом в уравнениях (9), (11), (12),

полагая значения τ_D = 0 и τ_ϕ = 0. В этом случае

уравнения (14) сводятся к виду:

Δx = x^S - x^L, τ_D = γM_D - время релаксации потока

примеси, M_D > 0 - коэффициент мобильности для



(∑

)



p^S x^S =-ϕ˙Δxp^′ + p^S∇

p^αD_α

∇x^α ,

B-атомов, τ_ϕ = βM_ϕ - время релаксации скорости



изменения фазового поля и M_ϕ > 0 - коэффициент

(∑

)

(

)

(15)



мобильности фазового поля. Получившееся выраже-

p^L

x^L =p^L∇

p^αD_α∇xα + ∇ MDΔµ∇p ,

ние (10) точно соответствует определению движущих

сил в работе [4].

где коэффициенты диффузии определены через

Для определения диссипативных источников за-

диффузионную мобильность:

метим, что в силу линейности, потоки есть частные

решения уравнения (8), соответствующие неоднород-

∂µ^α

D_α = M_D

ностям в правой части. Разделим полный поток на

∂x^α

диффузионный (фиковский) J_F и сегрегационный

Умножая последнее соотношение на p^α и суммируя

J_S :

по α, находим обратную связь:

J=J_F +J_S,

∑

определяя J_S как частное решение уравнения

p^αD_α

M_D =

∑

(16)

p^α(∂µ^α/∂x^α)

τ_D

J_S + J_S = -M_DΔµ∇p,

(11)

В силу того, что на границах диффузного меж-

соответствующее неоднородному источнику Δµ∇p,

фазного слоя уравнения

(15) плохо определены,

где Δµ = µ^S - µ^L. Для диффузионной части потока

удобно использовать среднюю по фазам концентра-

имеем уравнение

цию (6), уравнение для которой

∑

(

(∑

))

τ_D

J_F + J_F = -M_D p^α∇µ^α,

(12)

ξ = ∇ M_D∇ p^αµ^α(x^α)

(17)

связывающее фиковский поток J_F с неоднородно-

хорошо определено в обеих фазах.

стью химпотенциала. С учетом того, что при за-

В целях простоты рассмотрим одномерную за-

твердевании процесс вытеснения примеси (сегрега-

дачу затвердевания раствора Si-As на бесконечном

ция) может наблюдаться только в убывающей (бо-

интервале при заданном переохлаждении (при по-

лее неупорядоченной) фазе, определим диссипатив-

стоянной температуре T

= const, которая меньше

ные источники q^S , q^L в рамках ICM как

температуры ликвидуса при заданной концентрации

q^S = -q^L = -p^S∇ · J_S.

(13)

T < T_L(x^L)). Будем полагать, что направление дви-

жения фронта совпадает с положительным направ-

Такое определение диссипативных источников на

лением оси z.

границах фаз позволяет избавиться от сегрегации в

Из стационарного одномерного решения уравне-

твердой фазе и полностью замкнуть модель ICM для

ния (9) вида

бинарных растворов.

(

(z-z₀))

В итоге, уравнения диффузии (4) при затверде-

ϕ₀ =

1 - tanh

(18)

вании раствора в ICM примут вид

2δ



где δ - полуширина, z₀ - положение диффузной гра-

 p^S ˙x^S = - ϕΔxp^′ - p^S∇ · J_F,

ницы, находим соотношения [17]

(14)



p^L x^L = -p^L∇ · J_F - ∇ · J_S.

3σ

ε² = 6σδ,

W =

(19)

При условии положительности мобильностей,

уравнения (9)-(12), (14) гарантируют невозрастание

где σ - коэффициент поверхностного натяжения.

энергии Гиббса и выполнение закона сохранения

Для перехода к безразмерным переменным в

вещества для примеси.

уравнениях (9), (10), (15), (17) сделаем замену z → zδ

Письма в ЖЭТФ том 115 вып. 3 - 4

2022

Динамика перераспределения примеси на границах фаз растворов. . .

259

и t → tδ²/ν, выразив мобильность фазового поля че-

квазистационарный режим движения фронта (дви-

рез параметр ν [19]:

жение с определенной постоянной скоростью). Этот

случай отмечен как ΔT₁ на рис. 1. Во втором случае,

M_ϕ =

(20)

6σδ

В этом случае уравнение для фазового поля в од-

номерном случае будет иметь вид

(

)

ϕ=ν ϕ^′′ -

g^′ - ΔΩp^′

(21)

где ΔΩ = ΔΩδ/6σ выступает в роли безразмерной

движущей силы фазового перехода, ϕ^′′ = ∂²ϕ/∂z².

При переходе к безразмерным переменным вид

уравнений (15) и (17) изменится только за счет пе-

реопределения координат, времени и коэффициента

диффузии

D^α →

D^α = D^α/ν.

Численное моделирование уравнений (15), (17),

Рис. 1. (Цветной онлайн) Часть равновесной диаграм-

(21) выполнялось на конечном интервале длины L.

мы Si-As и переохлаждение системы ∆T2 внутри и ∆T1

Для имитации бесконечного интервала, при прибли-

вне двухфазной зоны

жении границы фронта к правой границе, выполнял-

ся сдвиг в противоположном направлении. На гра-

если переохлаждение соответствует области двух-

ницах расчетной области выставлялись однородные

фазной зоны, движение фронта (при заданной тем-

граничные условия фон Неймана. Для численного

пературе) происходит до установления концентраци-

интегрирования уравнений использован следующий

онного равновесия в растворе, случай ΔT₂, рис. 1.

алгоритм:

Для численного решения были использованы пара-

метры из табл. 1. Потенциалы Гиббса взяты из рабо-

• явно выполнялся шаг по времени для фазового

ты [21]. Характер движения фронта затвердевания

поля при заданных концентрациях;

полностью подтверждается численными расчетами.

На рисунке 2 представлены сформировавшиеся про-

• по явной схеме для (17) находились новые зна-

чения ξ во всей расчетной области;

• в интервале 0 ≤ p ≤ 1/2 решалось уравнение

для x^L, в интервале 1/2 ≤ p ≤ 1 определялось

x^S из системы (14);

• недостающая концентрация внутри диффузной

границы восстанавливалась по ξ из соотноше-

ния (6);

• в областях, где доли второй фазы малы (p(ϕ) <

< ε или p(ϕ) > 1 - ε, при ε < 10^-8), счита-

лось, что концентрация второй фазы совпадает

с концентрацией в пребладающей фазе.

Рис. 2. (Цветной онлайн) Профили средней концен-

Для верификации модели рассмотрим численное

трации при разных переохлаждениях ∆T в интервале

моделирование изотермического процесса направ-

40 K ≤ ∆T ≤ 200 K. Красная сплошная линия соот-

ленного затвердевания в расплаве Si-As. Из общей

ветствует профилю фазового поля. Черная штриховая

теории [20] затвердевания следует, что в этом случае

линия соответствует установившемуся профилю сред-

могут наблюдаться два сценария. В первом случае,

ней по фазам концентрации ξ при переохлаждении на

если переохлаждение (отклонение от температуры

величину ∆T = 41 K, зеленый пунктир - при переохла-

ликвидуса) задано вне двухфазной зоны диаграммы

ждении ∆T = 45 K, синий штрих-пунктир - при пере-

равновесия (ниже солидуса), должен наблюдаться

охлаждении ∆T = 60 K

Письма в ЖЭТФ том 115 вып. 3 - 4

2022

260

В.Г.Лебедев

фили примеси относительно положения фронта, за-

даваемого фазовым полем. Отметим, что движение

фронта происходит с полным захватом примеси, т.е.

концентрации примеси в жидкой и твердых фазах

совпадают. На рисунке 3 представлен график зави-

Рис. 4. (Цветной онлайн) Профили средней по фазам

концентрации в разные моменты времени относительно

положения фазового поля при переохлаждении ∆T =

= 200 K. Красная сплошная линия соответствует про-

филю фазового поля. Черная штриховая линия соот-

ветствует установившемуся профилю средней по фазам

концентрации ξ в момент времени t3, зеленый пунк-

тир - в момент t2, синий штрих-пунктир - в момент t1,

Рис. 3. Зависимость скорости движения фронта от пе-

такие, что t1 < t2 < t3

реохлаждения системы ∆T

Таким образом, в работе получены уравнения

симости скорости движения фронта от величины пе-

термодинамически согласованной модели затверде-

реохлаждения в растворе. Из приведенного графи-

вания бинарных концентрированных растворов. Осо-

ка видно, что при больших переохлаждениях наблю-

бенностью модели является учет захвата примеси

дается линейная зависимость V (ΔT ). Такая зависи-

растущей фазой и сегрегационного потока приме-

мость, в общем случае, не характерна для V (ΔT ),

си в неупорядоченной фазе. Численное моделирова-

являясь, в данном случае, артефактом вида потен-

ние затвердевания расплава Si-As и сравнение с из-

циалов Гиббса для Si-As. Кроме того, из графика

вестными ранее результатами показывает адекват-

V (ΔT ) видно, что при переохлаждениях T ≤ 40 K

ность модели и ее пригодность для использования

движение фронта происходит уже по второму вари-

в задачах моделирования микроструктуры материа-

анту.

лов.

Таблица 1. Параметры для расплава Si-0.1 ат. % As [19]

Величина

Значение

1. Р. А. Кончаков, А. С. Макаров, А. С. Аронин,

1.2 · 10^-5 м³/моль

Н. П. Кобелев, В. А. Хоник, Письма в ЖЭТФ 113(5),

D_L

1.5 · 10^-9 м²/с

341 (2021).

D_S

3 · 10^-13 м²/с

2. В. В. Бражкин, Т. И. Дюжева, И. П. Зибров, Письма

0.477 Дж/(м²)

в ЖЭТФ 114(8), 541 (2021).

Mϕ

8.777 м³/(Дж · с)

3. P. K. Galenko, V. Ankudinov, K. Reuther,

1.22 · 10^-8 м²/с

1.875 · 10^-9 м

Rettenmayr,

Salhoumi,

and

E. V. Kharanzhevskiy, Philos. Trans. R. Soc. A

Math. Phys. Eng. Sci. 377, 2143, 20180205 (2019).

Поведение фронта затвердевания ΔT = 200 K для

4. M. Hillert and V. Rettenmayr, Acta Mater. 51, 2803

второго варианта представлено на рис.4, где пока-

(2003).

заны профили концентрации относительно фазового

5. I. Singer-Loginova and H. M. Singer, Rep. Prog. Phys.

поля. Из рисунка 4 видно, что эти профили выравни-

71, 106501 (2008).

ваются, принимая в итоге вид ступеньки. Параллель-

6. M. Hillert, Phase Equilibia, Phase Diagrams and

но с процессом релаксации концентрации, граница

Phase Transformations: their thermodynamic basis,

фронта движется, постоянно замедляясь, до полной

Cambridge University Press, UK, Cambridge (2008).

остановки в равновесии.

7. http://web.micress.de/ from 10.11.2021.

Письма в ЖЭТФ том 115 вып. 3 - 4

2022