А.2. Формализм.94.
А.З. Реализация.96.
А.4. Утилиты командной строки .102.
А.5. Полезные ссылки .107.
1.
Введение
.
Для численного моделировании конденсированных сред в настоящее время широко используются методы молекулярной динамики. Их преимуществом является базирование на так называемых «первых принципах», когда задавая явно потенциалы межчастичного взаимодействия, можно изучать динамику различных процессов, получать уравнения состояния. Основной проблемой молекулярной динамики, ограничивающей ее применение для решения практических задач, является высокая вычислительная сложность: вообще говоря, необходимо решать систему из порядка 1023 уравнений движения, тогда как существующие вычислительные мощности позволяют решать систему из лишь порядка Ю74″ 10 уравнений. Проведение расчетов в некоторой ограниченной области с периодическими граничными условиями навязывает системе «сверхдальний порядок» и может приводить к возникновению различных побочных эффектов.
Для решения практических задач применяются, как правило, методы сплошной среды, адекватно описывающие самые разные физические процессы. Но серьезной проблемой является необходимость привлечения предположения о локальном термодинамическом равновесии, возникающего из-за использования уравнений состояния, которые могут быть получены экспериментально или численно. А между тем, в конденсированных средах существует большое количество различных фазовых переходов, для адекватного описания динамики которых необходимо учитывать нарушение термодинамического равновесия.
Молекулярная динамика оперирует с избыточным объемом информации (координаты и импульс каждой частицы), что и приводит в свою очередь к избыточной вычислительной сложности задачи. Методы сплошной среды, напротив, используют минимум информации (средние по физически бесконечно малому объему скорость, плотность и т. д.), что существенно уменьшает вычислительную сложность задачи, но позволяет описывать лишь процессы сохраняющие локальное термодинамическое равновесие.
Для газов и плазмы переход от гамильтоновой системы к уравнениям сплошной среды возможен при помощи кинетической теории. Кинетические уравнения, описывающие эволюцию функций распределения, могут быть получены из уравнений Гамильтона различными методами, например при помощи цепочки Боголюбова. Функция распределения как раз является тем минимумом информации, который необходим для описания неравновесных процессов с приемлемой вычислительной сложностью. Уравнения гидродинамики могут быть получены из кинетических уравнений методами Чепмена-Энскога или Боголюбова [1].
Для широкого круга физических проблем (фазовые переходы в магнетиках и сегнетоэлектриках, напыление тонких пленок, различные мар-тенситные переходы и т. д.) возникает необходимость решать эволюционную, динамическую неравновесную задачу, поскольку система имеет много различных равновесных и квазиравновесных решений, и без решения эволюционной задачи невозможно определить какое из них будет в итоге реализовано. В этом случае методы молекулярной динамики часто не позволяют получить решение из-за черезмерной вычислительной сложности. Представляется целесообразным использовать кинетические уравнения типа Фоккера-Планка, левая часть которых получена при помощи квазиклассического корреляционного метода песимметризованного самосогласованного поля (КНСП) [2]. Фактически такой подход означает расщепление по физическим процессам. Воздействие температуры рассматривается как воздействие окружающей среды на броуновские частицы. Уравнения Фок-кера-Планка являются традиционным инструментом для описания броуновского движения [3]. Использование уравнений типа Фоккера-Планка позволяет прозрачным образом, явно, задавать статистику описываемой системы (Больцмана, Бозе или Ферми). Кроме того для уравнений Фоккера-Планка хорошо разработаны численные методы решения [4].
1.1. Состояние вопроса.
Попытки применить кинетическую теорию для описания конденсированных сред предпринимались неоднократно, но не все из них были успешными. Так, использование уравнения Власова в статистической терии кристаллов [5, 6] приводит к серьезным противоречиям. Одновременное использование предположения мультипликативности и симметричности функции распределения означает, что каждый атом может с равной вероятностью находиться в любой ячейке кристаллической решетки, при этом существует значительная вероятность того, что несколько атомов могут находиться в одной ячейке. В результате, в частности переоценивается вклад отталкивания межатомных потенциалов во внутреннюю энергию кристалла, а при расхождении потенциала межатомного взаимодействия на малых расстояниях (как например для потенциала Леннарда-Джонса) возникает расходимость соответствующих интегралов.
Базарову [7] отчасти удалось преодолеть эти противоречия, введя обрезание одночастичной функции распределения на некотором расстоянии от узла кристаллической решетки, определяемом из вариационного принципа. Но основные проблемы остались, в частности атомы все равно могли «туннелировать» из ячейки в ячейку.
Первые успехи были достигнуты на основе идеи Терлецкого [8, 9], Стру-минского [10] и Кога [11] о необходимости использования несимметричной УУ-частичной функции распределения, что привело к созданию квазиклассического корреляционного метода нссимметризованного самосогласованного поля (КНСП). Получаемая в итоге система кинетических уравнений в некотором смысле эквивалентна исходной системе гамильтоновых уравнений движения частиц, но каждая частица описывается не координатой и импульсом, а своей одночастичной функцией распределения по координате и импульсу. Такая система значительно более сложна, чем исходная система уравнений движения, но она может быть существенно упрощена различными методами. В частности при введении дальнего порядка для идеального кристалла система из N кинетических уравнений, описывающих эволюцию N одночастичных функций распределения, сводится к системе из К уравнений в периодических граничных условиях, где К — число частиц в элементарной кристаллической ячейке. Так же дальний порядок может вводится лишь по некоторым направлениям, что позволяет переходить к квазиодномерным и квазидвумерным моделям, описывающим, например, тонкие пленки.
Левая часть кинетических уравнений, получаемых из КНСП (фактически уравнений Власова), непосредственно следует из уравнения Лиувилля и является консервативной, для нее гамильтониан системы является иитегралом движения, то есть любая функция от гамильтониана является стационарным решением 1 [1] (стр. 17). В левую часть температура системы не входит. Одна из возможностей нахождения равновесного решения полученной системы — разложение потенциалов межчастичного взаимодействия в ряд Тэйлора в окрестностях минимума и использование равновесного распределения Больцмана, что позволяет в итоге перейти к системе алгебраических трансцендентных уравнений для моментов функций распределения, решаемых численно. Такой подход позволяет ввести в модель температуру и получить уравнения состояния для различных идеальных кристаллов [2].
В то же время решение нестационарных задач, за исключением самых тривиальных случаев может быть получено лишь численно. Для численного решения уравнений типа Фоккера-Плаика оптимальным представляется использование метода стохастического аналога [4, 12−15].
Сегнетоэлектрики в настоящее время изучены очень хорошо [16], и являются традиционным объектом для апробации различных моделей фазовых переходов второго рода. Многие модели сегнетоэлектриков использует явную феноменологическую зависимость потенциала от температуры Т типа / Т^с, оо и (х, Т) = а — + Ь—. v ' ' Тс 2 4.
В частности, такой вид имеет термодипамичский потенциал Ф в феноменологической теории фазовых переходов второго рода Ландау [17]. Но больший интерес представляют модели построенные из «первых принципов», в которых фазовый переход возникает в результате совместного действия самосогласованных полей и шума (температуры).
При условии, что гамильтониан не зависит явно от времени, например за счет внешних периодических сил.
В настоящее время в магнетиках [18−20] обнаружен ряд различных фазовых переходов, имеющих большое прикладное значение. Для их описания и анализа широко применяется численное моделирование. Обычно динамика магнетиков рассматривается на основе уравнений Ландау-Лифшица (или Гильберта) для магнитных моментов частиц системы [21, 22], стационарные состояния рассчитываются из различных энергетических соображений.
Система магнитных моментов как правило имеет несколько устойчивых состояний. Релаксация системы к тому или иному состоянию существенно зависит от неоднородностей (дефектов, границ образца), начальных условий, моделей релаксации и температуры (шума). В частности, именно шум приводит к спонтанному перемагничиванию (переходу между устойчивыми состояниями). Подавляющее большинство современных моделей магнетиков при описании динамики не учитывают воздействия шума (уравнение Ландау-Лифшица), либо не учитывают собственно динамику (процессы релаксации).
Описание динамики магнетиков па основе уравнений Ландау-Лифшица или Гильберта с учетом температуры при помощи источника Ланжевена уже проводилось в некоторых работах. В работе [23] шум вводится как дополнительное случайное слагаемое ?(х, г) к эффективному полю Случайное поле задается в виде дельта-коррелированной случайной гауссовой величины с интенсивностью, определяемой из флуктуационно-диссипаци-онной теоремы.
В работе [24] аналогичное слагаемое добавлялось к эффективному полю лишь в диссипативный член уравнения Гильберта.
Однако, такой модельный, феноменологический подход представляется не очень корректным. Свойства шума могут существенно влиять на результаты расчетов. Для корректного определения свойств шума представляется целесообразным использовать кинетические уравнения и метод стохастического аналога [4, 15].
В данной работе предложена система уравнений типа Фоккера-План-ка для эволюции функций распределения по магнитному моменту частиц системы, что позволяет корректно учитывать действие температуры, вводить понятия ближнего и дальнего порядка (в том числе по избранным направлениям внутри идеальных областей). Описание на уровне функций распределения является промежуточным между молекулярно-динамиче-ским и термодинамическим подходами и оперирует с необходимым минимумом информации для моделирования сильно неравновесных процессов в неоднородных системах.
Аналогичное уравнение Фоккера-Планка описывающее распределение магнитных моментов в сферических координатах для ансамбля монодоменных невзаимодействующих микрочастиц ферромагнетика (так называемая суперпарамагнитная система) было выведено Брауном [25] на основе уравнений Гильберта и Ланжевена. В дальнейшем это уравнение упрощалось и сводилось к системе уравнений для моментов фунцкции распределения (сферических гармоник) [26−30], а затем в некоторых частных случаях решалось численно, но его применение не выходило за рамки анализа различных систем монодоменных микрочастиц.
Так в работе [31] рассматривалось уравнение Фоккера-Планка для динамики намагниченности однодоменной частицы в тепловом резервуаре с температурой Т, но при решении уравнение было существенно упрощено, сведено к системе уравнений моментов проекции намагниченности на направление задаваемое внешнем полем, которые были разложены в ряд Фурье и решены численно методом матричной прогонки. Полученное решение было рассмотрено лишь в общем виде с точки зрения возникновения шумового резонанса.
В работе [32] были доказаны теоремы о существовании и единственности обобщенного и и классического решений задачи Коши для системы кинетических уравнений среднего поля типа Власова-Пуассона, описывающей модель твердого магнетика. Но в исследованной системе отсутствовала диссипация и не учитывалось влияния шума (уравнения выводились для магнитных моментов, прецессирующих в среднем поле).
Напыление тонких магнитных пленок является актуальной задачей. Тонкие магнитные пленки широко используются в различных устройствах хранения и преобразования информации. В настоящее время на эту тему имеется большое количество экспериментальных (напр. [33, 34]) работ.
Для моделирования напыления тонких пленок традицонно используются два подхода. На микроскопическом уровне применяют методы молекулярной динамики [35−38], но моделирование систем подвижных магнитных частиц при помощи методов молекулярной динамики является в настоящее время скорее некоторой экзотикой.
На значительно больших пространственных масштабах используют различные кинетические уравнения описывающие эволюцию функции распределения размеров зародышей или слоев. Подробный обзор этого подхода содержиться в [39]. Кроме того, на эту тему существует большое количество работ Г. И. Змиевской с соавторами [14,40−43], а так же других авторов [44−47].
1.2. Целью работы являются.
1. Разработка численных кинетических моделей, описывающих при помощи уравнений Фоккера-Планка фазовые переходы в сегнетоэлек-триках и магнетикахпроцесс напыления тонких магнитных пленок.
2. Разработка численных методов для решения используемых в моделях уравнений Фоккера-Планка и их реализация в виде соответствующих программ.
3. Кинетическое моделирование фазовых переходов в сегнетоэлектри-ках и магнетиках.
4. Кинетическое моделирование процесса напыления тонких магнитных пленок с учетом поступательного движения и эволюции магнитных моментов частиц.
1.3. В работе решены следующие задачи.
В рамках единого подхода (квазиклассического корреляционного метода несимметризованного самосогласованного поля) получены кинетические уравнения Фоккера-Планка, описывающие динамику сегнетоэлектриков, магнетиков и процесса напыления тонких магнитных пленок находящихся в термостате. Для описания отдельных частиц используются функции распределения по пространству, скорости, магнитному моменту, трактующиеся как результат усреднения по малым временам. Полученные системы уравнений легко упрощаются в случае предположения о дальнем порядке, позволяют эффективно описывать процессы с существенно различными характерными временами.
Показано, что феноменологическая теория фазовых переходов второго рода Ландау для сегнетоэлектриков и магнетиков может быть получена из предложенных уравнений Фоккера-Планка в случае термодинамического равновесия.
Предложенные системы уравнений решаются численно, на массиве траекторий, с использованием расщепления по физическим процессам и метода стохастического аналога. Для численного моделирования пленок так же использован метод, близкий к методу PIC (частица в ячейке). Численные схемы реализованы в виде высокопроизводительного кода на языках С++ и Python. Для управления расчетами разработан оригинальный пакет RACS (Results к Algorithms Control System — система контроля результатов и алгоритмов), берущий на себя всю черную работу по реализации интерфейса, упорядоченного хранения и поиска результатов, хранения исходного кода алгоритмов, воспроизведения расчетов и т. д.
Уравнения Фоккера-Планка получены для пространственно-однородных одноосных ионных сегнетоэлектриков из уравнений движения ионов вдоль «легкой» оси. Численно решена 1D1V2 задача. Исследован эффект динамического температурного гистерезиса. Уравнения Фоккера-Планка обобщены на случай адиабатически изолированной системы, проведены.
2 В дальнейшем для обозначения размерности задачи будут исполыцоваться выражения вида /DmVnS, где I— размерность конфигурационного пространства, m — размерность пространства скоростей (импульсов), п — размерность пространства магнитного момента. расчеты по адибатическому охлаждению внешним полем.
Для магнетиков уравнения Фокксра-Планка получены из системы уравнений Ландау-Лифшица, описывающей динамику системы магнитных моментов атомов кристаллической решетки.3 Численно решена ЗБ задача. Результаты кинетического моделирования существенно отличаются от результатов для системы Ландау-Лифшица:
• изменяются времена релаксации;
• появляются новые механизмы перемагничивания;
• учет шума при решении динамической задачи приводит к появлению новых эффектов, таких как динамический температурный гистерезис.
Приведены результаты расчетов для перемагничивания ансамбля невзаимодействующих магнитных моментов, динамики возникновения спонтанной поляризации при фазовом переходе парамагнетик-ферромагнетик, движения доменной границы в ферромагнетике.
Для напыления тонких магнитных пленок использовалась модельная зависимость обменного интеграла от расстояния. Численно решены Ю1/35 и ЮЗБ задачи. Исследована зависимость эффективности напыления пленки от величины обменного интеграла, температуры и других параметров.
Аналогичное уравнение Фоккера-Планка было выведено Брауном в 1963 г. из уравнения Гильберта при помощи метода Ланжевена [25], но только для ансамбля монодоменных ферромагнитных микрочастиц и использовалось различными исследователями только для изучения явлений суперпарамагнетизма.
1.4. Научная новизна.
Впервые построена кинетическая модель магнетиков на базе системы уравнений Фоккера-Планка, полученных при помощи КНСП, описывающих эволюцию функций распределения по магнитным моментам отдельных атомов.
Впервые построена кинетическая модель процесса напыления тонких магнитных пленок на базе системы уравнений Фоккера-Планка, описывающих эволюции функций распределения отдельных частиц по пространству, импульсу и магнитному моменту.
Показано, что для сегнетоэлектриков и магнетиков в случае термодинамического равновесия от построенной модели возможен переход к теории фазовых переходов второго рода Ландау.
Использованные в модели кинетические уравнения впервые решены численно методом стохастического аналога.
Для управления численными расчетами создан оригинальный пакет РАСБ (система контроля результатов и алгоритмов).
Численно обнаружена существенная зависимость эффективности напыления тонкой магнитной пленки от величины обменного интеграла.
1.5. Практическая ценность.
Разработанные модели могут быть использованы для изучения различных неравновесных процессов в твердом теле, для описания которых необходимо учитывать действие шума (температуры). В частности, результаты, получаемые из моделирования кинетических уравнений могут служить основой для получения эффективных уравнений состояния, необходимых в моделях сплошной среды.
Численное решение предлагаемых уравнений типа Фоккера-Планка может быть использовано для изучения различных резонансов (стохастического, шумового, силового), в том числе в магнетиках. В настоящее время на основе этих эффектов ведется разработка новых типов детекторов.
Кинетическое моделирование напыления тонких магнитных пленок может использоваться для оптимизации параметров технологического процесса.
1.6. На защиту выносятся следующие положения.
1. С использованием квазиклассического корреляционного метода несим-метризованного самосогласованного поля на основе кинетических уравнений типа Фоккера-Планка построены математические модели, описывающие фазовые переходы в сегнетоэлектриках и магнетиках, а так же процесс напыления тонких магнитных пленок. Показано, что теория фазовых переходов второго рода Ландау для сегнетоэлектриков и магнетиков может быть получена из уравнений Фоккера-Планка в случае термодинамического равновесия.
2. Для полученных систем уравнений разработаны численные схемы и алгоритмы, основанные на методах стохастического аналога и расщепления по физическим процессам, адаптированных к специфике задачи. Для численного моделирования пленок дополнительно использован метод аналогичный методу PIC (частица в ячейке). Численные схемы реализованы в виде высокопроизводительного программного комплекса на языках С++ и Python. Для управления расчетами и анализа результатов разработан оригинальный пакет RACS (система контроля результатов и алгоритмов).
3. Показана адекватность построенных моделей. В квазистационарных случаях результаты расчетов совпадают с известными аналитическими решениями. В сильно неравновесных случаях изменяются времена релаксации, появляются динамический температурный гистерезис и новые механизмы перемагничивания. Обнаружена существенная зависимость эффективности напыления от величины обменного интеграла при моделировании процесса напыления тонких магнитных пленок.
1.7. Личный вклад автора.
Уравнение Фоккера-Планка для сегнетоэлектриков было получено Ю. Л. Климонтовичем [3], но автор использовал несколько иное уравнение полученное из других соображений.
Численная схема для уравнения Фоккера-Планка описывающего сегне-тоэлектрики была разработана совместно с Г. И. Змиевской.
Остальные результаты работы получены самостоятельно.
1.8. Достоверность и обоснованность научных результатов.
Уравнения Фоккера-Планка получены из теоремы Лиувилля с использованием цепочки Боголюбова. Использовано приближение мультипликативности, традиционное для твердого тела. Добавленная феноменологически правая часть (фоккер-планковский «интеграл столкновений») необходима для обеспечения релаксации системы к равновесному решению Больц-мана.
Корректность полученных уравнений подтверждается для сегнетоэлек-триков и магнетиков возможностью перехода в равновесном случае к теории фазовых переходов второго рода Ландау.
Достоверность и обоснованность результатов численного моделирования достигается использованием метода расщепления по физическим процессам, метода стохастического аналога, схемы Рунге-Кутта и метода, аналогичного PIC (частица в ячейке), а так же многочисленными тестами и сравнениями с известными аналитическими решениями для различных частных случаев.
1.9. Апробация работы.
Результаты, представленные в диссертации, докладывались и обсуждались на Всероссийских и Международных научных конференциях: XXIVth Inetrnational Conference of Phenomena in Ionized Gases «ICPIG-26» (Варшава, 1999), «Звенигородская конференция по физике плазмы и УТС» (Звенигород, 1999, 2006, 2007), «Уравнения состояния вещества» (пос. Эльбрус, 2004, 2006), «Физика Экстремальных Состояний Вещества» (пос.Эльбрус, 2005), «Математические модели нелинейных возбуждений, переноса, динамики, управления в конденсированных системах и других средах» (Тверь, 1998; Москва, 2000), 25-th International Symposium of Rarefied Gas Dynamics «RGD25» (Санкт-Петербург, Репино, 2006). На 49-й научной конференции.
МФТИ «Современные проблемы фундаментальных и прикладных наук» (г.Долгопрудный, 2006).
1.10. Публикации.
Результаты, представленные в диссертации, опубликованы в работах: [48−57].
Заключение
.
В данной работе в рамках единого подхода (квазиклассического корреляционного метода несимметризованного самосогласованного поля) получены кинетические уравнения Фоккера-Планка, описывающие динамику сегнетоэлектриков, магнетиков и процесса напыления тонких магнитных пленок находящихся в термостате. Для описания отдельных частиц используются функции распределения по пространству, скорости, магнитному моменту, трактующиеся как результат усреднения по малым временам. Полученные системы уравнений легко упрощаются в случае предположения о дальнем порядке, позволяют эффективно описывать процессы с существенно различными характерными временами.
Предложенные системы уравнений решаются численно, на массиве траекторий, с использованием расщепления по физическим процессам и метода стохастического аналога. Для численного моделирования пленок так же использован метод, близкий к методу PIC (частица в ячейке). Численные схемы реализованы в виде высокнроизводительного кода на языках С++ и Python с использованием пакета RACS.
Показано, что феноменологическая теория фазовых переходов второго рода Ландау для сегнетоэлектриков и магнетиков может быть получена из предложенных уравнений Фоккера-Планка в случае термодинамического равновесия.
Проведено тестирование предложенных моделей для постановок имеющих аналитическое решение, проведено сравнение результатов традиционных подходов и кинетического моделирования.
Предлагемый подход можеть быть обобщен на случай систем со статистиками Феми и Бозе, что открывает возможность моделирования динамики различных квантовых эффектов с минимальными вычислительными затратами.
Хочется выразить глубокую благодарность Г. И. Змиевской за ценные замечания, многочисленные консультации и активное обсуждение различных аспектов данной работы на всем ее протяжении.