CC BY
24
УДК 536.421.4, 536.421.48
ЧИСЛЕННОЕ РЕШЕНИЕ УРАВНЕНИЯ АЛЛЕНА-КАНА НА ОСНОВЕ ВАРИАЦИОННОГО ПРИНЦИПА
ОБУХОВ А.А., *ЛЕБЕДЕВ В.Г., ЛАДЬЯНОВ В.И., *НОВИКОВА ТА.
Физико-технический институт УрО РАН, 426000, г. Ижевск, ул. Кирова, 132
*Удмуртский государственный университет, 426034, г. Ижевск, ул. Университетская, 1
АННОТАЦИЯ. Сформулирован вариационный принцип для уравнения Аллена-Кана. На его основе проведены численные расчеты в тестовой задаче моделирования процесса затвердевания расплава для одномерного случая.
КЛЮЧЕВЫЕ СЛОВА: фазовый переход, микроструктура, вариационный принцип, компьютерное моделирование.
ВВЕДЕНИЕ
Согласно общим представлениям вблизи температуры кристаллизации, при замораживании физической системы в неупорядоченной фазе протекают процессы релаксации замороженных степеней свободы, происходящие в течение некоторого характерного времени [1, 2]. Размеры упорядоченных областей, возникающих в исходной неупорядоченной (симметричной) фазе, могут расти с различными скоростями, из-за неоднородности состава, неоднородности поля локальных температур и т.п. Это приводит к конкуренции роста различных фаз с нарушенной симметрией, в результате чего в процессе кристаллизации реализуется неоднородная упорядоченная структура. Из-за сложности процесса, одним из перспективных методов теоретического исследования процесса затвердевания является компьютерное моделирование [3 - 5]. Оно позволяет исследовать и прогнозировать характеристики неравновесных физических процессов, происходящих в реальных условиях.
Важным шагом в исследовании является выбор модели для описания процесса. Например, квазиравновесная теория двухфазной зоны [6], рассматривает гетерогенную двухфазную зону как дисперсную бесструктурную среду [7, 8]. Неравновесная модель кристаллизации позволяет рассчитывать дендритную и зеренную структуру благодаря введению зависимости скорости кристаллизации от переохлаждения, для которой необходимо знать кинетический коэффициент роста фаз, определяемый эмпирически [9].
В этой статье рассматривается вариационный метод для численного моделирования процесса затвердевания, используя фазово-полевую модель [12]. В отличие от неравновесной, фазово-полевая модель основана на динамическом подходе, позволяющем описывать движение границы раздела фаз.
Вариационный подход к решению уравнений критической динамики, описанный в [10], дает возможность получать численные алгоритмы для моделирования процессов с любым заданным порядком аппроксимации по пространству и времени. В частности, его использование, позволило получить численный алгоритм для моделирования процесса спинодального распада, описываемое уравнением Кана-Хилларда, с высоким порядком аппроксимации по пространству и времени [11].
УРАВНЕНИЕ АЛЛЕНА-КАНА
Теоретическое описание процессов затвердевания на мезоскопическом уровне основано на модели фазового поля, которая описывается уравнением Аллена-Кана (АК) (модель движения границ зерен [12])
р=м^ар-/;(р,т)], (1)
где Т - температура; р - параметр порядка; ер - корреляционный коэффициент; Мр > 0 - кинетический коэффициент мобильности для несохраняющегося параметра
порядка, зависящие в общем случае от параметра порядка и температуры. Обозначение /' = д//др соответствует производной по аргументу равновесной плотности свободной энергии.
В основу модели заложено предположение о существовании локального параметра порядка р , называемого фазовым полем. Параметр порядка позволяет описать диффузионную границу между твёрдой и жидкой фазами непрерывным образом с условием, что фазовое поле может изменяться достаточно резко, оставаясь непрерывной и гладко дифференцируемой функцией в интервале значений от 0 до 1.
Рассмотрим вывод уравнения АК из условия уменьшения свободной энергии. В предположении изотермичности динамика процессов релаксации определяется свободной энергией системы, которая может быть выбрана в виде суммы вкладов от ее равновесной /(р), зависящей от локального параметра порядка р , и локально-равновесной части, учитывающей пространственную неоднородность в системе [13]. Плотность свободной энергии (потенциал Гиббса) может быть восстановлена из эксперимента [14]. Коэффициенты в локально-равновесной части с характеристиками границы раздела фаз обычно полагаются константами [10,11] для уравнений Аллена-Кана (АК) [12]
Ъ [р]=1 \ 1 е2 (V р)2 + / (р)
йУ. (2)
Рассматривая процесс для системы с несохраняющимся параметром порядка в объеме У с границей £, и дифференцируя функционал ^ = ^ [ р] по времени, получаем
= \у рр (/р -е2А р )йУ + е2&рру-рй£п * 0 , (3)
где /р - первая производная равновесной плотности свободной энергии по параметру порядка. Обозначение / ' = д/(р)/др соответствует производной по аргументу равновесной
плотности свободной энергии из соотношения (1).
Чтобы исключить перенос свободной энергии через границу, описываемый вторым интегралом, для изолированной системы в качестве естественных граничных условий удобно выбрать в виде
V- р £ = 0. (4)
Тогда утверждение о монотонном убывании ^ в простейшем варианте гарантируется выбором
/р - е2А р = -мМ-р, (5)
-
где Мп > 0 - кинетический коэффициент мобильности для несохраняющегося параметра порядка, в общем случае зависящий как от параметра порядка, так и от температуры.
В рассматриваемом изотермическом случае далее будем предполагать лишь зависимость вида Мп = Мп (р ). Полученное уравнение является уравнением АК
Ф = М-(е2А р - /р) (6)
и используется для описания движения межфазной границы при фазовых переходах первого рода [12].
Теперь рассмотрим получение уравнения АК, используя вариационный принцип.
ВАРИАЦИОННЫЙ ПРИНЦИП ДЛЯ АК
В рамках локально-равновесного подхода уравнение АК является уравнением первого порядка по времени. Для построения вариационного принципа, приводящего к уравнению первого порядка по времени, рассмотрим функционал вида:
Ф = JdQJdz - + 2 - У + + LpP,w)
(7)
где введенная функция " играет роль ''потоковой'' переменной для первой производной по времени. Варьируя (7), находим
t
ЗФ = J dQ J dr\dwiw ~ф - " + Зр(у/ - 2" + l() + Sju(p-y + 2ф)]
0 (8)
+ J dQ(2"-w)3(0,
где Lw и L p - частные производные плотности L по функциям у и p соответственно.
Из экстремальности функционала (8), исключая " и считая, что при т = 0 задана функция ((0, r ) = p0 (г ), приходим к уравнениям
Ф = - L ,
W (9)
У = 4( + L^,
с дополнительным условием в момент времени т = t вида y(t )= 4p (t), играющим вид ''начального'' условия для функции у.
Таким образом, экстремальность функционала (7) приводит к уравнениям первого порядка по времени (9), описывающим эволюцию поля р(т, r) с начальным условием p = p0 (r) в момент времени т = 0, к моменту времени т = t, и эволюцию поля у(т, r) с начальным условием у = 4p (t, r ) в момент времени т = t к моменту времени т = 0 .
Чтобы завершить построение вариационного принципа для уравнения АК, конкретизируем вид пространственной части для функционала (7), выбирая ее в виде:
L = - JV(M (У )+ M (p){(J + е2Vp)+ yf( }, (10)
где е2 > 0 - const; M(p)> 0 - коэффициент мобильности, обозначение f ( = df (pp)jdp соответствует производной по аргументу равновесной плотности свободной энергии из соотношения (2), а J и { являются неизвестными функциями (потоками). Соответствующие вариации по p, у, J, { приводят к выражению:
+
~ * +
sJLdQdt = J§J (m{ - V(My))dQdt ■
+ J MSy(yJ + fp)dQdt + J MSI(J + e2V p(dQdt JSp^Mp({(J + е2Vp)+ yVJ)+ - е2V(м{)}dQdt
+ J dt J [(е2ЛЫ - J МУЗ - MJnSy]iSn,
(11)
где индекс р соответствует производной по параметру порядка р; индекс п обозначает
нормальную компоненту вектора. Вариация по потоковым переменным J и Я приводит к определению потоков
Я = М-'V (Мщ) или Л = -а2У р, (12)
а вариации по р и щ совместно с уравнениями (9) дают уравнения для параметра порядка р
р = М (е2 А р - /;) (13)
и для ''сопряженной'' функции \
\ = 4 рр-е2 А (му)+ м\/^. (14)
Естественные граничные условия появляются при занулении поверхностного вклада, определяемого интегралом:
$[(е2ЛпМ - ЗпМ\)5р - ЗпМ8ц К = 0, (15)
из которого следует, что требование произвольности изменения функции 5\ на границе
области приводит к обращению в ноль потока J на границе:
5ц\£ * 0 ^ Jn\s = 0 ^У- р\£ = 0, (16)
а произвольность вариации 5 на границе приводит к условию:
£2Л-13 = 0 ^У£ = 0. (17)
Таким образом, естественные граничные условия для выбранного функционала соответствуют граничным условиям, полученным из условия монотонного уменьшения свободной энергии.
Важно отметить особенности уравнения АК, которые влияют на выбор численного метода при решении. Уравнение является первого порядка по времени, описывает нестационарный процесс, имеет нелинейный член в виде свободной энергии. По этим причинам существующими численными методами довольно сложно получить правильное и устойчивое численное решение.
Полученное уравнение АК из вариационного принципа является более общим и имеет ряд преимуществ, в отличие от требования монотонного уменьшения свободной энергии. Вариационная постановка задачи позволяет ослабить математические ограничения, накладываемые на разыскиваемые решения, самосогласует дифференциальное уравнение с граничными условиями, а также позволяет строить устойчивые численные алгоритмы для их численной реализации [16].
ФУНКЦИЯ ПЛОТНОСТИ СВОБОДНОЙ ЭНЕРГИИ
Рассмотрим член уравнения АК, который приводит данное уравнение к нелинейному виду. Динамика диффузионной границы в рамках неравновесной термодинамики определяется плотностью свободной энергии системы / (р, Т) [17]. Функция
g(р) = р2 (1 - р)2 и функция интерполяции р(р)= р3 (бр2 -15 р +10) необходимы для определения функции плотности свободной энергии.
Для однокомпонентного расплава А функция плотности свободная энергии определяется как
/А(р,Т)= (1 -р(р))/А(Т) + р(р)/АА (т) + ^(р). (18)
Эта функция объединяет свободные энергии жидкой и твердой фаз с функцией интерполяции р(р) и добавляет барьер энергии ЖА между ними.
Для однокомпонентного расплава, определим вид свободной энергии для жидкой и твердой фазы.
/А (т ) = 0,
А (ТМА - т) (19)
/аа (Т )- /А (Т ) =
А
Подставляя (19) в (18) , получаем свободную энергию для однокомпонентного расплава в таком виде
ТА — т
fA (р, т ) = wAg (р) + ь,™— Р (р) • (20)
1 м
Подставляя в уравнение АК производную функции плотности свободной энергии по не сохраняющемуся параметру порядку р, получаем необходимое уравнение, для численного решения
р = м
С тА — т
2 Л „ ТЛ , Т Тм 1
с1 Ар — WAg'(p+ ЬА^—Р\
V Тм У
(21)
где g , р производные функций g(р), р(р) по р.
РЕЗУЛЬТАТЫ И ОБСУЖДЕНИЯ
В работе получено численное решение уравнения Аллена-Кана с модельными безразмерными коэффициентами: б2 = 0,03, WA = 10, ЬА = 500, Т = 1580, ТМ = 1600, М = 0,1. Результатом решения является зависимость фазового поля р от координаты х и времени. Фрагменты динамики движения фазового поля представлены на рис. 1 в различные моменты времени.
На рис. 2 представлено решение уравнения Аллена-Кана с параметром Ьл = 0 в различные моменты времени. На рис. 2, а изображены два решения в момент времени t = 0, численное - это заданное начальное распределение р, и аналитическое
Р(х) = 2
(
1 - tanh
2г.
Л"
(22)
соответствующее равновесному состоянию двух фаз. На рис. 2, б изображены решения в момент времени t = 5, из которого видно, что численное решение достаточно точно совпадает с аналитическим в классе линейных функции.
1 - аналитическое решение; 2 - численное решение Рис. 2. Решение уравнения Аллена-Кана с параметром Ьл = 0 при t = 0 (а) и t = 5 (б)
Моделирование проводилось на области в 20 ячеек с шагом по времени t = 0,03. На рис. 1 и 2 видны изломы, это связано с тем, что была использована линейная аппроксимация функций в ячейке р, у, 3, Я , /и . Что бы получить более гладкий вид решения, необходимо использовать аппроксимацию по пространственным переменным второго и более порядка.
ВЫВОД
В работе показано, что вариационный принцип [10] может быть использован для построения эффективных численных алгоритмов с высоким порядком аппроксимации по пространственным и временным переменным для компьютерного моделирования уравнения Аллена-Кана. В работе построен алгоритм второго порядка по времени и первого порядка по пространству, проведены численные эксперименты по моделированию процесса затвердевания расплава для одномерной задачи, качественно правильно описывающие физическое поведение системы.
Работа выполнена при частичной финансовой поддержке Программы Президиума РАН (проект 12-П-2-1044).
СПИСОК ЛИТЕРАТУРЫ
1. Bray A.J., Rutenberg A.D. Growth laws for phase ordering // Phys. Rev. E. 1994. V. 49, № 1. P. R27.
2. Bray A.J. Theory of phase ordering kinetics // Adv. Phys. 1994. V. 43, P. 357.
3. Gunton J.D., Toral R., Chakrabarti A. Numerical Studies of Phase Separation // Phys. Scripta. 1990. T. 33. P.12-19.
4. Emmerich H. Advances of and by phase-field modelling in condensed-matter physics // Advances in Physics. 2008. V. 57, № 1. P. 1-87.
5. Singer-Loginova I., Singer H.M. The phase field technique for modeling multiphase materials // Rep. Prog. Phys. 2008. V. 71. P. 106501.
6. Журавлев В.А. Затвердевание и кристаллизация сплавов с гетеропереходами (физические основы, теория, эксперименты, практика). М. - Ижевск : НИЦ «Регулярная и хаотическая динамика», Институт компьютерных исследований, 2006. 560 с.
7. Борисов В.Т., Журавлев В.А. Развитие теории двухфазной зоны и применение к промышленным проблемам // Теплообмен при кристаллизации и конденсации металлов. Новосибирск : Изд-во ИТФ СО АН СССР, 1981. С. 51-55.
8. Журавлев В.А., Борисов В.Т. Развитие теории двухфазной зоны металлических сплавов и ее приложение к проблемам литья // Кристаллизация. Теория и эксперимент / Ижевск : Изд-во УдГУ, 1987. С. 5-15.
9. Уманцев А.Р., Виноградов В.В., Борисов В.Т. Математическое моделирование роста дендритов в переохлажденном расплаве // Кристаллография. 1985. Т. 30, вып. 3. С. 455-460.
10. Новикова Т.А., Обухов А.А., Обухов А.В. и др. Вариационные принципы для уравнений критической динамики // Химическая физика и мезоскопия. 2010. Т. 12, № 3. С. 338-347.
11 Обухов А.В., Обухов А.А., Лебедев В.Г. и др. Численное моделирование спинодального распада на основе вариационного подхода // Вестник Удмуртского университета. Сер. Физика и химия. 2011. Вып. 1. С. 31-39.
12. Allen S.M., Cahn J.W. A microscopic theory for antiphase boundary motion and its application to antiphase domain coarsening // Acta Metal. 1979. V. 27. P. 1085-1095.
13. Cahn J.W., Hilliard J.E. Free energy of a nonuniform system. I. Interfacial free energy // J. Chem. Phys. 1958. V. 28. P. 258-267.
14. Cahn J.W., Hilliard J.E. Free Energy of a Nonuniform System. III. Nucleation of a Two-Component Incompressible Fluid // J. Chem. Phys. 1959. V. 31. P. 688-699.
15. Kessler D. Sharp interface limits of a thermodynamically consistent solutal phase field model // J. Cryst. Growth. 2001. V. 224. P. 175.
16. Ректорис А.П. Вариационные принципы в математической физике и технике. М. : Мир, 1985. 590 с.
17. Boettinger W. J., Warren J. A., Beckermann C., Karma A. Phase-field simulation of solidification // Annu. Rev. Mater. Res. 2002. V. 32. P. 163-194.
NUMERICAL SOLUTION OF EQUATIONS OF ALLEN-KAHN ON VARIATIONAL PRINCIPLES
Obuhov A.A., *Lebedev V.G., Ladiyanov V.I., *Novikova T.A.
Physical-Technical Institute, Ural Branch of the Russian Academy of Sciences, Izhevsk, Russia *Udmurt State University, Izhevsk, Russia
SUMMARY. Were carried out numerical simulations of solidification of the melt for one-dimensional problem, based on the variational principle for the equations of Allen-Cahn.
KEY WORDS: phase transition, microstructure, variational principles, computer simulation.
Обухов Александр Андреевич, аспирант ФТИ УрО РАН ФТИ, e-mail: Obuh87@bk.ru
Лебедев Владимир Геннадьевич, кандидат физико-математических наук, доцент кафедры теоретической физики УдГУ, тел. (3412) 91-61-30, e-mail: lvg@udsu.ru
Ладьянов Владимир Иванович, доктор физико-математических наук, директор ФТИ УрО РАН, тел. (3412) 21-89-11, 43-03-02, e-mail: las@pti.udm.ru / direct@fti.udm.ru
Новикова Татьяна Алевтиновна, кандидат педагогических наук, доцент кафедры общей физики УдГУ, тел. (3412) 916-139, e-mail: Novikovata72@mail.ru
