НАНОСИСТЕМЫ
DOI: 10.17725/rensit.2021.13.141
Корреляция поступательного и вращательного движения молекул воды в молекулярно-динамических моделях Волошин В.П.
Институт химической кинетики и горения им. В.В. Воеводского, Сибирское отделение РАН,
http://www.kinetics.nsc.ru/
Новосибирск 630090, Российская Федерация
Е-mail: [email protected]
Поступила 30.04.2021,рецензирована 20.05.2021, принята 1.05.2021
Аннотация: Движение жёсткой молекулы в компьютерной модели можно описать как сумму перемещения центра масс молекулы и вращения молекулы вокруг оси, проходящей через этот центр. Взаимодействие молекулы со своим окружением приводит к возникновению корреляции между этими движениями. Ранее в [1] нами были изучены распределения углов между направлениями векторов перемещения и осей вращения молекул воды, а также между ними и внутренними векторами этих молекул. В данной работе описывается корреляция величин этих перемещений, то есть длин векторов перемещений и углов сопутствующих вращений. Рассчитаны коэффициенты корреляции этих характеристик на интервалах различной продолжительности при разных температурах и давлениях, определены характерные времена сохранения этих корреляций. Мы полагаем, что найденные нами времена представляют собой времена жизни локальных окружений молекул. Для молекул воды изменение локального окружения сопровождается изменением структуры ближайшего участка сетки водородных связей, а потому времена жизни локальных окружений близки к временам жизни этих связей.
Ключевые слова: молекулярно-динамическое моделирование, вода, пространственное смещение, угловое смещение, корреляция перемещения и вращения
УДК 532.74
Для цитирования: Волошин В.П. Корреляция поступательного и вращательного движения молекул воды в молекулярно-динамических моделях. РЭНСИТ, 2021, 13(2):141-148. DOI: 10.17725/ rensit.2021.13.141._
Correlation of Translational and Rotational Motion of Water Molecules in Molecular Dynamic Models
Vladimir P. Voloshin
Voevodsky Institute of Chemical Kinetics and Combustion SB RAS, http://www.kinetics.nsc.ru/ Novosibirsk 630090, Russian Federation Е-mail: [email protected]
Received 27 April, 2021, peer-reviewed 10 May, 2021, accepted 24 May, 2021
Abstract: The motion of a rigid molecule in a computer model can be described as the sum of the displacement of the center of mass of the molecule and the rotation of the molecule around an axis passing through this center. The interaction of a molecule with its environment leads to a correlation between these movements. Earlier in [1], we studied the distributions of angles between the directions of the displacement vectors and the axes of rotation of water molecules, as well as between them and the internal vectors of these molecules. This paper describes the correlation of the magnitudes of these displacements, that is, the lengths of the displacement vectors and angles of accompanying rotations. The correlation coefficients of these characteristics are calculated for intervals of different durations at different temperatures and pressures, and the characteristic
ВОЛОШИН В.П.
НАНОСИСТЕМЫ
times of preservation of these correlations are determined. We believe that the times found by us represent the lifetimes of the local environments of molecules. For water molecules, a change in the local environment is accompanied by a change in the structure of the nearest section of the network of hydrogen bonds, and therefore the lifetimes of local environments are close to the lifetimes of these bonds.
Keywords: molecular dynamics simulation, water, translational displacement, angular displacement, correlation of displacement and rotation UDC 532.74
For citation: Vladimir P. Voloshin. Correlation of Translational and Rotational Motion of Water Molecules in Molecular Dynamic Models. RENSIT, 2021, 13(2):141-148. DOI: 10.17725/rensit.2021.13.141._
Содержание
1. 2. 3.
5.
6.
Введение (142) модели (143)
поступательные и угловые перемещения (143)
корреляция поступательных и угловых перемещений (145)
Времена жизни водородных связей (147) заключение (147)
Литература (148) 1. ВВЕДЕНИЕ
Движение любой молекулы можно разложить на два формально независимых движения — поступательное и вращательное. В качестве поступательного можно рассмотреть
перемещение центра масс молекулы. Тогда оставшееся движение будет представлять собой вращение молекулы вокруг оси, проходящей через этот центр. Если молекула жёсткая, то есть изменением расстояний между её атомами можно пренебречь, такое описание будет однозначно предсказывать положение и движение каждого атома молекулы.
Результатом молекулярно-динамического моделирования являются координаты центра каждого атома модели в любой момент времени относительно исходно заданной системы координат. Положение атома относительно начала координат можно изобразить с помощью радиус-вектора, проведённого из этого начала к текущему положению центра атома. Перемещение атома описывается с помощью вектора перемещения, равного разнице радиус-векторов его конечного и начального положений. То есть, и положение, и перемещение атома описывается в одной и той же форме — в виде вектора, проведённого
к текущему положению атома либо из начала координат, либо из начального положения на интервале движения.
Мы предлагаем свести описание вращения молекул к столь же простой схеме. Для этого для молекул каждого вида определим собственную систему координат, поместив её начало в центр масс молекулы, и выбрав направления осей по положению её атомов. Ориентацию молекулы, в которой направления осей собственной системы совпадают с направлениями осей системы координат модельного бокса, будем называть базовой или нулевой. Любую текущую ориентацию опишем параметрами вращения, переводящего молекулу из базовой ориентации в текущую. Таким образом, мы зададим пространство ориентаций выбранного вида молекул, в котором роль начала координат будет играть базовая ориентация этих молекул, а параметры вращения, переводящего её в текущие ориентации — координатами этих ориентаций. Вращение, переводящее молекулу из некоторой начальной ориентации в конечную, будут сопровождаться изменением координат в этом пространстве, то есть представлять собой перемещение в нём. То есть, при таком подходе и текущие ориентации, и вращения между ними будут записываться в одной и той же форме — параметрами вращения.
Такой подход можно реализовать для любого способа описания вращений, например, с помощью углов Эйлера или матриц вращения. Однако наиболее удобной для этого оказалась алгебра кватернионов. Кватернион представляет собой гиперкомплексное число с одной действительной и тремя разными мнимыми компонентами. Каждому вращению на угол ф вокруг единичной оси и соответствует
НАНОСИСТЕМЫ
кватернион, действительная часть которого равна cos(ф/2), а три мнимые соответствуют вектору направленному вдоль
оси вращения. То есть, в отличие от других способов описания вращения, кватернион содержит параметры вращения практически в явном виде. В Приложении к работе [1] приведено детальное изложение процедур по использованию алгебры кватернионов для описания текущих ориентаций, вращений между ними и преобразованию векторов при переходе от внутренней системы координат молекулы к системе координат модели и обратно.
Поступательное и вращательное движения являются независимыми только у изолированной молекулы. Взаимодействие молекулы с окружением приводит к возникновению корреляции между этими видами движений. В работе [1] мы показали, что в модели жидкой воды наиболее часто вектор перемещения перпендикулярен оси вращения. Одновременно с этим вектор перемещения часто параллелен вектору нормали, а ось вращения параллельна отрезку, соединяющему центры атомов водорода данной молекулы. Эти закономерности соответствуют движению, при котором молекула воды раскачивается на двух наиболее сильных водородныхсвязях, обычно таких,вкоторыхданная молекула участвует как донор протона. Впрочем, все распределения достаточно широкие, и потому представленная картина описывает лишь наиболее вероятные виды движения молекул воды, и не охватывает всё их разнообразие. В данной работе мы изучаем корреляции не направлений, но величин пространственных перемещений и углов вращения молекул, а также определяем характерные времена сохранения этих корреляций.
2. МОДЕЛИ
В работе использованы модели воды, приготовленные с помощью пакета молекулярно-динамического моделирования LAMMPS
[2]. Каждая модель содержит 8000 молекул с потенциалом взаимодействия Т1Р4Р/2005
[3] в кубическом боксе с периодическими граничными условиями, шаг моделирования 2 фс. При давлении 1 бар были приготовлены модели с температурой 260, 280, 300, 330 и 360
К. Кроме того, при температуре 300 К были рассчитаны модели с давлением 2, 5, 10 и 15 Кбар. Предварительная релаксация длилась не менее 1 нс. Для анализа перемещений и вращений использовались по 20 000 мгновенных конфигураций через каждые 200 фс. Для каждой молекулы каждой записанной конфигурации были определены координаты центров масс молекул, а также кватернионы их текущих ориентаций [1]. Расчёт пространственных и угловых перемещений производился для интервалов времени 0.2, 0.4, 0.6, 1, 2, 4 и 6 пс. Для каждого интервала в качестве начальной и конечной конфигураций использовались все возможные пары, разделённые таким интервалом; таким образом, для наиболее короткого интервала было использовано по 20 000 наборов смещений для каждой из 8000 молекул, а для наиболее длинного по 19 970 наборов.
Для вычисления автокорреляционных функций поступательных и вращательных скоростей, а также для определения времён жизни водородных связей, дополнительно использовалось по 60 000 мгновенных конфигурации для тех же моделей, но на этот раз на каждом шагу моделирования, то есть через 2 фс.
3. ПОСТУПАТЕЛЬНЫЕ И УГЛОВЫЕ ПЕРЕМЕЩЕНИЯ
Рассмотрим вначале поступательные и вращательные смещения по отдельности. На Рис. 1 представлены распределения поступательных (а) и вращательных (Ь) смещений молекул воды за разные интервалы времени в модели при температуре 300 К и давлении 1 бар. На врезках показаны зависимости средних квадратов соответствующих смещений от длительности движения. Для поступательных смещений эта зависимость уже после 1 пс выглядит как прямая линия в соответствии с формулой Эйнштейна для трёхмерного пространства. Вращательные смещения можно рассматривать как смещения на поверхности сферы. В модели случайных независимых вращений при малых углах полного поворота, где поверхность сферы практически плоская (косинус угла очень близок к 1), средний квадрат
ВОЛОШИН В.П.
НАНОСИСТЕМЫ
Рис. 1. Распределения поступательных (а) иугловых (Ь) перемещений молекул воды при P = 1 бар, T = 300 заразные промежутки времени (указаны около кривых). На врезках — зависимости средних квадратов соответствующих
смещений от времени.
конфигурациях.
На Рис. 2 приведены автокорреляционные функции поступательных (пунктир) и вращательных (сплошная линия) скоростей молекул воды, нормированные на величину среднего квадрата соответствующих скоростей. Нормированные автокорреляционные
функции начинаются со значения 1 (идеальная корреляция начальных скоростей самих с собой) и вначале монотонно уменьшаются, поскольку со временем модули и направления скоростей молекул постепенно изменяются. Положение первого минимума отмечает время, при котором
угла поворота также пропорционален времени согласно формуле Эйнштейна для случайного движения на плоскости. Однако даже в такой модели по мере увеличения угла наклон зависимости постепенно понижается, вначале из-за нарастающего отклонения от плоскости, а затем из-за ограничения достоверно определяемого угла поворота значением п. В модели воды помимо этого существует корреляция между последовательными молекулярно-динамическими шагами, а потому даже на коротких временах зависимость среднего квадрата угла поворота отклоняется от прямой, аналогично тому, как это проявляется на зависимости среднего квадрата длины перемещения центра масс. В результате зависимость среднего квадрата угла поворота в модели воды (Рис. 1^) нигде не описывается прямой линией.
Корреляция между скоростями молекул в последовательных шагах моделирования может быть показана с помощью автокорреляционной функции. Для расчёта данной функции мы использовали в качестве начальных 50 000 конфигураций из набора конфигураций через 2 фс, и вычисляли средние скалярные произведения скоростей в начальных и последовательных конфигурациях. Для поступательного движения использовали скорости атомов кислородов, а угловые скорости извлекали из кватернионов вращения, рассчитанных сравнением кватернионов ориентаций молекул в начальной и конечной
Рис. 2.
поступательных и вращательных скоростей молекул воды при Р = 1 бар и Т = 300 К Стрелками показаны положения вторых максимумов, соответствующих периодам основных трансляционных и либрационных колебаний.
НАНОСИСТЕМЫ
наибольшее количество молекул изменило свои скорости на противоположные. Следующий за ним второй максимум представляет время, при котором наибольшее количество молекул оказалось вновь движущимся в первоначальном направлении. Положения этих максимумов (показаны вертикальными стрелками)
соответствуют периодам основных колебаний. Для трансляционных колебаний период таких колебаний составляет примерно 130 фс (от 117.1 фс в модели с наибольшим давлением до 134.8 фс в модели с наибольшей температурой), для либрационных — около 47 фс (от 45.3 фс до 49.3 фс в тех же моделях). Это означает, что в течение одного трансляционного колебания молекула воды успевает совершить несколько либрационных колебаний. Сложная форма автокорреляционной функции угловых скоростей наводит на мысль, что молекула воды участвует одновременно в нескольких либрационных колебаниях с близкими частотами.
4. КОРРЕЛЯЦИЯ ПОСТУПАТЕЛЬНЫХ И УГЛОВЫХ ПЕРЕМЕЩЕНИЙ
Корреляцию между модулями поступательных и вращательных перемещений можно изобразить с помощью двухмерной диаграммы, показанной на Рис. 3. Здесь использовались смещения молекул воды для интервала 0.6 пс в той же модели, что на предыдущих рисунках. Данная диаграмма показывает, что корреляция между этими модулями явно присутствует, и она
ф И
СИ
<L> Т>
30 80-
70
Ф
Е m и Ч CL и
30
60-1
50
40
2010
о Сс
Т1Р4РГ2005 P = 1 bar
T = 30ÛK
Ss
¡Ai)
w
it = 0 6 ps
положительная, хотя довольно слабая.
Численно корреляцию между
произвольными характеристиками х и у одной и той же системы из п объектов можно описать с помощью коэффициента корреляции Пирсона:
(X -(*)) • (у -( у) )
xy
=1
Vi >-< x) ,:,(y, -< y )2
(1)
где (х) = (Е "х )/п - среднее значение характеристики х. Ошибка вычисления этого коэффициента определяется по формуле
ту =7 (1 - гХУ )/(п - 2).
В наших расчётах мы использовали идентичную, но более удобную формулу:
хУ) -(х) •( У)_ , (2)
r =
xy
< ■ <
= j (x2) -( x)2.
о
00 05 1.0 1.5 20 Transi atiûnal dis placement, А
Рис. 3. Совместное распределение длины перемещения центра масс молекулы воды (по горизонтали) и угла её поворота (по вфтикали) за время движения 0.6 пс.
где \хУ) = (Е"=1 Х- • У')/п' ах ,
^х= (Ехг2)/п. Коэффициент корреляции принимает значения от максимального +1 (идеальная корреляция, когда характеристики растут и уменьшаются синхронно и в равной степени) до минимального -1 (идеальная антикорреляция, когда ростодной характеристики сопровождается синхронным понижением другой характеристики). Для полностью независимых характеристик коэффициент равен 0. Коэффициент корреляций смещений, изображённых на Рис. 3, равен гу = 0.303, а ошибка его определения т = 0.00024. То есть, между величинами поступательных и вращательных смещений в данной модели достоверно присутствует слабая положительная корреляция.
Величина коэффициента корреляции поступательных и вращательных смещений молекул воды зависит от термодинамических условий модели, а также от длительности движения Д/. Результаты представлены на Рис. 4 для моделей с разными температурами (а), и с разными давлениями (В). Повышение температуры или давления при любой фиксированной длительности движения уменьшает значение коэффициента
корреляции.
При изменении длительности движения коэффициент корреляции в любой модели
1=1
ВОЛОШИН В.П.
НАНОСИСТЕМЫ
04
с
и 0 3
ЗГ:
01 п
о с о г
о
73
£ 1_ 0.1 ■
Г)
О
00
а ^-гп* **"*■ 2ИЖ
' 280К
»ЭКЖ
Т1Р4Р/2006 \ ' Э30К
Р = 1 бар
» эеок
0 25 0.5 1 2 4 М, рз
0.4
с ш
3 0.3 Е ш
|о,
О
¡0.1
о
О
00
1 Ьзг ь
2 КЬаг
10 кьаг
16 КЬаг
Т1Р4Р/2005
Т = 300 К
0.125 0.25 0.5 1 2 4 р5
Рис. 4. Зависимость коэффициента корреляции поступательных и вращательных смещений от длительности движения Дt для моделей с разной температурой (а) и разным давлением (Ь).
проходит через максимум. Длительность, при котором коэффициент принимает максимальное значение, очень сильно зависит от температуры (от 0.5 до 2.7 пс) и почти не зависит от давления (от 0.6 до 0.9 пс) (см. также Рис. 5). Эти времена существенно больше, чем периоды трансляционных колебаний (порядка 130 фс) и тем более либрационных колебаний (порядка 50 фс), которые мы определили по автокорреляционным функциям скоростей.
Таким образом, за время сохранения корреляции между поступательным и вращательным перемещением молекула воды успевает совершить множество трансляционных колебаний, и ещё больше либрационных. Сохранение корреляции не взирая на колебания означает, что все они происходят в неизменном локальном
окружении. Рост коэффициента корреляции на малых временах связан с нивелированием случайных отклонений текущей формы этого окружения от его средней равновесной формы. Понижение коэффициента при больших временах происходит из-за роста числа молекул, окружение у которых изменилось. В результате всё чаще начальные и конечные положения оказываются принадлежащими движениям в разных условиях, между которыми корреляция практически отсутствует. Таким образом, длительность интервала времени, при которой корреляция поступательного и вращательного движения проявляется в наибольшей степени, является характерным временем сохранения локального окружения молекул воды в почти неизменном виде, то есть временем жизни этого локального окружения.
Рис. 5. Сравнение характерных времён жизни локальных окружений молекул воды ДШах (я) со средними временами жизни водородных связей, определённых согласно чистому геометрическому критерию (о), а также
комбинации геометрического и динамического критерия (▲)
НАНОСИСТЕМЫ
5. ВРЕМЕНА ЖИЗНИ ВОДОРОДНЫХ СВЯЗЕЙ
Для молекулы воды смена локального окружения сопровождается разрывом каких-либо старых водородных связей и заменой их на новые. На Рис. 5 произведено сравнение времён жизни локальных окружений А^ и средних времён жизни водородных связей, рассчитанных согласно двум разным критериям: на Рис. 5а в зависимости от температуры, на Рис. 5В от давления. Времена жизни локальных окружений А^ показаны чёрными квадратиками. Пустые окружности представляют времена жизни водородных связей, определённых согласно простому геометрическому критерию, часто используемому для их расчёта. Связь между молекулами считается существующей, если расстояние между центрами их кислородов Яоо < 3.3 А, а минимальное расстояние между кислородом одной молекулы и одним из водородом другой &он < 2.45 А. Ранее в работах [4-5] мы отмечали, что при использовании подобного критерия существует множество связей, претерпевающих кратковременные разрывы с последующим восстановлением. Поскольку восстановление водородной связи означает, что данная пара молекул осталась в окружении друг друга, мы посчитали, что наиболее кратковременные разрывы можно игнорировать. Это соображение мы реализовали во втором критерии в виде динамической добавки к первому: связи, возникшие согласно первому критерию, считались по второму критерию разорванными только в том случае, если длительность разрыва превышала половину периода основных трансляционных колебаний, определённых по положению второго максимума автокорреляционной функции скорости, как это было сделано на Рис. 4. В разных моделях длительность игнорируемого разрыва составляла от 58 до 68 фс. Времена жизни связей, определённые согласно этому комбинированному критерию, на Рис. 5 показаны треугольниками.
Как видно из Рис. 5, времена жизни по геометрическому критерию оказались в 1.5-2 раза меньше, чем времена жизни локального окружения. Однако времена жизни по комбинированному критерию соответствуют
этому времени гораздо лучше. Наибольшие отклонения наблюдаются при наименьшей температуре (260 К) и большом давлении (от 5 Кбар и выше). Возможно, эти отклонения связаны с тем, что при этих условиях молекулы воды согласно комбинированному критерию имеют самое большое среднее количество водородных связей (около 4 при Т = 360 К и около 4.5 при Р = 15 Кбар). То есть, для изменения локального окружения требуется разорвать больше связей, причём одновременно, поскольку иначе ранее разорванная связь может восстановиться и продлить продолжение существования текущего окружения.
6. ЗАКЛЮЧЕНИЕ
Взаимодействие молекул воды с окружением приводит к возникновению корреляции между поступательным и вращательным движением. Коэффициент корреляции Пирсона между длиной пространственного перемещения молекулы воды и углом её поворота за это же время не превышает 0.4, что означает присутствие слабой положительной корреляции между этими характеристиками. С ростом температуры и давления коэффициент уменьшается. Однако его величина зависит также от длительности движения. Достигнув максимума при движении в течение примерно пикосекунды, далее коэффициент корреляции монотонно падает. Мы полагаем, что уменьшение корреляции происходит после изменения локального окружения молекулы, и потому предлагаем считать время достижения максимального значения коэффициента корреляции средним временем жизни локальных окружений. Полученное время существенно превышает периоды основных трансляционных и либрационных колебаний, а значит за время жизни локального окружения молекулы успевает совершить множество таких колебаний.
В рамках картины теплового движения молекул в жидкости, предложенной Я.И. Френкелем, каждая молекула жидкости заметную часть времени колеблется внутри «клетки», сформированной молекулами локального окружения. Выход из клетки происходит только после её разрушения. Возможно, найденное
ВОЛОШИН В.П.
НАНОСИСТЕМЫ
нами время представляет собой время оседлой жизни молекулы воды, и мы реализовали метод по его определению. В работе [6] было высказано мнение, что именно наличие водородных связей приводит к столь существенному различию в поведении молекулы воды во время оседлой жизни внутри клетки и во время переходов между клетками. Поэтому неудивительно, что найденное нами время жизни локального окружения в воде близко ко времени жизни водородных связей в ней. Более того, эти времена демонстрируют практически одинаковые зависимости от температуры и давления.
ЛИТЕРАТУРА
1. Волошин ВП, Наберухин ЮИ. Описание вращательных движений молекул в компьютерных моделях воды с помощью кватернионов. РЭНСИТ, 2020, 12(1):69-80; DOL10.17725/rensit.2020.12.069.
2. Plimpton S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comp. Phys., 1995, 117:1-19. URL: http://lammps.sandia.gov.
3. Abascala JLF, Vega C. A general purpose model for the condensed phases of water: TIP4P/2005. J. Chem. Phys, 2005, 123:234505. DOI: 10.106з/1.2121687.
4. Волошин ВП, Желиговская ЕА, Маленков ГГ, Наберухин ЮИ, Тытик ДЛ. Структуры сеток водородных связей и динамика молекул воды в конденсированных водных системах. Росс. хим. ж., 2001, XLV(3):31-37.
5. Волошин ВП, Наберухин ЮИ. Распределение времени жизни водородных связей в компьютерных моделях воды. Журн. структ. химии, 2009, 50(1):84-95.
6. Родникова МН. Особенности растворителей с пространственной сеткой H-связей. Журн. физ. химии,, 1993, 67(2):275-280.
Волошин Владимир Петрович
к.ф.-м.н, с.н.с.
Институт химической кинетики и горения им. В.В. Воеводского СО РАН
3, ул. Институтская, Новосибирск 630090, Россия [email protected].