Метод интегралов по траекториям был предложен Фейнманом в 1948 году. Позднее была предложена идея использовать такой формализм для учета ядерных квантовых эффектов в методе молекулярной динамики, которые, как было показано в работе [1], дают значимый вклад даже при комнатной температуре.
Основная идея метода сводится к представлению ансамбля атомов в виде большого количества систем с разными микроконфигурациями, которые связаны между собой пружинками, образующими циклическую цепь. При таком подходе атом представляет собой так называемый кольцевой полимер, каждая частица которого имеет связи с другими атомами своей микроконфигурации (рис. 1 [2]).
Рис. 1. Взаимодействие двух частиц в соответствующих репликах в представлении подхода молекулярной динамики на интегралах по траекториям
В работе [3] было показано, что метод молекулярной динамики на интегралах по траекториям позволяет получить лучшую сходимость с экспериментальными данными, нежели классическая молекулярная динамика, для жидкой воды. Нас интересовал потенциал данного подхода в контексте задачи о стабильности различных протонных водных катионов в зависимости от температуры, рН, а также приложенного внешнего электрического поля.
Для проверки была выбрана модель воды tip4p/2005 [4], которая хорошо воспроизводит результаты (с небольшими поправками) для всей фазовой диаграммы воды. В качестве источника протонов была выбрана соляная кислота. Метод молекулярной динамики на интегралах по траекториям реализовывался посредством программных пакетов lammps [5] для расчета сил и i-pi [6] для расчётов координат каждого звена кольцевого полимера.
Для всех расчетов количество звеньев в цепи было равно шести. Вычисления производились в периодических граничных условиях. Частицы взаимодействовали посредством кулоновского потенциала и потенциала Леннарда — Джонса (6-12), который обрезался на расстоянии 9 Å. Для расчёта куловского взаимодействия использовалось суммирование Эвальда. Параметры потенциала 6-12 для ионов Cl– были взяты из статьи [7]. Параметризация такого взаимодействия для протонов с другими типами частиц описана в работе [8]. В свете отсутствия у протонов электронной оболочки для получения коэффициентов взаимодействия использовалось арифметическое правило смешивания [9].
На первом шаге решалась задача получения корректной плотности воды. Для этого 216 молекул H2O релаксировали при постоянных давлении и температуре (NpT-ансамбль) до достижения сходимости. Для расчета объема области моделирования были использованы экспериментальные данные [10] о плотности соляной кислоты. Далее в данной области вода релаксировалась, а затем равномерно разыгрывались координаты пятидесяти пар ионов Cl– и H+, что соответствует 13-молярному раствору соляной кислоты.
На рис. 2 представлены радиальные функции распределения (РФР) вновь внесенных частиц относительно остальных. Итоговые РФР получены посредством усреднения РФР, полученных от каждого звена кольцевого полимера.
Полученная система релаксировала для достижения равновесной плотности в NpT-ансамбле в течение 10 пикосекунд с временным шагом 0.25 фемтосекунд. Затем моделирование производилось при постоянных объеме и температуре (NVT-ансамбль) в течение 50 пс с тем же временным шагом. Все расчеты производились в нормальных условиях.
В результате были получены РФР тем же образом, что и ранее (рис. 3). Можно видеть наличие пика на длинах, соответствующих характерным расстояниям между соответствующими частицами в катионе Цунделя. Также был проведен анализ углов полученных катионов, что позволило судить о наличии в растворе не только катиона Цунделя, но и гироксония.
Кроме того, была проведена оценка количества и длины водородных связей, возникающих как в воде, так и в водном растворе соляной кислоты. С этой целью вода моделировалась в своем равновесном объеме в NVT-ансамбле. Среднее количество водородных связей на молекулу воды составило 3,50, для исследуемого раствора соляной кислоты: 1.99. «Обрезка» водородной связи проводилась по длине 3,5 Å и углу 120º.
На рис. 4 представлено изображения катиона Цунделя, полученное посредством моделирования выбранным методом. Красные шарики соответствуют атомам кислорода, белые — протонам.