Атомистические механизмы и кинетика распада перегретых и переохлажденных жидкостей тема автореферата и диссертации по физике, 01.04.07 ВАК РФ

На правах рукописи УДК 538.9+539.42

Писарев Василий Вячеславович

Атомистические механизмы и кинетика распада перегретых и переохлажденных жидкостей

Специальность 01.04.07 - Физика конденсированного состояния

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

1 З ДЕК 2012

Долгопрудный - 2012

005056957

005056957

Работа выполнена в Федеральном государственном автономном образовательном учреждении высшего профессионального образования Московский физико-технический институт (государственный университет) на кафедре физики высокотемпературных процессов (Базовый институт Федеральное государственное бюджетное учреждение науки Объединённый институт высоких температур РАН).

Научный руководитель: доктор физико-математических наук, зав. отд.

профессор Норман Генри Эдгарович (Объединенный институт высоких температур РАН)

Официальные оппоненты: доктор физико-математических наук,

заместитель директора по науке, профессор Рыжов Валентин Николаевич (Институт физики высоких давлений РАН им. Л.Ф. Верещагина, Москва)

доктор технических наук, профессор

Белащенко Давид Кириллович

(Национальный исследовательский технологический

университет «Московский институт стали и сплавов

Москва)

Ведущая организация: Лаборатория теоретической физики им. Н. Н. Боголюбова

Объединенного института ядерных исследований, Дубна

Защита состоится 26 декабря 2012 г. в 14.00 часов на заседании диссертационного совета Д 212.156.06 при Московском физико-техническом институте (государственном университете), расположенном по адресу: 117393, Москва, ул. Профсоюзная, д. 84/32, корпус В-2

Отзывы направлять по адресу: 141700, Московская обл., г. Долгопрудный, Институтский пер., д. 9, МФТИ, Диссертационный совет Д 212.156.06.

С диссертацией можно ознакомиться в библиотеке МФТИ. Автореферат разослан 26 ноября 2012 г.

Ученый секретарь диссертационного совета кандидат технических наук Н. П. Чубинский

Общая характеристика работы

Диссертационная работа посвящена изучению механизмов распада перегретых и переохлажденных жидкостей. Исследование проведено методом молекулярной динамики (МД). Исследованы процессы зародышеобразования (нуклеации) полостей в растянутой жидкости вблизи спинодали и кристаллической фазы в переохлажденном расплаве, процесс роста полостей и кинетика формирования околокритических кристаллических зародышей. Полученные результаты сопоставлены с классической теорией нуклеации (КТН). На основе данных МД построена модель разрыва жидкости при высокоскоростном растяжении.

Актуальность работы Развитие экспериментальной техники, освоение нано-, пико- и фем-тосекундных диапазонов открывает возможности получения метастабильных состояний конденсированных веществ со значительными степенями как перегрева (растяжения), так и переохлаждения (сжатия) по сравнению с равновесным состоянием. Следует упомянуть сильные ударные волны, нагрев фемтосекундным лазером, наносекундный электровзрыв проводников ударное и изэнтропическое сжатие жидкостей и др. Такое развитие науки ставит перед теорией задачи о предсказании свойств метастабильных состояний и описания их распада.

При теоретическом описании и численном моделировании процессов разрушения используется два основных подхода. Первый — рассмотрение разрушения за счет зарождения полостей [1-3]. Скорость этого процесса обычно вычисляется по формулам КТН. Второй подход — использование критериев на основе баланса энергии [4] или импульса [5]. Эти два подхода могут применяться как по отдельности [1—3], так и вместе [б, 7]. Эти подходы не лишены проблем, поскольку эти модели откола не позволяют предсказать, например, распределение размеров пор, определяющих откол. Требует более внимательного рассмотрения и применимость классической теории нуклеации (КТН) для импульсных процессов: при достигаемых степенях метзстэбильности КТН предсказывает размер зародыша всего в несколько раз больше межмолекулярного расстояния.

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

В связи с этим МД моделирование метастабильных сред является плодотворным методом исследования, поскольку позволяет изучить микроскопические механизмы фазовых и структурных переходов в метастабильных средах [9-12].

Цели диссертационной работы

1) Исследование гомогенной нуклеации вблизи спинодали в перегретых жидкостях и в сильно переохлажденных жидкостях .

2) Анализ применимости КТН для распада при высоких степенях метастабильности.

3) Развитие многомасштабной модели вскипания жидкости при импульсном создании метастабиль-ности.

4) Определение зависимости свободной энергии поверхности кристалл-расплав от температуры.

5) Исследование кинетики формирования кристаллического зародыша в переохлажденном расплаве.

Научная новизна.

Исследована нуклеация полостей в растянутой леннард-джонсовской жидкости. Показано, что зависимость частоты нуклеации от температуры и давления аппроксимируется в форме КТН при введении зависимости поверхностного натяжения от размера критического зародыша. Предельная скорость нуклеации. определенная из МД моделирования, значительно превышает предсказание КТН.

Исследован рост полостей нэнометровых размеров в растянутой жидкости. Показано соответствие роста полостей, получаемого из МД моделирования, уравнению Рэлея-Плессе.

Предложена модель разрыва жидкости на основе подхода «нуклеация и рост». В модели используется зависимость частоты нуклеации от давления и температуры, рассчитанная методом МД, и уравнение Рэлея-Плессе для описания роста. Данная модель позволяет детально рассчитать кинетику разрушения.

Вычислена зависимость свободной энергии межфазной поверхности жидкость-кристалл от температуры для алюминия методом МД с использованием многочастичного потенциала взаимодействия. Показано увеличение свободной энергии поверхности с ростом температуры. Рост свободной энергии поверхности жидкость-кристалл с температурой для потенциала Леннарда-Джонса, обнаруженный в работе [13] вместе с полученным в настоящей работе результатом позволяет предположить, что это является общим свойством для всех веществ.

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

Практическая значимость. Результаты диссертации могут быть использованы для описания фазовых переходов при гидродинамическом моделировании импульсных воздействий (ударного, лазерного, ионного и т.п.) на вещество.

На защиту выносятся следующие основные положения:

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

2. Модель разрыва жидкости на основе подхода «нуклеация и рост», позволяющая по данным МД рассчитать кинетику разрушения жидкости на экспериментально достижимых скоростях растяжения.

3. Полученная зависимость частоты нуклеации при кристаллизации расплава алюминия при сильном переохлаждении от давления и температуры.

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

5. Зависимость среднего времени первого появления от размера зародыша, соответствующая зависимости свободной энергии зародыша от размера, имеющей несколько максимумов.

Апробация работы. Результаты диссертации докладывались на конференциях; «Современные проблемы фундаментальных и прикладных наук» (МФТИ 2007 - 2011); «Воздействие интенсивных потоков энергии на вещество» и «Уравнения состояния вещества» (Эльбрус 2009 - 2010, 2012); «Проблемы физики ультракоротких процессов в сильнонеравновесных средах» (Новый Афон 2008 -2009, 2011 - 2012); «Аэрофизика и физическая механика классических и квантовых систем» (Москва 2008-2010); Nudeation Theory and Applications (Дубна 2008 - 2012); «Физика прочности и пластичности материалов» (Самара 2009); APS Topical Meeting on Shock Compression of Condensed Matter (Нэшвилл, США 2009); New Models and Hydrocodes for Shock Wave Processes in Condensed Matter (Париж, Франция 2010); 2010 Annual Conference of the British Association for Crystal Growth (Манчестер 2010); Advanced Problems in Mechanics (Санкт-Петербург 2012)

Публикации. Материалы диссертации опубликованы в 6 статьях в рецензируемых журналах, статьях в сборниках научных трудов и в тезисах докладов российских и международных конференций. Личный вклад автора. Все представленные в диссертации результаты получены автором

ЛИЧНО-

Структура и объем диссертации. Диссертация состоит из введения, пяти глав и заключения, изложена на 99 страницах, включает 19 рисунков, 7 таблиц, библиографию из 128 наименований.

Содержание работы

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

Первая глава носит обзорный характер.

В разделе 1.1 приведены общие сведения о понятии метастабильного состояния. Приводится описание метастабильных состояний и нуклеации с точки зрения термодинамики и статистики.

В разделе 1.2 дается введение в метод МД: его формулировка, понятие МД моделирования, основные методы диагностики МД системы. Даются понятия равновесного и неравновесного МД расчета, приведены примеры расчета свойств системы в равновесных МД расчетах.

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

где "гегаашСО - число МД траекторий, в которых к моменту времени £ не произошла нуклеация, пил<А — полное число МД траекторий, соответствующих заданным давлению и температуре, £о -минимальное время, необходимое на образование зародыша критического размера. Результаты МД моделирования нуклеации в перегретой и переохлажденной жидкости хорошо соответствуют этой зависимости, что показывает стационарность нуклеации, наблюдаемой в МД моделировании, и независимость частоты нуклеации от предыстории. Среднее время жизни метастабильной фазы вычислялось по наклону зависимости в координатах 1п(пгета;„) — Ь.

В третьей главе представлены результаты МД моделирования элементарных процессов при распаде перегретой (растянутой) жидкости — нуклеации и роста пор. Эти процессы рассмотрены на примере леннард-джонсовской жидкости с потенциалом взаимодействия:

Здесь и(г) — потенциальная энергия взаимодействия двух частиц, находящихся на расстоянии г друг от друга, а сг и е — параметры, определяющиеся выбором конкретного вещества. В дальнейшем температура, давление, плотность и другие величины будут указываться в обезразмеренных единицах: Т* = кдТ/е, Р* = Рсг3/е,р* = р<т3/т, где т - масса частицы, к в - постоянная Больцмана. Время

Для сравнения результатов МД моделирования с классическими теориями проведено независимое вычисление величин поверхностного натяжения и вязкости жидкости.

"гетат(0 = ' ехр (-(« - ¿о)/ (0) >

(1)

будет указываться в единицах т = (та^/е)1^2.

Раздел 3.1 посвящен исследованию гомогенной нуклеации полостей в жидкости. Проведено МД моделирование зарождения полостей в леннзрд-джонсовской жидкости. Получена зависимость частоты нуклеации от температуры и давления. Для сравнения с классической теорией нуклеации (КТН) проведено определение поверхностного натяжения полостей независимым способом.

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

тМ = 00 1(1 1 + S/г

где то - поверхностное натяжение на плоской поверхности раздела. Зависимость частоты нуклеации от давления при постоянной температуре в КТН имеет вид

J = Joexp("Ë)' (2)

167Г73

где IV' = ^ ^ 2 - работа образования критического зародыша, 7 - поверхностное натяжение, Р - давление в жидкости, Р' - давление пара в полости.

При сравнении результатов моделирования с (2) полагалось, что давление в жидкой фазе Р равно среднему давлению на равновесном участке МД траектории, а давлением газа в пузырьке Р' можно пренебречь по сравнению с Р. Действительно, критическая температура для леннард-джонсовской системы составляет Тс* и 1.3 [14, 15]. Рассматриваемые температуры заметао ниже критической, поэтому давление насыщенного пара мало, и небольшие пузырьки являются практически пустыми.

Сравнение результатов МД с предсказанием КТН представлено на рис. 16. КТН недооценивает частоты нуклеации на несколько порядков величины, что происходит по двум причинам. Во-первых, в классической теории нуклеации величина поверхностного натяжения 7 обычно считается постоянной и принимается равной поверхностному натяжению на плоской границе 70- Учет уменьшения поверхностного натяжения для малых пузырьков приводит к уменьшению работы образования критического зародыша ( W* ос 73) и существенному увеличению величины J, поскольку работа образования зародыша входит в экспоненту в выражении для J. Во-вторых, КТН дает величины J0 ~ 1 в интервале температур 0.5 < Т* < 0.8. Результаты МД расчетов для всех температур и давлений аппроксимируются выражением (2), если учесть зависимость поверхностного натяжения от давления и принять величину предэкспоненциального множителя Jq = 3 • 105.

s

(a)

0.4 0.6 1.6 1.8 2.0 2.42.6 i.2 3.4 3.6 3.8 4.0 42 4.4 4.6 l/p'2

1.5 2 2.5 3 3.5 4 4.5 1/P'2

Рис. 1. Зависимость частоты нуклеации от давления и температуры в крупном (а) и мелком (6) масштабах. Кружки — Т* = 0.50; квадраты — Т* = 0.71; ромбы — Т* = 0.75; треугольники — Т* = 0.80. Сплошные линии — вычисления по (2) с Зц = 3 • 105 и поверхностным натяжением, зависящим от давления. Штриховые линии — вычисления по (2) с из КТН и величиной поверхностного натяжения для плоской границы раздела.

В разделе 3.2 приводятся результаты МД моделирования роста полостей малого радиуса. Получены кривые роста полости при различных давлениях. Для сопоставления с теоретической моделью проведено вычисление вязкости жидкости.

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

Вязкость вычислена методом обращенной неравновесной молекулярной динамики [16]. Получена зависимость вязкости при нулевом давлении от температуры, которая хорошо аппроксимируется формулой Френкеля-Андраде (аррениусовская зависимость). Также получена зависимость вязкости от плотности жидкости, которая соответствует формуле Бачинского

vm-w

где vm — молярный объем жидкости. В этой модели нет явной зависимости от температуры, но при постоянном давлении она возникает за счет теплового расширения. При дальнейшем анализе использовалась зависимость вязкости от плотности (в безразмерных единицах) l/rj ~ 2.2/р 2.1, которая соответствует результатам МД моделирования при температурах Т* от 0.5 до 0.8.

Полученный в МД расчетах рост полости сопоставлен с уравнением Рэлея-Плессе, описывающим рост полости в несжимаемой жидкости;

2 р R pR р

42-

(а)

12-q 10: 8-

1.2 1.4 1.6 1.8 2.0 Vt

0.6-

(b)

л -о-jA

□ -1

о -2

л -3

о -л

------5

1/P

Рис. 2. Зависимость вязкости от температуры при Р = 0 (а) и от плотности при различных температурах (Ь). (а): квадраты — МД моделирование, штриховая линия — аппроксимация т) = 0.145•е212/"г. (b) Т*: 1 - 0.5, 2 - 0.6, 3 - 0.7, 4 - 0.8, штриховая линия (5) — единая аппроксимация l/rj — 2.2/р— 2.1.

Показано хорошее соответствие решения уравнения (3) с величинами вязкости и поверхностного натяжения, вычисленными методом МД, с кинетикой роста полостей в прямом МД моделировании (рис.З). На временах МД расчетов наблюдался выход на почти линейный рост радиуса полости, что соответствует асимптотическому решению уравнения Рэлея-Плессе. Часто используемая в моделях роста полостей в металлах упрощенная модель «вязкого роста», в данном случае неприменима, т.к. дает слишком быстрый рост.

R-Ro, а

Рис. 3. Сравнение роста полости в МД расчетах с уравнением Рэлея-Плессе (а) и законом вязкого роста (б), (а) Сплошные линии — МД моделирование, штриховые линии — вычисления по (3). (б): 1 — вязкий рост, 2 — МД моделирование, 3 — уравнение Рэлея-Плессе (3)

Четвертая глава посвящена моделированию процесса разрыва жидкости при растяжении. Рассмотрены прямое МД моделирование и многоуровневый подход к описанию разрушения. Анализируются ограничения прямого МД моделирования с точки зрения вычислительной сложности задачи.

В разделе 4.1 представлены результаты прямого численного моделирования растяжения жидкости с постоянной скоростью £ = V/Уо методом МД для достаточно больших систем.

Моделирование растяжения проведено для режима постоянной температуры Т* = 0.707. При растяжении по мере роста степени метастабильности в системе происходит зарождение и рост полостей (рис.4). При этом рассмотрены размеры системы, когда появление первого зародыша не приводит к сильному снижению метастабильности, препятствующему дальнейшему зародышеобразо-ванию, так что разрушение обусловлено зарождением и одновременным ростом множества полостей. Получены величины максимальных растягивающих напряжений в зависимости от скорости растяжения. Показано, что наблюдается рост прочности с увеличением скорости растяжения с выходом при сверхвысоких скоростях растяжения на постоянное значение, равное спинодальной прочности.

Рис. 4. Визуализация зарождения и роста полостей в расчетной ячейке. Показаны только частицы, лежащие на поверхности полостей. Снимок в момент £* = 120 соответствует моменту, когда достигнуто максимальное напряжение.

Приведена оценка роста вычислительной сложности прямого МД моделирования разрушения при уменьшении скорости растяжения. Рост сложности определяется, во-первых, увеличением характерного времени процесса, а во-вторых, ростом характерного масштаба — среднего расстояния между полостями. За счет этих факторов вычислительная сложность прямого МД моделирования разрушения растет с уменьшением скорости растяжения как £ . Для прямого моделирования скоростей растяжения ~ 10б с-1, близких к максимальным для классических ударно-вол новых экспериментов с использованием химических взрывчатых веществ, необходим расчет с 10 триллионами частиц, при этом необходимо около 10' шагов численного интегрирования. Решение такой задачи за разумное время лежит за пределами возможностей даже мощнейших из современных суперкомпьютеров.

В разделе 4,2 предложена двухуровневая модель на базе подхода «нуклеация и рост полостей» (nudeation and growth — NAG) для вычисления динамической прочности жидкости и описания кинетики разрушения при меньших скоростях деформирования. Первый уровень составляют данные о кинетике зарождения и роста, полученными для единичных полостей в МД расчетах. Второй уровень представляет собой уравнения, связывающие кинетику зарождения и роста полостей и историю нагружения.

Суммарный объем полостей в момент времени t можно вычислить как

t

Vcav(t) = J Vcav(t - T)n{r)d,T, (4)

0

где vcav(i — т) - объем полости, возникшей в момент г, к моменту t, п(г) - скорость появления полостей в образце объема V.

Объем единичной полости и скорость появления полостей выражаются через известные частоты нуклеации и скорости роста:

4

Vcav(t - г) = -w(R(t - т))\ h{r) = J{P{t))V, (5)

где R(t — r)- радиус полости, получаемый решением уравнения (3) с начальными условиями Л(т) = 27(Р(г))/Р(г),Д(т) = 0.

Таким образом, при известных уравнении состояния материала P(V,T) и истории нагружения выражения (4) - (5) позволяют вычислить объем полостей и их распределение по размерам в произвольный момент времени. Далее можно, например, вычислить максимально достижимое напряжение в материале. При объединении предложенной модели с методами механики сплошной среды можно рассчитать профили плотности и скорости свободной поверхности при ударно-волновом нагружении [17].

В качестве примера применения модели рассмотрен частный случай растяжения с постоянной скоростью V/VQ = è = const. Введем момент откола tsp как момент, когда скорость нарастания объе-

I dVcav

ма пор сравнивается со скоростью растяжения: --—-—{t,p) = é. Зная из уравнения состояния ЗЗВИ-

Vo «f

симость P(t), по предложенной модели можно рассчитать достигнутое напряжение в момент откола, т.е. откольную прочность. Такой критерий откола позволяет оценить максимальные растягивающие напряжения, достижимые в материале, и соответствует критерию, использованному в прямом МД моделировании.

Данные МД моделирования, полученные для леннард-джонсовской жидкости, перемасштабированы для описания жидкого гексана с использованием следующих параметров потенциала: с/кд = 413 К, и = 5.909 нм [18]. Зависимость откольной прочности гексана от скорости растяжения при тем-

пературе Т = 292 К, вычисленная по разработанной модели, показана на рис.5а. Хорошее согласие результатов модели NAG с экспериментом и прямым МД моделированием доказывает ее применимость для описания разрушения при скоростях растяжения, достигаемых в эксперименте.

(б)

25

20

15

Г Ps, МПа

(а)

Спинодаль

NAG Прямой МД расчет

е, с'

?0з

ю4

105 106 107

108

109

0.1

О 2.4 4.8 7.2 9.6 12 14.4 Радиус полости (нм)

Рис. 5. (а): Зависимость откольной прочности гексана от скорости растяжения при Т — 292 К. 1 — вычисление по модели NAG, 2 - вычисление по нуклеационной модели (см. текст), 3 - энергетический критерий откола, 4 — экспериментальные данные [1], 5 — прямое МД моделирование, б - спинодаль. (6): Распределение полостей в жидком гексане по размерам в момент откола t = tsp при скорости растяжения Е = 3 ■ 10® с"1.

Для сравнения на рис. 5а показана откольная прочность гексана, рассчитанная по энергетическому критерию откола [4]. В этом случае откольная прочность дается выражением:

Psp = (6p4V)1/3 (6)

где Iо - плотность жидкости, cs — скорость звука, 7 - поверхностное натяжение.

Как видно, откольная прочность гексана, вычисленная по энергетическому критерию, в несколько раз ниже экспериментально наблюдаемой величины. Это, по-видимому, связано с тем, что при выводе формулы (б) рассматривается предельный случай, когда работа, совершаемая при растяжении жидкости, в основном тратится на образование новой поверхности. Это предположение может быть при наличии в жидкости большого количества гетерогенных центров нуклеации (микропузырьков газа, твердых частиц и т.д.) с самого начала. Тогда даже при небольшом растяжении происходит активное зародышеобразование и рост поверхности полостей. В нашем случае лимитирующей стадией разрушения является ожидание появления полостей по гомогенному механизму. При этом работа, совершаемая при растяжении, в основном не идет на образование поверхности, а запасается в виде

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

Модель NAG предсказывает достижение откольной прочностью теоретического предела прочности леннард-джонсовской жидкости при скорости деформирования é яз 3 • 108 с-1. Это означает, что при более высоких скоростях растяжения откольная прочность должна оставаться постоянной. Это предсказание модели подтверждается при прямом МД расчете разрушения (см. рис. 5а).

Прямой проверкой модели NAG является сравнение предсказываемого распределения пор по размерам в некоторый момент времени с распределением, которое получается в прямых МД расчетах. Такая проверка была сделана для момента откола при скорости растяжения é = 6 ■ 107 с-1. Сравнение полученных распределений размеров пор показано на рис. 5Ь. Вычисления по модели NAG демонстрируют полное согласие с результатами прямого МД моделирования. Как видно, в момент откола в системе имеется большое количество малых полостей и несколько более крупных. Радиус полости определяется временем ее появления. Так как величина растягивающих напряжений растет со временем, то растет и частота нуклеации, вследствие этого большая часть полостей — это те, которые образовались на более поздних этапах растяжения. Следует заметить, однако, что общий объем малого количества крупных полостей, которые образовались на более ранних этапах, превышает общий объем множества недавно зародившихся, малых полостей. Поэтому в момент откола основной вклад в объем пустот и скорость его нарастания вносит не зарождение новых пор, а рост тех, что образовались несколько раньше момента откола. Для иллюстрации этого на рис.5а показана откольная прочность, вычисленная по нуклеационной модели, не учитывающей рост полостей (т.е. в (5) полагается i?(f — г) = RCt(t) = const). Такое упрощение модели NAG приводит к переоценке откольной прочности на ~20% и расхождении предсказаний модели с результатами прямого МД моделирования и экспериментов.

В разделе 4.3 рассматривается приближенная аналитическая зависимость откольной прочности от скорости растяжения.

Зависимость, предложенная в работе [1] на основе нуклеационной модели разрушения без рассмотрения роста пор, имеет вид:

Откольная прочность в модели NAG также приближенно описывается этим выражением. Для обоснования применимости этой формулы проведен анализ модели NAG, где для упрощения полагается, что полости, определяющие скорость разрушения в момент откола, зарождаются в интервале длиной At на временах trn ~ o-t9p. В этом случае зависимость откольной прочности от скорости растяжения имеет тот же вид, что и в нуклеационной модели, при этом по величине коэффициента Л можно оценить величину а. Для режима, рассмотренного в работе, получена величина а г» 0.8.

В пятой главе представлены результаты исследования нуклеации при кристаллизации и стеклования. Рассмотрен пример жидкого алюминия. Для МД моделирования использовался многочастичный потенциал погруженного атома (ЕАМ) в параметризациях, предложенных в [19, 20].

В разделе 5.1 приводится общее описание двух возможных режимов затвердевания жидкостей при высокоскоростном охлаждении: кристаллизации и стеклования.

Раздел 5.2 посвящен определению частоты нуклеации кристалла в переохлажденном расплаве. Вычислены частоты нуклеации в зависимости от давления вдоль изотерм. Проведен анализ образующихся кристаллических кластеров.

Частоты нуклеации вычислены для изотерм 670, 750 и 900 К. Проведено сравнение полученной зависимости частоты нуклеации от температуры и давления с моделью классической теории нуклеации. В КТН частота нуклеации представляется в виде:

ДР,Т) = МР,Т)ехр(-^у (7)

где ДС — свободная энергия образования критического зародыша. В КТН используется следующее выражение для энергетического барьера:

= <»>

где 7 — свободная энергия единицы поверхности кристалл-расплав, и — атомный объем в кристалле, А/г — разность химических потенциалов расплава и кристалла.

В этом выражении ■у и Д^ являются неизвестными величинами. Однако разность химических потенциалов можно оценить как

Д^-Д^-^-Г)/^, (9)

где Дкт — теплота плавления (на один атом) на кривой сосуществования, Тт — температура плавления. Величина теплоты плавления взята из [19]. Температура плавления в зависимости от давления оценивалась по формуле Клапейрона-Клаузиуса:

Тт(Р) -Т,„о~ Т^Р&Ут/Акп, (10)

где Тт(Р) — температура плавления при давлении Р, I— температура плавления при атмосферном давлении, Аут — изменение атомного объема при плавлении. Величины Тто и Avm взяты из [19].

Таким образом, в формулах (7) - (9) остаются неизвестными 'только свободная энергия единицы поверхности кристалл-расплав у и кинетический множитель 7о- Для аппроксимации результатов МД моделирования использовались формулы (7) — (10), где величины ^ и 7 рассматривались как свободные параметры. Результат аппроксимации результатов МД моделирования в форме (7) — (10)

Рис. б. Зависимость частоты нуклеации от давления. Кружки — Т = 670 К; треугольники — Т = 750 К; ромбы — Т = 900 К; сплошные линии — аппроксимация результатов моделирования формулой КТН.

с показаны сплошными линиями на рис. 6- Из этой аппроксимации получены следующие оценки величины свободной энергии поверхности расплав-кристалл:

7(670 К) ~ 97 мДж/м2

7(750 К) ~ 114 мДж/м2

7(900 К) ~ 155 мДж/м2

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

Для двух запусков (при Г=670 К, р=0.0568 ат./А3, Р=18.3 кбар) проведен анализ размеров образующихся кристаллических зародышей. Для выделения атомов, имеющих кристаллическое локальное окружение, использовался метод анализа общих соседей (common neighbor analysis — CNA) [21]. На основе CNA выделялись атомы, имеющие кристаллическое окружение. Атом считался принадлежащим кристаллическому зародышу, если на расстоянии З.б А от него находился хотя бы еще один «кристаллический» атом. Анализ показывает, что в жидкости появляются зародыши размером 3-15 атомов, которые являются докритическими и быстро распадаются. По результатам расчетов можно оценить размер критического зародыша в 20 - 30 атомов, поскольку на метастабильном участке траектории размеры кристаллита не превышают величину 20 атомов, а зародыш размером более 30 атомов растет необратимо. В расчетах, стартующих из макроскопически идентичных состояний, образование критического зародыша происходит в различных местах расчетной ячейки. При

исследованных переохлаждениях фазовый переход инициируется образованием компактного зародыша (рис.7), в отличие от более глубоких переохлаждений, для которых в работе [11] был обнаружен режим распада через «коагуляцию» атомов с кристаллическим окружением.

Рис. 7. Визуализация околокритических зародышей. Показаны только те атомы, локальное окружение которых имеет симметрию ГЦК решетки.

Проведено сопоставление параметров нуклеации, вычисляемых по формулам КТН и полученных при МД моделировании. Показано, что в КТН переоценивается размер критического зародыша (в 1,5 - 2 раза) и отличается величина кинетического множителя </о, а также неверно воспроизводится ее температурный ход.

Таблица 1. Сравнение результатов МД и КТН

Величина КТН МД

Число частиц в критическом зародыше при 40 20-30

Г=670 К, Р=18,3 кбар

30, (пс-1) при 670 К 0,14 0,011

.¡о, (пс-1) при 750 К 0,24 0,117

■10, (пс-1) при 900 К 0,44 1,48-Ю3

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

Вычисление проводилось методом капиллярных флуктуации. В основе метода лежит исследование спектра капиллярных флуктуаций в квазидвумерной двухфазной системе [22, 23]. Квазидвумер-ность означает, что размер системы 6 вдоль направления Оу гораздо меньше размеров системы вдоль

осей Ох и Ог. В случае, когда среднее направление нормали к поверхности раздела параллельно оси Ог, форма поверхности определяется зависимостью > ^ спектр капиллярных флуктуаций

связан со свободной энергией поверхности следующим соотношением:

квт

(11)

ЪЧ? {у+ -/"(&))№'

где - преобразование Фурье функции га„г^асе(з;), IV - размер системы вдоль оси Ох, кд — постоянная Больцмана, Т - температура, 7 — свободная энергия поверхности, к - волновое число, угловые скобки означают усреднение по времени. Из этого выражения находится величина 7-н-у"(0), называемая поверхностной жесткостью. Будем предполагать, что зависимость свободной энергии поверхности от ориентации 7(11) имеет вид [13, 23]

7(") = Ю

1 + МЕ

(12)

где щ - компоненты вектора нормали к поверхности в координатах, связанных с кристаллографическими осями. Тогда для нахождения 7(п) необходимо знать величины 70, ех и 62, которые могут быть получены при вычислении поверхностной жесткости для трех различных ориентации межфазной границы.

Таблица 2. Вычисленные значения параметров 70, 61, 62, а также величин анизотропии свободной

Величина Г = 935 К Т = 1027 К Т = 1105 К

7о (мДж/м2) 82(2) 89(3) 104(4)

«1 0.061(7) 0,08(1) 0,06(1)

€2 -0.003(1) -0.002(2) 0,002(2)

7юо - Tiio (мДж/м2) 2,1(3) 3,4(5) 3,3(7)

7мо - 7111 (мДж/м2) 3,5(4) 4,9(6) 3,9(8)

Методом молекулярной динамики вычислены величины поверхностной жесткости для ориен-таций межфазной поверхности в плоскостях (100), (110) и (111), а по ее величинам рассчитаны величины свободной энергии поверхности кристалл-расплав различных ориентации. Обнаружена анизотропия свободной энергии поверхности: уюо > 7110 > Tin Также обнаружено увеличение средней величины свободной энергии межфазной поверхности с ростом температуры.

Раздел 5.4 посвящен анализу кинетики формирования кристаллического зародыша в расплаве. Рассмотрено среднее время первого появления зародыша в зависимости от размера.

Метод анализа средних времен появления предложен и обоснован в работах Бартепла и By [24], а также Ведекинда и др. [25]. Для процесса, имеющего единственный максимум на профиле

О 10 20 30 40 SO

Число атомов в кластере

Рис- 8. Вычисленные зависимости среднего времени первого появления от размера зародыша для различных давлений при Т = 900 К. Линии - аппроксимация выражениями (13) и (14).

свободной энергии, получено следующее соотношение:

т(х) = Tj/2 ■ [1 + erf (Ъ (х - ж*))] (13)

где т(х) - среднее время первого пересечения (mean first-passage time - MFPT) координаты x, x* -координата переходного состояния, tj — величина, обратная к частоте перехода барьера, коэффициент Ъ связан с множителем Зельдовича Z = {d2G(x)/дх2)„*/2тгкТ выражением Ъ = Zy/ж. Функция (13) имеет точку перегиба при х = х*, при этом т{х") = tj/2. Таким образом, нахождение средних времен первого пересечения координаты реакции и их аппроксимация в виде (13) позволяет сразу получить и частоту перехода барьера, и координату переходного состояния.

Для выделения кристаллической фазы использован метод CNA. Атомы с кристаллическим окружением, находящиеся на расстоянии, меньшем 3.70 А, друг от друга, считались принадлежащими одному и тому же кластеру (плотность жидкости соответствовала среднему межатомному расстоянию 2.49 А). По атомным конфигурациям анализировались размеры и количество кластеров и вычислялось среднее время первого появления кластера в зависимости от размера. Полученные зависимости MFPT от размера п при различных давлениях при температуре Т = 900 К приведены на рис.8. Кривые, соответствующие высоким давлениям (или высоким переохлаждениям) хорошо соответствуют формуле (13). При снижении степени переохлаждения расплава зависимость MFPT(ra) начинает отклоняться от формы (13), на ней появляются «ступеньки» и дополнительные точки перегиба. Небольшие «ступеньки» видны уже при Р = 12.8 ГПа, при понижении давления до Р = 12.1 ГПа они становятся более выраженными. Вместо формулы (13) зависимость MFPT(n) представляется в виде:

= ~ [1 + erf(bi (п - «!))] + T-f [1 + erf(b2 (n - П2))} + ... (14)

т.е. появляются дополнительные «критические» зародыши с размерами п*. Заметим, что положение образующихся «ступенек» (т.е. дополнительных «критических» размеров) практически не меняется при изменении давления, но при этом меняется их относительная высота. При рассмотрении кинетики зародышеобразования в рамках подхода Дёрингз-Фольмера-Зельдовича показано, что полученная ступенчатая зависимость среднего времени первого появления от размера может быть интерпретирована как наличие более сложной зависимости свободной энергии зародыша от размера, чем полагается в КТН, и наличия нескольких максимумов на этой зависимости.

В разделе 5.5 приведены результаты исследования процесса стеклования расплава при высокоскоростном охлаждении в методе МД. Рассмотрены критерии перехода жидкость-стекло в МД расчетах.

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

о

Основные результаты и выводы работы

Проведено исследование распада перегретых и переохлажденных жидкостей методом классической молекулярной динамики. Исследованы режимы стационарной нуклеации, когда частота нуклеа-ции не зависит от времени наблюдения и, следовательно, истории создания метастабильного состояния. Обнаружены заметные отличия от классической теории нуклеации (КТН). Получены следующие результаты.

1. Исследована кинетика нуклеации полостей в простой жидкости на примере леннард-джонсов-ской системы частиц. Вычислены частоты нуклеации в интервале температур 0.5 < кТ¡с < 0.8 в широком диапазоне отрицательных давлений. Зависимость частоты нуклеации от давления и температуры можно аппроксимировать формулой КТН только после введения поправок на уменьшение поверхностного натяжения малых зародышей и изменения величины предэкспо-ненциального множителя в аррениусовской зависимости.

2. Построена модель типа «нуклеация и рост» для описания кинетики разрыва жидкости при растяжении с использованием данных о механизмах зарождения и роста полостей, полученных методом молеклярной динамики.

3. Получена зависимость частоты нуклеации при кристаллизации расплава алюминия от температуры и давления. Зависимость мажет быть аппроксимирована в виде КТН с двумя подгоночными параметрами: кинетический множитель и свободная энергия поверхности жидкость-кристалл. Полученная зависимость кинетического множителя от температуры отличается от предсказываемой КТН. Наблюдаемый при молеклярно-динамическом моделировании размер критического зародыша в 1,5 — 2 раза меньше, чем размер, оцененный в КТН.

4. Рассчитана свободная энергия поверхности кристалл-расплав для алюминия на кривой плавления при температурах 935 К - 1110 К. Оценена анизотропия поверхностной энергии. Энергии основных поверхностей расположены в порядке 71оо > 7U0 > 7ш- Наблюдается увеличение свободной энергии поверхности с ростом температуры и давления.

5. Методом среднего времени первого появления (MFPT) исследована кинетика формирования кристаллического зародыша. Полученная зависимость MFPT от размера может быть интерпретирована в рамках подхода Дёринга-Фольмера-Зельдовича как наличие нескольких максимумов на зависимости свободной энергии кристаллического зародыша от размера, в отличие от классической теории нуклеации, в которой предполагается единственный максимум.

Публикации автора по теме диссертации

В реферируемых журналах:

В. В. Писарев. Определение свободной энергии поверхности кристалл-расплав // Теплофизика высоких температур. — 2012. — Т. 50, № 4. — С. 769-774.

Г Э. Норман, В. В. Писарев. Молекул я рно-динамический анализ кристаллизации переохлажденного расплава алюминия // Журнал физической химии.— 2012. — Т. 86, № 9. С. 1578-1583. /. V. Morozov, А. М. Kazennov, R. G. Bystryi, G. Е, Norman, V. V. Pisarev, V. V. Stegailov. Molecular dynamics simulations of the relaxation processes in the condensed matter on CPUs // Computer Physics Communications. — 2011. —Vol. 182, no. 9. — Pp. 1974-1978.

А. Ю. Куксин, Г. Э. Норман, В. В. Писарев, В. В. Стегайлов, А. В. Янилкин. Кинетическая модель разрушения простых жидкостей // Теплофизика высоких температур. — 2010. — Т. 48, № 4, — С. 536-543. A. Yu. Kuksin, G. Е. Norman, V. V. Pisarev, V. V. Stegailov, A. V. Yanilkin. Theory and molecular dynamics modeling of spall fracture in liquids // Phys. Rev. B. — 2010. — Vol. 82, no. 17. — 174101 [1-10]. А. В. Янилкин, П. А. Жиляев, А. Ю. Куксин, Г. Э. Норман, В. В. Писарев, В. В. Стегайлов, Применение суперкомпьютеров для молекулярно-динамического моделирования процессов в конденсированных средах // Вычислительные методы и программирование. — 2010.— Т. 11, № 1. — С. 111-116.

В реферируемых сборниках:

А. Ю. Куксин, Г. Э. Норман, В. В. Писарев, В. В. Стегайлов, А. В. Янилкин. Кинетическая модель откола в простых жидкостях / Физико-химическая кинетика в газовой динамике. 2010. Т. 9. http://www.chemphys.edu.ru/medi3/files/2010-01-12-014.pdf V. V. Pisarev, A. Yu. Kuksin, G. Е. Norman, V. V. Stegailov, A. V. Yanilkin. Microscopic theory and kinetic model of spall in liquids / In Shock Compression of Condensed Matter — 2009. Ed. by M.D. Furnish et al. American Institute of Physics, New York: 2011. Pp. 801-804.

В сборниках статей:

G. E. Norman, V. V. Pisarev. Homogeneous crystal nucleation in a supercooled aluminum melt. A molecular dynamics study / Physics of Extreme States of Matter - 2011. Ed. by ас. V. E. Fortov et al. — Chernogolovka: 2011. —Pp. 121-124.

A. Yu. Kuksin, P. R. Levashov, V. V. Pisarev, M. E. Povarnitsyn, A. V. Yanilkin, A. S. Zakharenkov. Model of fracture of liquid aluminum based on atomistic simulations / ibid. — Pp. 57-59.

Г. Э. Норман, В. В. Писарев, кинетика гомогенной кристаллизации переохлажденного расплава алюминия / В сб. трудов Всероссийской школы-семинара «Аэрофизика и физическая механика классических и квантовых систем — 2010» М.: 2011.

А. Ю. Куксин, Г. Э. Норман, В. В. Писарев, В. В. Стегайлов, А. В. Янилкин Атомистическое моделирование гомогенной нуклеации / ОИВТ РАН: Итоги и перспективы. М.: 2010. С. 360-383.

G. Е. Norman, V. V. Pisarev, V. V. Stegailov, А. V. Yanilkin. Spall kinetics in liquids at high strain rates / In Physics of Extreme States of Matter — 2009. Ed. by ас. V. E. Fortov et al. Chernogolovka: 2009. Pp. 105-108.

Список литературы

[1] А. В. Уткин, В. А. Сосиков, and А. А. Богач, ПМТФ 44, 27 (2003).

[2] А. В. Уткин, В. А. Сосиков, ПМТФ 46, 29 (2005).

[3] В. Е. Виноградов, Метзстзбильные состояния и фазовые переходы. Вып. 10 (УрО РАН, Екатеринбург. 2007), С. 39.

[4] D. Е. Grady, J. Mech. Phys. Solids 36, 353 (1988).

[5] А. А. Груздков, Ю. В. Петров, ЖТФ 78, 6 (2008).

[6] М. Е. Povarnitsyn, К. V. Khishchenko, P. R. Levashov, Int. J. Impact. Eng. 38, 1723 (2005).

[7J M. E. Povarnitsyn, Т. E. Itina, M. Sentis, P. R. Levashov, К. V. Khishchenko, Phys. Rev. В 75, 235414 (2007).

[8] А. И. Фунтиков, TBT 49, 454 (2011).

[9] Yu. D. Fomin, N. V. Gribova, V. N. Ryzhov, S. M. Stishov, D. Frenkel, J, Chem. Phys. 129, 064512 (2008).

[10] Д. К. Белащенко, ЖФХ 80, 2207 (2006).

[11] Д. К. Белащенко. О. И. Островский, ЖФХ 82, 443 (2008).

[12] К. Т. Kholmurodov, К. Yasuoka, X. С. Zeng, J. Chem. Phys. 114, 9578 (2001).

[13] R. L. Davidchack, В. B. Laird, J. Chem. Phys. 118, 7651 (2003).

[14] А. Ю. Куксин, Г. Э. Норман, В. В. Стегайлов, ТВТ 45, 43 (2007).

[15] V. G Baidakov, S. P. Protsenko, Phys. Rev. Letters 95, 015701 (2005).

[16] F. MCiller-Plathe, Phys. Rev. E 59, 4894 (1999).

[17] A. Yu. Kuksin, P. R. Levashov, V. V. Pisarev, M. E. Povarnitsyn, A. V. Yanilkin, A. S. Zakharenkov, in Physics of Extreme States of Matter - 2011, edited by ас. V. E. Fortov et al. (Chernogolovka. 2011), pp. 57-59.

[18] Справочник химика. Т. 1 (M., Л.:Химия, 1966), Под ред. Никольского Б. П.

[19] М. I. Mendelev, М. J. Kramer, С. A. Becker. М. Asta, Philosophical Magazine 88,1723 (2008).

[20] X. Liu, F. Ercolessi, J. B. Adams, Model. Simul. Mater. Sci. Eng. 12, 665 (2004).

[21] D. Faken, H. Jonsson, Comput. Mater. Sci. 2, 279 (1994).

[22] J. Hoyt, M. Asta, A. Karma, Phys. Rev. Lett. 86, 5530 (2001).

[23] J. R. Morris. X. Song, J. Chem. Phys. 119, 3920 (2003).

[24] L. S. Bartell, D. T. Wu, J. Chem. Phys. 125, 194503 (2006).

[25] J. Wedekind, R. Strey, D. Reguera, J. Chem. Phys. 126, 134103 (2007).

[26] J. W. P. Schmelzer, J. Chem. Phys. 136, 074512 (2012).

Подписано в печать 19 ноября 2012 г. Объем 1,2 п.л. Тираж 100 экз. Заказ № 651 Отпечатано в Центре оперативной полиграфии ООО «Ол Би Принт» Москва, Ленинский пр-т, д.37