Презентация

Цель: исследовать явление термодиффузии во флюиде Леннарда-Джонса методом молекулярной динамики, получить зависимости фактора термодиффузии бинарной системы от разности масс компонентов, плотности и температуры системы.
Чтобы это осуществить, нужно решить несколько задач:
1. Исследовать необходимую литературу о методе молекулярной динамики и алгоритмах, которые он использует
2. Исследовать литературу о явлении термодиффузии, методах расчета и экспериментального определения фактора и коэффициента термодиффузии
3. Разработать алгоритмы, позволяющие вычислить фактор термодиффузии бинарной сиситемы частиц в рамках метода молекулярной динамики
4. Проверить работоспособность программы молекулярной динамики на примере системы, находящейся в термодинамическом равновесии.
5. Провести вычислительный эксперимент с помощью программы молекулярной динамики по исследованию явления термодиффузии многочастичной системы, проанализировать полученные результаты.
Цели и задачи проекта

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

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

Потенциал межмолекулярного взаимодействия

ММД основывается на численном решении систем дифференциальных уравнений движения вида:
Для численного решения необходимо представить уравнение (2), с учетом равенства (3), в виде системы уравнений первого порядка.
Численное решение уравнений движения

Как можно увидеть из (3), основная вычислительная нагрузка приходится на расчет результирующей силы парного взаимодействия, следовательно алгоритмы для решения данной системы не должны производить более одного расчета правой части на итерацию. Для итерационной схемы верны следующие равенства:
Далее, разлагая (5) и (6) в ряд Тейлора, можно получить выражения для различных схем интегрирования, от простейших методов Эйлера и Эйлера-Кромера до алгоритма Верле и Бимана. Последние два алгоритма наиболее востребованы из за своей высокой точности.
Алгоритм Бимана обладает наибольшим порядком точности. Он может применяться в случаях, когда весомое значение имеют фазовые траектории отдельных частиц системы. Главным недостатком данной схемы является то, что она не обладает самостартуемостью, т.к. для вычисления состояния в последующий момент времени используются значения в предыдущий момент. Для решения данной проблемы первый шаг интегрирования обычно производят менее точными алгоритмами, например Эйлера-Кромера
Метод Верле, является компромиссом между точностью процедуры и скоростью её реализации. Этот метод является, пожалуй, одним из самых точных в классе методов интегрирования первого порядка

При проведении МД-расчетов для свободных нанокластеров при умеренной температуре не возникает необходимости в ограничении моделируемого пространства. Однако если система содержит жидкие или газообразные фракции, то со временем она будет расширяться во всех направлениях. Поэтому иногда приходится ограничивать моделируемую область.
Простейшим вариантом такой границы является ящик с отражающими стенками. Всякий раз, когда частица в процессе моделирования достигает внутренней поверхности ящика, ее перпендикулярная компонента скорости меняет знак. При таком граничном условии сохраняются энергия и объем системы, а также остается постоянным число частиц. Иными словами, имеется в виду изолированная система. В другом варианте отражение может зависеть от импульса частицы. Если перпендикулярная компонента скорости частицы превышает некоторое значение, она может покинуть ящик и, таким образом, исключается из дальнейшего рассмотрения. Есть ещё более оптимальный вариант: ящик заменяется сферическим потенциальным барьером, который частица может преодолеть, только обладая достаточной кинетической энергией.
Граничные условия

Граничные условия определяют термодинамическое окружение системы, поэтому они составляют существенную часть молекулярно-динамических моделей. Наиболее важным методом является применение периодических граничных условий. Данная идея возникла из проблемы, касающейся моделирования объемного кристалла или жидкости. Используемый прием состоит в том, чтобы ограничить моделируемое пространство кубическим ящиком, который частицы не могут покинуть. Но вместо полного отражения предполагается, что частица, вышедшая из ящика через одну из его граней, тут же возвращается в ящик через противоположную грань без изменения вектора скорости. Подобный подход обеспечивает постоянство энергии, объема и числа частиц. Однако учет взаимодействия между частицами производится с учетом новых граничных условий.
Периодические граничные условия позволяют рассматривать сравнительно небольшой «кубик» пространства, в котором расположена изучаемая молекула. Молекулы внутри кубика со временем совершают конформационные движения и перемещаются в пространстве, причем могут пересечь границы кубика. Суть метода заключается в том, что пространство разбивается на одинаковые кубики и предполагается, что их содержимое одинаково, а границы соприкасаются. При пересечении молекулой границы одного кубика она оказывается в другом, но при этом в первый кубик с противоположной стороны попадает такая же молекула. При таком подходе достаточно смоделировать динамику лишь одного такого кубика. Естественно, его размер должен быть достаточно большим, чтобы исключить возможность краевых эффектов.

В реальных экспериментах температура системы поддерживается за счёт энергообмена с внешней средой. Детальный учёт взаимодействия молекулы с внешней средой часто невозможен. Для учёта эффектов энергообмена с внешней средой используются специальные алгоритмы – термостаты.
Алгоритм (термостат) Берендсена корректирует производную температуры по времени согласно уравнению
т.е. отклонения температуры (T-T0) экспоненциально убывают c временной постоянной равной τT . Поскольку температура определяется кинетической энергией, которая пропорциональна квадрату скорости частиц, то корректируется не сама температура, а скорость частиц. Для этого после определения скоростей частиц их значения умножаются на фактор λ :
Термостат Берендсена

Этот метод имеет то преимущество, что интенсивность теплообмена может изменяться и адаптироваться под конкретные цели: для начального уравновешивания системы величина τT может быть довольно малой (например, 0,01 пс), но для надежного моделирования в равновесии она должна быть достаточно большой (например, 0,5 пс). Термостат Берендсена подавляет флуктуации кинетической энергии. Это значит, что он не обеспечивается правильный канонический ансамбль, поэтому, очевидно, выборка будет неточной. Погрешность уменьшается как 1/N, так что для очень больших систем большинство средних по ансамблю не будут существенно затронуты, за исключением собственно распределения кинетической энергии. Однако свойства, которые определяют величину флуктуаций (например, теплоемкость) будут затронуты.

Макроскопические свойства системы рассчитываются только после достижения системой стабильных статистических характеристик. Ниже приведены способы вычисления некоторых свойств.
Энергия. Потенциальная энергия системы ???? и кинетическая энергия системы ???? вычисляются как суммы, соответственно, потенциальных и кинетических энергий всех частиц:
где ? — число атомов в системе. В случае многоатомных жестких молекул к ???? добавляется кинетическая энергия вращательного движения
Полная энергия системы ? представляет собой сумму потенциальной и кинетической энергий
Расчет макроскопических свойств

Температура может быть вычислена через кинетическую энергию. Так как каждая вращательная степень свободы вносит вклад размера ???/2 в кинетическую энергию, то температура равна:
где ?? — константа Больцмана, ? — число степеней свободы. При трехмерном моделировании простых веществ, атомы которых взаимодействуют при помощи потенциала Леннард—Джонса, ? = 3.
Давление может быть вычислено при помощи уравнения состояния, если оно имеется для рассматриваемой системы, или при помощи вириальной теоремы по следующей формуле
где ? — объем области, ⟨. . .⟩ — обозначает осреднение по времени.
Плотность. Локальная плотность может быть вычислена при помощи разбиения области на ячейки и определения распределения частиц по ячейкам. Тогда локальная плотность в некоторой ячейке объемом ? , содержащей ? частиц массой ?, определяется по следующей формуле

В изолированной системе, содержащей смесь n компонентов с неоднородным первоначальным распределением концентраций, возникает процесс переноса массы компонентов, стремящийся к выравниванию концентраций. Процесс переноса может вызываться и другими физическими величинами – разностью температур и давлений.
Виды молекулярной диффузии (по причине возникновения):
1. концентрационная – за счёт разности концентраций;
2. бародиффузия – за счёт разности парциальных давлений;
3. термодиффузия – за счёт разности температур (температурный градиент).
Термодиффузия (термическая, или тепловая, диффузия) - перенос компонент газовых смесей или растворов под влиянием градиента температуры. Если разность температур поддерживается постоянной, то вследствие термодиффузии в объёме смеси возникает градиент концентрации, что вызывает также и обычную диффузию. В стационарных условиях при отсутствии потока вещества термодиффузия уравновешивается с обычной диффузией. Термодиффузия в газах была теоретически предсказана английским учёным С. Чепменом и шведским учёным Д. Энскогом (1911—17) на основе кинетической теории газов и экспериментально обнаружена английскими учёными С. Чепменом и Ф. Дутсоном в 1917.
Термодиффузия

В бинарной смеси при постоянном давлении в отсутствии внешних сил полный диффузионный поток вещества равен
ji = – nD12 gradci – n (DT/T) grad Т,
где D12 — коэффициент диффузии, DT — коэффициент термодиффузии, n — число частиц смеси в единице объёма, ci = ni/n — концентрация частиц i-й компоненты (i = 1,2).
с и (с-1) – концентрации компонентов смеси двух газов соответственно
Это приближенное аналитическое выражение для фактора термодиффузии лишь для случая идеального газа. Величина αT идеального газа зависит от разности масс компонент и показателя степени n, отвечающего за отталкивание в выражении для потенциала взаимодействия