Молекулярная динамика процессов зарождения и протекания фазовых переходов в неупорядоченных системах Галимзянов Булат Наилевич

Содержание к диссертации

Введение

1 Метод моделирования молекулярной динамики 12

1.1 Вводные замечания 12

1.2 Моделирование равновесной молекулярной динамики 14

1.3 Реализация термостатирования и баростатирования 21

1.4 Неравновесная молекулярная динамика 26

1.5 Потенциалы межчастичного взаимодействия 30

1.6 Эффективные потенциалы. Методы параметризации 36

1.7 Монте-Карло моделирование 40

1.8 Выводы 42

2 Расчет параметров и характеристик фазовых переходов на основе результатов моделирования молекулярной динамики 44

2.1 Введение 44

2.2 Кластерный и структурный анализ на основе конфигурационных данных 45

2.2.1 Многогранники Вороного 46

2.2.2 Симплексы Делоне 48

2.2.3 Параметры порядка 51

2.2.4 Характеризация формы структурных образований 54

2.3 Методы определения характеристик нуклеации 55

2.3.1 Скорость стационарного зародышеобразования

2.3.2 Скорость пристегивания частиц к зародышу 58

2.3.3 Расчет времени задержки нуклеации и времени индукции 60

2.4 Анализ среднего времени первого появления зародыша новой фазы 62

2.5 Термодинамическое интегрирование в оценке свободной энергии 66

2.6 Термодинамическое интегрирование в оценке поверхностного натяжения 70

2.7 Структурные характеристики 72

2.8 Трехчастичный структурный анализ 75

2.9 Динамика трехчастичных корреляций 82

2.10Выводы 86

3 Кристаллическая нуклеация в модельных переохлажденных жидкостях и стеклах 88

3.1 Введение 88

3.2 Приведенная температурная шкала 90

3.3 Моделирование молекулярной динамики кристаллизации модельных стекол 94

3.4 Процессы структурного упорядочения и СВПП-анализ 99

3.5 Характеристики кристаллической нуклеации 103

3.6 Кинетические аспекты нуклеации 108

3.7 Выводы 112

4 Кристаллическая нуклеация в модельных системах в условиях устойчивого неравновесия. Сдвиг Куэтта 114

4.1 Вводные замечания 114

4.2 Моделирование молекулярной динамики кристаллизации стекол под сдвигом Куэтта 115

4.2.1 Детали моделирования молекулярной динамики 115

4.2.2 Результаты кластерного анализа 116

4.3 Механизмы структурного упорядочения 121

4.4 Внутреннее давление в аморфной системе под однородным сдвигом 124

4.5 Термодинамические и кинетические характеристики аморфной системы под однородным сдвигом 127

4.6 Выводы 132

5 Молекулярно-динамические расчеты конденсации паров воды 134

5.1 Введение 134

5.2 Детали моделирования капельной нуклеации паров воды 136

5.3 Расчет характеристик капельной нуклеации воды 138

5.4 Расчет поверхностных характеристик капель 144

5.5 Скорость нестационарной капельной нуклеации 146

5.6 Кинетические характеристики нуклеации и роста капель 148

5.6.1 Законы роста капель 148

5.6.2 Кинетика конденсации паров воды 149

5.7 Выводы 152

Заключение 154

Приложение 156

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

• Реализация термостатирования и баростатирования
• Анализ среднего времени первого появления зародыша новой фазы
• Процессы структурного упорядочения и СВПП-анализ
• Внутреннее давление в аморфной системе под однородным сдвигом

Введение к работе

Актуальность темы исследования. В настоящее время задача, связанная с расчетом характеристик зарождения и протекания фазовых переходов в конденсированных системах, является одной из наиболее актуальных. Согласно классической теории нуклеации фазовый переход начинается, когда кластеры формирующейся новой фазы, в результате случайных флуктуации достигают так называемого критического размера. Идентифицировать эти кластеры с помощью традиционных экспериментальных методов структурного анализа крайне сложно. Причиной этому является малый размер зародышей формирующейся новой фазы (например, кристаллитов, формирующихся в переохлажденных жидкостях и в стеклах). Также отсутствуют детальные исследования нуклеационных процессов, происходящих в аморфных системах при глубоких переохлаждениях. Здесь возникает необходимость обнаружения универсальных закономерностей в температурном поведении основных параметров, характеризующих кинетику процесса зародышеобразования (например, скорости стационарной нуклеации, времени индукции, скорости роста кластеров). Численные значения этих параметров определить экспериментальными методами практически не представляется возможным. В связи с этим, значительную роль приобретает численный метод моделирования молекулярной динамики. Основываясь на наработках теоретической физики (прежде всего статистической физики и термодинамики) метод моделирования молекулярной динамики позволяет достаточно точно рассчитать практически все характеристики фазовых переходов в конденсированных системах.¹

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

^:D. Kashchiev Nucleation: Basic Theory with Applications (Butterworth-Heinemann, Frenkel Oxford, 2000).

ции эффективных потенциалов для описания физических свойств и структурных характеристик таких систем как вода, растворы солей, многокомпонентные материалы, сплавы металлов и т. д.² Например, до недавнего времени отсутствовали исследования кинетических и термодинамических характеристик конденсации насыщенного водяного пара с использованием крупнозернистых моделей. Как правило в большинстве работ, посвященных исследованию свойств воды методами моделирования молекулярной динамики, за основу берутся атомистические потенциалы, такие как TIP4P, SPC/E и др.

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

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

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

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

2. Развит оригинальный подход в рамках метода термодинами-

²D. Frenkel, and В. Smit Understanding Molecular Simulation (San Diego: Academic Press, 2001).

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

3. Введена оригинальная приведенная температурная шкала для температурной области 0 < Т < Т_т (Т_т - температура плавления), которая учитывает как особенности плавления, так и стеклования конкретной системы. Данную температурную шкалу можно использовать в качестве альтернативы, например, переохлаждению и приведенной температуре, используемой для построения так называемого графика Энжелла (Angell plot), применяемой для оценки хрупкости стеклообразующих систем.³

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

5. Методом моделирования неравновесной молекулярной динамики рассчитаны характеристики процесса кристаллизации стекольной системы, индуцируемого внешним однородным сдвиговым воздействием.

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

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

³С.А. Angell // Science. - 1995. - Vol. 267, P. 1924.

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

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

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

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

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

Апробация работы. Основные результаты и выводы настоящей работы докладывались на Международной молодежной научной школе "Современная нейтронография" (г. Дубна, ОИЯИ, 2012 г.); Международной научной Интернет-конференции "Математическое и компьютерное моделирование в биологии и химии. Перспек-

тивы развития" (г. Казань, К(П)ФУ, 2012 г.); Всероссийской конференции "Необратимые процессы в природе и технике" (г. Москва, МГТУ им. Н.Э. Баумана, 2013, 2015); Международной Интернет-конференции "На стыке наук. Физико-химическая серия" (г. Казань, К(П)ФУ, 2013, 2015); Всероссийской школе-конференции студентов, аспирантов и молодых ученых "Материалы и технологии XXI века" (г. Казань, К(П)ФУ, 2014 г.); XVIIIth Research Workshop "Nucleation Theory and Applications" (г. Дубна, ОИЯИ, 2014 г.); Международном молодежном научном форуме "Ломоносов-2015" (г. Москва, МГУ им. М.В. Ломоносова, 2015 г.); XIV Международной школе-конференции "Проблемы физики твердого тела и высоких давлений" (г. Сочи, ИФ-ВД РАН, ФИАН, 2015 г.); I Международной школе-конференции молодых ученыхг"Биомедицина, материалы и технологии XXI века"ґ(г. Казань, КФУ, 2015 г.), а также на научных семинарах кафедры вычислительной физики К(П)ФУ (2011-2015 гг.).

Публикации. По результатам выполненных исследований опубликовано 19 печатных работ, из них 6 в научных изданиях, включенных в базу цитирования Scopus, Web of Science и РИНЦ, 3 статьи в сборниках научных работ, 10 тезисов докладов на между народных и всероссийских научных конференциях. Получено 6 свидетельств о государственной регистрации расчетных и программных комплексов для ЭВМ.

Структура диссертации. Диссертация состоит из пяти глав, заключения и списка цитируемой литературы. Работа изложена на 179 страницах, которые включают в себя 46 рисунков, 15 таблиц и список литературы из 180 наименований.

Реализация термостатирования и баростатирования

В современной физике метод моделирования молекулярной динамики имеет весьма широкую область применения [32,40,41]. В некоторых случаях этот метод выступает в качестве замены традиционному эксперименту, применимость которого иногда может быть ограниченной (например, при исследовании фазовых переходов в переохлажденных жидкостях и в стеклах) [3,26,32,41,53]. Также метод молекулярной динамики успешно применяется для воспроизведения структурно-механических свойств в атомарных и молекулярных системах при условиях, наиболее приближенных к реальным [42]. Этот метод позволяет рассчитать транспортные свойства, релаксационные параметры, временные корреляционные функции, включая структурно-динамические характеристики таких материалов как полимеры, жидкие кристаллы, металлические сплавы [14,54,55].

В последние десятилетия активное развитие компьютерного моделирования и его активное внедрение в науку открыли новые возможности в исследовании разнообразных фундаментальных проблем современной физики, например, связанных с фазовыми переходами в неупорядоченных систе 13 мах [1,7,18,32,39]. На сегодняшний день наиболее популярными методами компьютерного моделирования являются: метод моделирования молекулярной динамики, основанный на законах классической механики [32,40,56,57]; метод ab-initio моделирования (моделирование из первых принципов), где силы межчастичного взаимодействия определяются на основе решения кван-тово-механических уравнений [58]; метод стохастической динамики, который преимущественно используется для моделирования броуновского движения частиц (в коллоидном растворе) на основе стохастических дифференциальных уравнений [59,60]; метод Монте-Карло, где моделирование осуществляется на основе стохастических численных методов и вероятностных распределений [61,62].

Указанные методы позволяют выполнить расчет траекторий движения частиц (или подсистем), формирующих некоторый физический объект, и отслеживать временную эволюцию многочастичной системы [32]. При этом можно проследить за движением всех частиц системы, что экспериментально выполнить практически не представляется возможным. Траектории частиц могут быть детерминированными или определяться на основе различных вероятностных параметров (или распределений) [59,62]. В настоящей работе под термином "частицы" подразумеваются атомы, молекулы или ионы, из которых состоит система.

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

Сложности при реализации моделирования молекулярной динамики преимущественно связаны с выбором, подходящего для исследуемой системы, модельного потенциала межчастичного взаимодействия. Большинство реальных систем состоят из огромного числа частиц различных сортов, которые могут образовывать соединения со сложной структурой и геометрией, свойства которых сложно описать с помощью сферических парных или трехча-стичных потенциалов. При этом большинство потенциалов позволяют рассматривать только ограниченную область фазовой диаграммы исследуемой системы. Поэтому задача, связанная с реконструкцией эффективных потенциалов межчастичного взаимодействия, адаптированных для изучения структурно-механических свойств конденсированных систем в широкой области фазовой диаграммы, является крайне важной [33,34,51].

Развитие классического метода моделирования равновесной молекулярной динамики непосредственно связано с работами Б. Дж. Олдера и Т. Э. Уэй-нрайта [63], которые моделировали свойства системы твердых сфер в различных агрегатных состояниях. В 1964 году А. Рахманом были выполнены молекулярно-динамические расчеты жидкого аргона на основе потенциала Леннард-Джонса с 864 атомами [56]. На основе этого же потенциала в 1967 году Л. Верле рассчитал фазовую диаграмму для аргона [57]. В известной монографии М. П. Аллена и Д. Дж. Тилдеслия [40] были системно изложены базовые алгоритмы метода моделирования молекулярной динамики, которые позволили наглядно продемонстрировать применимость этого метода как для самостоятельных исследований, так и применительно к интерпретации экспериментальных результатов. Несмотря на то, что первые шаги к использованию метода молекулярной динамики были сделаны полвека назад, этот метод на сегодняшний день является весьма популярным [32].

Рассмотрим реализацию метода моделирования равновесной молекулярной динамики на примере простой модельной системы. В качестве примера рассмотрим систему, состоящую из N идентичных частиц с одинаковой массой т, расположенных внутри симуляционной ячейки. В зависимости от условий решаемой задачи симуляционная ячейка может быть задана произвольной формой (например, кубической, сферической, ромбической). В случае трехмерной прямоугольной ячейки с длинами ребер Lx, Ly, Lz и объемом V = Lx Ly Lz система будет характеризоваться численной плотностью рп = N/V. На начальном этапе положения частиц внутри этой ячейки определяются либо произвольным образом, например, с помощью генератора случайных чисел, либо задаются в соответствии с геометрией кристаллической решетки (гранецентрированной кубической, гексагональной плотно упакованной, объемноцентрированной кубической и т.д.) [40]. Скорости Vi = Vi(vx,Vy,vz) и направление движения каждой частицы задаются случайным образом или при условии, что полная кинетическая энергия системы Е соответствует ее внутренней температуре Т (Т ос Е /кв), которая является изначально заданной [32,40].

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

При рассмотрении жидкости или аморфного материала интерес представляет локальное окружение атомов (или молекул), что весьма важно при описании коллективных эффектов, возникающих на этапе зарождения фазовых переходов в многочастичных системах [1, 2, 6, 12, 53]. Например, в процессе фазового перехода корректное выявление очагов новой фазы позволяет достаточно аккуратно определить структурные характеристики, размеры и форму образовавшихся зародышей [2, 34, 53, 90] и т. д. При этом основное внимание уделяется к выявлению наличия периодичности в расположении частиц, формирующей трансляционную и вращательную симметрии, свойственные кристаллам [115]. Для того чтобы определить тип кристаллической решетки, возникающей в процессе нуклеации, нужно учитывать взаимное расположение в пространстве атомов или молекул, образующих зародыш. В молекулярно-динамических расчетах это возможно осуществить с помощью специальных алгоритмов кластерного и структурного анализа: например, метода многогранников Вороного [116,117], симплексов (триангуляции) Делоне [117-119], метода расчета параметров ориентационного порядка Стейнхардта-Нельсона [43,44].

Метод многогранников (или диаграмм) Вороного представляет собой разбиение плоскости некоторого конечного множества точек на многоугольники. Иногда многогранники Вороного сопоставляют с многоугольниками (или ячейками) Дирихле и ячейками Вигнера-Зейтца [117]. Метод построения многоугольников впервые был изучен русским математиком Вороновым в конце XIX века и на сегодняшний день используется в различных областях науки и техники (например, в вычислительном материаловедении, в геодезии, в численном моделировании и в структурном анализе) [116]. Этот метод используется для разбиения моделируемой системы на отдельные участки, а также применяется для определения ближайшего окружения атомов (или молекул) многочастичных систем. Также построение многогранников Вороного успешно применяется для выявления упорядоченных (кристаллических) областей в аморфных системах, например, в стеклах.

Рассмотрим более подробно алгоритм построения многогранников Вороного на двумерной плоскости. Предположим, что множество S содержит некоторое количество N частиц (см. рисунок 2.2.1). Точками будут называться вершины многоугольников, расположенных в множестве S. Изначально, относительно каждой частицы pi (г = 1, 2,...,7V) необходимо определить расположение (координаты и(х,у)) точек вершин многоугольников. Если имеются две частицы pi и pj, то множество точек, более близких К Рі, чем К Pj, представляет собой полуплоскость h(j)i,Pj), ограниченную прямой. Данная прямая перпендикулярна отрезку р Щ и делит его пополам. Множество точек V(i), которое является наиболее близким к частице р,;, можно определить через соотношение V(i) = f)h(pl}pJ)}j = l} 2,..., Ж (2.2.1)

Здесь множество V{i) называется многоугольником Вороного для частицы Pi. Каждая из N частиц может принадлежать лишь одному многоугольнику Вороного. Поэтому частица pi может рассматриваться в качестве ближайшего соседа некоторой точки и(х,у) при условии и Є V(i). Таким образом, количество точек и(х,у), входящих в множество V(i), представляет собой число ближайших соседей Рі-ой частицы.

Основные свойства многогранников Вороного: каждая вершина многогранников Вороного является точкой пересечения в точности трех ребер; вершины диаграммы являются центрами окружностей С (и), определяемыми тремя частицами исходного множества S; каждый ближайший сосед частицы pi определяет ребро в многоугольнике ПО; многоугольник V(і) является неограниченным тогда и только тогда, когда точка pi лежит на границе выпуклой оболочки множества S; если точка РІ не лежит на границе выпуклой оболочки множества S, то она является внутренней точкой некоторого треугольника; диаграмма Вороного, образованная N частицами, имеет не более (2N — 5) вершин и (3N — 6) ребер. Рис. 2.2.1: Схематичное изображение многогранников Вороного в двумерной плоскости. Через pi, р2, ...,Рп обозначены частицы некоторого множества S. Малыми кружками щ, U2,---,un обозначены вершины многоугольников Vi, V2 и Уз. Через С(щ) и C(v,2) обозначены окружности с центрами в точках и\ и «2 2.2.2 Симплексы Делоне

Впервые задача построения симплексов Делоне была сформирована в 1934 году и является одной из основных в вычислительной геометрии [119]. Этот метод используется при решении пространственных задач, например, при нахождении менее и более плотных областей системы, а также в структурном анализе для расчета рельефа и контура моделируемых многочастичных систем [117,119]. Симплексы Делоне, по определению, являются тетраэдрами, в вершинах которых располагаются атомы (или молекулы). Полное разбиение пространства на симплексы показывает все возможные связи частиц системы. Благодаря этому можно выделить межчастичные области, поры различных размеров и форм, а также определить области, имеющие упорядоченную структуру (например, в твердых телах, в жидкостях) [117,118].

В настоящее время существует множество алгоритмов построения симплексов Делоне. Например, итеративные алгоритмы (с индексированием и кэшированием поиска треугольников), алгоритмы слияния (с разрезанием по диаметру, выпуклого полосового слияния), алгоритмы прямого построения (клеточный пошаговый, с ускорением поиска соседей), двухпроходные алгоритмы (ленточный, веерный, рекурсивного расщепления). Однако наиболее эффективными по скорости построения являются итеративный алгоритм, алгоритм динамического кэширования, а также ленточный алгоритм [119].

Процессы структурного упорядочения и СВПП-анализ

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

Эта функция определяет вероятность нахождения пары частиц на расстоянии г = \г\ друг от друга. В численном моделировании выражение (2.7.1) заменяется соотношением =5 \5 / (2-7-2) которое подразумевает усреднение по различным парам частиц [126]. Здесь Ащ(г) - число частиц в сферическом слое толщиной Аг на расстоянии г от j-той частицы.

На основе функции радиального распределения частиц можно рассчитать средние значения различных физических величин (/): Здесь рп - численная плотность системы, U(г) - потенциальная энергия межчастичного взаимодействия, гс определяет расстояние до первого минимума в функции д{г). Следует отметить, что все вышеперечисленные выражения представлены для трехмерной системы.

В случае жидкостей и аморфных тел функция д(г) на больших расстояниях стремится к единице, а на малых расстояниях к нулю, так как частицы не могут сблизиться на дистанцию меньше их условного диаметра: Г— 0 Кроме того, функция д(г) для жидкости имеет осциллирующее поведение, где первый максимум д(г) соответствует первой координационной сфере (с координационным числом z), второй максимум определяет вторую координационную сферу и т. д.

Статический структурный фактор. Функция радиального распределения частиц связана со статическим структурным фактором S(k) через трехмерное Фурье-преобразование: где к = \к\ - волновое число [42]. Величина S(k) непосредственно измеряется из экспериментов по рассеянию нейтронов и упругому рассеянию рентгеновских лучей 1 /N - N - \ s = (Е е гк0 Ее%к0) (2-7-10) \j=i і=і I Для статического структурного фактора S{k) можно записать следующие свойства, которые аналогичны свойствам функции д(г) (при рассмотрении жидкостей и аморфных тел):

В настоящее время большинство методов кластерного и структурного анализа позволяют идентифицировать области структурного упорядочения в многочастичных системах на основе различных корреляционных функций или распределений (например, функция радиального распределения частиц, параметры ориентационного порядка) [1,42,43,54,55]. При этом преимущественно учитываются лишь парные корреляции частиц. Однако рассмотрение лишь парных корреляций не позволяет в полной мере объяснить и описать многие эффекты, наблюдаемые в многочастичных системах, например, образования дальнего порядка (как ориентационного, так и трансляционного) в структурно-упорядоченных системах (например, в квазикристаллах) [42, 47-49], возникновение динамической неоднородности в жидкостях [127,128]. Для описания этих процессов, наряду с парными корреляциями, необходимо рассматривать трехчастичные корреляции [42].

Нами предложен оригинальный метод трехчастичного структурного анализа многочастичных систем, суть которого заключается в рассмотрении расположения любых трех частиц системы и расчете площади S, образуемого ими треугольника (триплета) [48, 49]. Преимущество метода в том, что для описания трехчастичной корреляции используется лишь один параметр - площадь S. Величина S характеризует взаимное расположение трех частиц, которые образуют вершины условного треугольника. Это, в свою очередь, позволяет провести структурный анализ результатов молекулярно-динамических расчетов, а также количественно охарактеризовать структурные трансформации, возникающие в процессе фазовых переходов. Предположим, что площадь г-го треугольника Si меняется в зависимости от взаимной ориентации его вершин (т.е. частиц) и определяется формулой Герона [47-49] Si = [ВД r12))№ " rf]){k - rfl))]1/2 , (2.8.1) , / (12) . (23) . (31)x /0 . ,лт где Ц = [гl +т\ +т\ )/2 - полупериметр г-го треугольника [1\т - количество всевозможных треугольников, і Є {1, 2,..., Л г}), Т І - расстояние между вершинами г-го треугольника с условными номерами 1, 2 и 3 (т.е. длина сторон треугольника, образованного тремя частицами). Очевидно, что система может содержать различное количество комбинаций треугольников (в зависимости от числа частиц). Это, в свою очередь, может привести к большим вычислительным затратам. Поэтому величина N? задается с учетом размера исследуемой системы. Также в качестве вершин можно рассматривать частицы одного или разных сортов (например, в случае однокомпонентных или полидисперсных систем).

К примеру, на рисунке 2.8.1 показана временная зависимость величины S, полученная для жидкого и аморфного алюминия. Результаты получены в ходе моделирования динамики N = 864 атомов алюминия, взаимодействующих друг с другом через "ЕАМ-потенциал" [84]. Система приготовлена при температурах Т = 1000 К (жидкость) и Т = 100 К (аморфное тело) при постоянном давлении Р = 1 атм. Аморфный образец приготовлен через быстрое охлаждение жидкого образца со скоростью dT/dt = 1012 К/с также при постоянном давлении Р = 1 атм.

Из рисунка 2.8.1а видно, что для жидкого алюминий характерны кратковременные трехчастичные корреляции, где изменение площади триплета происходит в интервале от S — 0.25 нм до S — 1 нм . Однако в случае аморфного алюминия отсутствуют ярко выраженные осцилляции величины 5 , присущие жидкости (см. рисунок 2.8.lb). Поэтому изменение площади триплета происходит в гораздо меньшем интервале S Є [0.31; 0.34] нм . Это связано с тем, что в аморфном состоянии мобильность частиц системы очень низкая, где они в основном совершают только тепловые колебательные движения. Следует отметить, что результаты, представленные на рисунке 2.8.1, получены лишь для одного триплета. Это выполнено для того, чтобы продемонстрировать динамику S для различных агрегатных состояний.

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

Здесь S определяет площадь некоторой воображаемой окружности внутри сферы радиусом rs, где располагаются триплеты, An(S) - распределение по площадям S (по аналогии с величиной An (г) в функции д{г)). Скобки (...) обозначают усреднение по ансамблю триплетов. По аналогии сд(г) функция g(S,rs) определяет вероятность появления триплета с площадью S. Одна из вершин триплета рассматривается в качестве основной - фиксированной, относительно которой оценивается положение остальных вершин.

Внутреннее давление в аморфной системе под однородным сдвигом

Метод моделирования молекулярной динамики играет важную роль в исследовании микроскопических особенностей и свойств воды в широкой области фазовой диаграммы. Например, расчет характеристик конденсации перенасыщенного водяного пара в рамках молекулярно-динамических расчетов преимущественно выполняется с помощью различных атомистических моделей, например, SPC/E, ТІРЗР, ТІР4Р и ТІР5Р, основанных на электростатической природе межмолекулярных взаимодействий [36,167-169]. В настоящее время популярными являются идеи разработки крупнозернистых моделей, где взаимодействие между атомами кислорода и водорода учитывается эффективным образом [2,51,90].

Межмолекулярные взаимодействия в крупнозернистых моделях задаются с помощью некоторого эффективного многочастичного потенциала [89,90]. К примеру, была предпринята попытка описать структурно-динамические характеристики воды на основе одноатомной (mW) модели с модифицированным потенциалом Стиллинжера-Вебера, где взаимодействие между атомами водорода и кислорода не учитывается [33,90]. Каждая молекула воды рассматривается как отдельная частица, способная образовать с другими частицами тетраэдрические структуры [33,34,90]. Далее нами было показано, что в этой крупнозернистой модели потенциала используются оптимальные значения параметров (см. 1.6). Следует отметить, что изначально потенциал Стиллинжера-Вебера был предложен применительно к жидкому (и кристаллическому) кремнию и германию, а также для воспроизведения свойств жидкой воды, включая некоторые аномальные свойства. Наличие сходства между германием, кремнием и водой, которое выражается, в частности, в аномальном температурном поведении вязкости и плотности этих систем, позволяет рассматривать молекулы воды как частицы с короткодействующими силами взаимодействия [34,35,89,90].

Моделирование молекулярной динамики фазовых переходов, в частности, процесса конденсации водяного пара было выполнено К. Ясуоко и М. Мат-сумото на основе парного потенциала Леннард-Джонса, где рассчитывались значения скорости зародышеобразования, критического размера и нуклеаци-онного барьера при температуре Т = 350 К и давлении Р = 1 атм [37]. Также на основе mW-модели рассматривались процессы кристаллизации переохлажденной воды, выполнялась оценка значений скорости зародышеобразования Jst и критического размера пс для температурной области 150К Т 298К при давлении Р = 1 атм [34,90]. Также, в работе Ф. Зиполи [35] указывается на то, что mW-модель дает завышенные значения (порядка Jst 10 см с ) для скорости нуклеации капель по сравнению с экспериментом [140-142], к примеру, при температуре Т = 353 К и перенасыщении S = 1.2. Однако, результаты молекулярно-динамических расчетов через SPC/E модель, полученные при тех же условиях, также показали значения порядка Jst 10 см с [38]. Это, в свою очередь, указывает на необходимость проведения более детальных расчетов характеристик процесса зародышеобразования с помощью крупнозернистых моделей [2,34,141,142].

В настоящей главе представлены результаты молекулярно-динамических расчетов, где впервые была использована крупнозернистая модель воды применительно к рассмотрению процесса конденсации перенасыщенного пара в температурной области от 0С до 100С. Эти результаты опубликованы в работах [2,50,51,105,106,114,170].

Моделирование молекулярной динамики капельной нуклеации выполняется для одноатомной модели воды, состоящей из N = 8000 молекул, которые взаимодействуют через эффективный mW-потенциал (см. раздел 1.5) [33,34]. Моделируемая система представляет собой кубическую ячейку с периодическими граничными условиями по всем направлениям. Интегрирование уравнений движения осуществляется через скоростной алгоритм Верле с временным шагом t = 1 фс, где координаты и скорости движения молекул рассчитываются через уравнения (1.2.4)-(1.2.7) (см. 1.2). Расчеты выполняются в изобарическом-изотермическом ансамбле при постоянном давлении Р = 1 атм. Давление и температура системы контролируются через баростат и термостат Нозе-Гувера, определяемые выражениями (1.3.16) и (1.3.6), с параметрами Qp = 10єт2 и QT = 10єт2 соответственно. Исходная конфигурация системы, находящаяся в газообразной фазе при температуре Т = 900 К, изображена на рисунке 5.2.1, где молекулы водяного пара распределены равномерно и занимают весь объем моделируемой ячейки. Далее образцы охлаждаются со скоростью dT/dt = 1010 К/с при температурах из отрезка Т Є [273; 373] К при постоянном давлении Р = 1 атм.

Расчет траекторий молекул выполняется методом моделирования молекулярной динамики. При этом фиксируется каждое нуклеационное собы 137

Исходная трехмерная ячейка моделирования, где изображены молекулы водяного пара при температуре Т = 900 К и давлении Р = 1 атм. тие для каждой конфигурации системы при заданной температуре. Здесь под нуклеационным событием подразумевается формирование водяной капли, идентификация которой осуществляется на основе критерия Стиллинже-ра. Согласно этому критерию две молекулы принадлежат одному зародышу (т.е. одной капли воды), если расстояние между ними составляет г 3.56 А [33,34,89,171]. При этом зародышем считается образование, состоящее из четырех и более молекул.

При рассмотрении нуклеации основное внимание уделяется временным зависимостям размера капель воды n{t). На рисунке 5.2.2а изображены кривые n(t) для пяти независимых молекулярно-динамических итераций, полученные при температуре Т = 273 К и давлении Р = 1атм. Структурные изменения, происходящие вследствие перехода системы в жидкую фазу, определяются через функцию радиального распределения молекул д(г) (2.7.2) (см. 2.8). Кроме того, при фазовом переходе пар-жидкость во временной зависимости потенциальная энергия системы Ep(t) быстро уменьшается, как это показано на рисунке 5.2.2Ь. Такое поведение величины Ep(t) связано с образованием капель критического размера пс, при котором система достаточно быстро переходит в жидкую фазу