CC BY
98
ВЕСТНИК ТГГПУ. 2010. №4(22)
УДК 538.91 + 532.74
ВЛИЯНИЕ СЕТКИ ВОДОРОДНЫХ СВЯЗЕЙ В ЖИДКОЙ ВОДЕ НА ЕЕ СТРУКТУРНЫЕ И ДИНАМИЧЕСКИЕ СВОЙСТВА
© Р.Р.Зарипов
В работе представлены результаты компьютерного моделирования молекулярной динамики воды в широком диапазоне температуры и давления. На основе результатов компьютерного моделирования проведено исследование структурных характеристик и сетки водородных связей. Показана связь аномалии вязкости воды с особенностями сетки водородных связей.
Ключевые слова: аномальные свойства воды, Т1Р5Р, водородные связи, молекулярная динамика.
Введение
Вода - ценнейший природный ресурс, играющий исключительную роль в функционировании различных физических, химических и биологических процессов. Все населяющие нашу планету живые организмы более чем на 50% состоят из воды, а сама Земля покрыта водой в твердой и жидкой фазах на 70%. Кроме того, вода является основой для большинства существующих технологий в сельском хозяйстве, химической и целлюлозно-бумажной промышленности, черной и цветной металлургии. При этом большая часть используемой в различных сферах человеческой деятельности воды, получив изрядную долю загрязняющих веществ, возвращается в реки в виде стоков.
Современные химические и биологические методы очистки сточных вод такие, как биофильтр и аэротенк, позволяют существенно уменьшить выбросы загрязнений. Однако рост городов, развитие промышленности, сельского хозяйства и связанная с этим возрастающая потребность в воде приводят к тому, что использование существующих технологий начинает требовать непомерных материальных затрат, делая их, таким образом, экономически невыгодными. Дефицит пресной воды хорошего качества уже сейчас становится проблемой мирового масштаба.
На качественно новый уровень решения данной проблемы позволит перейти чрезвычайно популярное на сегодняшний день направление науки - нанотехнологии. Однако здесь мы сталкиваемся с весьма неожиданными трудностями. Несмотря на свою распространенность и повсеместность, а также сотни тысяч научных трудов, посвященных исследованию воды, свойства этого простого вещества до сих пор непонятны до конца. Это связано, прежде всего, с известными аномалиями воды [1], по поводу которых у ученых всего мира по сей день не сложилось единого мнения.
К примеру, вода имеет аномально высокую температуру плавления и кипения. Плотность воды в интервале от 100 до 4°С нормально возрастает, как и у большинства других жидкостей. Однако, достигнув максимального значения при температуре 3,98°С, при дальнейшем охлаждении уменьшается, а при замерзании скачкообразно падает, тогда как почти у всех остальных веществ кристаллизация сопровождается увеличением плотности. Удельная теплоемкость, удельная теплота плавления и кипения воды аномально высоки по сравнению с другими веществами, причем удельная теплоемкость воды минимальна при 40°С. Вязкость воды с ростом давления уменьшается, а не повышается, как следовало бы ожидать по аналогии с другими жидкостями. Вода способна к значительному переохлаждению, т.е. может оставаться в жидком состоянии ниже температуры плавления.
Исследование воды методом молекулярной динамики
Существует большое количество разнообразных экспериментальных методов исследования воды. Среди них инфракрасная спектроскопия, дифракция рентгеновских лучей, дифракция нейтронов, спектроскопия протонного магнитного резонанса и др. Однако с помощью экспериментальных методов можно получить только усредненные значения характеристик воды. Изучить детали строения и характер движения молекул воды можно только при помощи компьютерного моделирования (численного эксперимента).
Впервые в качестве объекта молекулярного моделирования вода была использована в 1971 году в работах Рахмана и Стилленджера [2]. В последующем исследование воды стало одной из самых популярных задач численных экспериментов. При этом используют либо моделирование "из первых принципов" на основе законов квантовой механики, либо основанные на классической механике методы Монте-Карло и молекулярной динамики. С помощью первого метода
можно получить весьма высокоточные результаты, однако затраты больших вычислительных ресурсов, необходимые для определения межатомных сил и энергий, накладывают ограничения на размеры исследуемой системы. Метод Монте-Карло не позволяет получать последовательную траекторию движения молекул, затрудняя, таким образом, анализ динамических свойств системы. Метод молекулярной динамики лишен недостатков перечисленных методов, поэтому он и был взят за основу в нашем исследовании.
Суть данного метода состоит в численном решении уравнений движения частиц, энергия взаимодействия которых описывается некоторым заданным потенциалом. Для молекулярных систем потенциал взаимодействия, как правило, состоит из вкладов, которые отвечают как за меж-молекулярные взаимодействия, так и за взаимодействия атомов внутри одной молекулы. Выражение для полной энергии взаимодействия частиц в системе имеет следующий вид:
и = Ти« +Е иь
Ua = I К (Г - Г0 )2
иаЬ =I
qiqte
aibj
-4s0
I кв(0-в0)
angles
- V2 f ~
(1)
V Га0Ь0 J
Y
V ra0b0 J
Здесь а и Ь - индексы молекул воды, а / и ] -индексы атомов в молекуле. Вклад иа отвечает за поддержание валентных длин и углов между атомами водорода и кислорода (ОН) в молекуле воды. Вклад иаЬ представляет собой сумму потенциалов Леннарда-Джонса и Кулона, отвечающих за Ван-дер-Ваальсовые и электростатические взаимодействия соответственно. В уравнение (1) входят параметры потенциала Леннар-да-Джонса: £0 - глубина потенциальной ямы, а0 -среднее межатомарное расстояние; величина г0 определяет равновесную валентную длину, г -текущая длина связи, кг и кв - силовые константы; в0 и в - соответственно, равновесные и текущие значения углов.
Число частиц (центров взаимодействия) в компьютерной модели может быть больше, чем атомов в реальной молекуле воды. За счет увеличения числа центров взаимодействия добиваются большего согласия с экспериментальными данными. К примеру, в 1984 году Йоргенсен предложил четырехточечную модель Т1Р4Р [3]. Новый виртуальный центр взаимодействия был добавлен на биссектрисе валентного угла между атомами водорода. Отрицательный заряд кислорода был смещен в этот центр. Продолжением этой модели стала разработанная в 2000 году
Махоуни и Йоргенсеном пятиточечная модель Т1Р5Р [4]. Отрицательный заряд молекулы в данной модели был разделен на две части и расположен на углах тетраэдра. При этом в ван-дер-ваальсовых взаимодействиях участвовал электрически нейтральный центр тетраэдра. Данная модель значительно лучше всех прочих воспроизводит структурные характеристики и температурную зависимость плотности реальной воды.
Детали моделирования В настоящей работе приводятся результаты исследования структурных особенностей и сетки водородных связей воды в различных термодинамических состояниях, выполненного численно с помощью компьютерного моделирования молекулярной динамики.
Исследуемая система состояла из N = 4000 молекул воды, расположенных в кубической ячейке с периодическими граничными условиями. Взаимодействие частиц в системе осуществлялось на основе пятиточечного модельного потенциала Т1Р5Р [4]. Для численного интегрирования уравнений движения частиц в системе был использован алгоритм Бимона [5] с временным шагом 10-15 с. Все моделирование было проведено в КрТ-ансамбле с термостатом и баростатом Берендсена [6].
Структурные характеристики системы, а также конфигурация сетки водородных связей были получены при значениях внешнего давления от 1 до 5000 атмосфер и при температурах от 200 до 300 К.
Структура и динамика сетки водородных связей Как правило, аномалии физических свойств воды связывают со структурой ее молекулы, а также свойственным воде специфическим типом межмолекулярного взаимодействия - водородными связями (ВС). ВС возникают по причине неравномерного распределения электронной плотности в молекуле воды. Электронное облако молекулы стягивается в сторону более тяжелого атома кислорода, таким образом, молекула становится полярной.
В рамках численного эксперимента для идентификации частиц, связанных водородной связью используют специальные критерии ВС. Принципиально они делятся на два типа: энергетические и геометрические. В первом случае молекулы считаются связанными, если их суммарная потенциальная энергия меньше некоторого порогового значения. В случае геометрического критерия учитываются расстояния и углы между молекулами:
bonds
Т < ТУС
ОіО] — ОО ->
Т < Тс
оіи] —^ои?
(2)
структурных свойств используют парциальные радиальные функции распределения частиц [8]:
§аР{г ) = -
I3
І-
1зР
{г )
Значения геометрического критерия, как правило, подбирают исходя из первых минимумов соответствующих парциальных радиальных функций распределения.
Радиальная функция распределения является одной из основных характеристик структуры. В случае многокомпонентных систем для анализа
NaNp\ Л4ПГ "Аг/ (3)
Здесь gCф(г) - вероятность расположения двух атомов в диапазоне от г до г + Аг, Ь -длина ребра моделируемой ячейки. Na и Nр
число атомов типа а ив соответственно. Величина пв (г) определяет число в - частиц в
сферическом слое толщиной Аг на расстоянии г от у -той частицы.
і(А)
і'(А)
г{А)
Рис.1. Парциальные радиальные функции распределения частиц для различных значений температуры
при внешнем давлении р=1 атм.
Как видно на графиках рис.1, высота первого и второго пиков радиальной функции распределения имеет сильную зависимость от температуры исследуемой системы. С уменьшением температуры первый и второй пики для всех типов распределения становятся более выраженными.
В системе происходят структурные изменения, в результате которых увеличивается число частиц, находящихся в первой и второй координационных сферах относительно друг друга. Однако можно также заметить, что изменения температуры не влияют на положения минимумов и максимумов радиальной функции распределения, что позволяет использовать одни и те же значения геометрического критерия для исследования систем при различных температурах. Кроме того, разработанный нами алгоритм кластеризации позволяет учитывать тот факт, что одна молекула воды может одновременно участвовать не более чем в четырех водородных связях, выбирая на каждом шаге оптимальную конфигурацию ВС.
Рис.2. Трехмерная сетка водородных связей в жидкой воде
На рис.2 изображена визуализация конфигурации перколяционного кластера, полученного с помощью указанных критериев существования ВС. Все молекулы, представленные на рисунке, связаны единой сетью водородных связей, что говорит о полном соответствии исследуемой модельной системы современным представлениям о структуре воды на молекулярном уровне [9].
Однако данный перколяционный кластер является лишь мгновенной конфигурацией исследуемой системы и не позволяет судить о динамических особенностях сетки водородных связей. Основной характеристикой, с помощью которой можно оценить динамическое поведение сетки ВС, является ее время жизни (существования). Различают два типа времени жизни водородной связи [8] - непрерывное и мерцающее. Непрерывное время жизни длится до первого нарушения в данном димере выбранного критерия водородной связи. Если данная связь разорвалась, а потом восстановилась, то она рассматривается как новая вне зависимости от периода разрыва. Мерцающее время жизни водородной связи длится до тех пор, пока один из участников ВС не найдет себе партнера для образования новой связи. В этом случае промежуточные нарушения критерия существования связи не учитываются вне зависимости от их длительности.
В рамках данной работы для исследования времени жизни водородных связей была рассмотрена зависимость вида:
N (>)
(4)
N (0)'
где N(0) и N(1) - число водородных связей в начальный момент и в момент времени I соответственно [10]. Вновь образовавшиеся связи при этом не учитываются (непрерывное время жизни ВС). Понятно, что через определенный промежуток времени все существовавшие в начале связи разрушатся. Аппроксимируя функцию (4) экспонентой, можно найти среднее время жизни ВС:
( , >
* ехр
I
‘‘ЕЕ /
(5)
На основе данной зависимости было получено среднее значение времени жизни водородной связи (тЕЕ ) для всех исследуемых термодинамических состояний системы.
Т(К)
Рис.3. Температурная зависимость среднего времени жизни ВС
Среднее значение времени жизни водородной связи при "нормальных условиях" вполне соответствует аналогичным результатам численных экспериментов других авторов [11; 12]. Однако в большинстве изученных литературных источниках не уделяется никакого внимания зависимости тЕЕ от термодинамического состояния системы. Такая зависимость была исследована в рамках данной работы. На рис.3 видно, что при понижении температуры значение тЕЕ увеличивается на порядок. Кроме того, численный эксперимент показал, что тЕЕ имеет обратно пропорциональную зависимость от давления, т. е. среднее время жизни ВС с увеличением давления уменьшается. Эта особенность сетки водородных связей, вероятно, является причиной известной аномалии воды - вязкость этой жидкости уменьшается с увеличением давления.
Заключение
В рамках данной работы на основе компьютерного моделирования молекулярной динамики выполнено исследование равновесных структурных свойств жидкой воды, а также структурных и динамических особенностей сетки водородных связей в различных термодинамических состояниях. Установлено, что среднее время жизни водородной связи существенно зависит от термодинамического состояния системы. Эта величина пропорционально возрастает с увеличением температуры, однако имеет обратно пропорциональную зависимость от внешнего давления. При повышении давления наблюдается уменьшение среднего времени жизни ВС, что может говорить о связи аномалии вязкости воды с особенностями сетки водородных связей. Таким образом, с развитием компьютерной техники молекулярнодинамическое моделирование становится реальным инструментом как для научных исследова-
ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ
ний, так и для решения прикладных задач нанотехнологий.
1. Eisenberg D., Kauzmann W. The Structure and Properties of Water. - London: Oxford, 1969.
2. Rahman A., Stillinger F.H. Molecular Dynamics Study of Liquid Water // J. Chem. Phys. - 1971. -Vol.55. - P.3336-3359.
3. Jorgensen W.L. and Madura J.D. Temperature and Size Dependence for Monte Carlo Simulations of TIP4P Water // J. Mol. Phys. - 1985. - Vol.56. -P.1381.
4. Mahoney M. W. and Jorgensen W.L. A five-site model for liquid water and the reproduction of the density anomaly by rigid, nonpolarizable potential functions // J. Chem. Phys. - 2000. - Vol.112, 20. -P.8910-8922.
5. Beeman D. Some Multistep Methods for Use in Molecular Dynamics Calculations // J. Comp. Phys. -1976. - Vol.20, 2. - P.130-139.
6. Berendsen H.J.C., Postma J.P.M., et.al. Molecular Dynamics with Coupling to an External Bath // J. Chem. Phys. - 1984. - Vol.81. - P.3684-3690.
7. Rapaport D.C. Hydrogen bonds in water // J. Mol. Phys. - 1983. - Vol.50. - P.1151-1162.
8. Allen M.P., Tildesley D.J. Computer Simulation of Liquids. - Oxford: Clarendon Press, 1987. - 385 p.
9. Маленков Г.Г. Структура и динамика жидкой воды // Журнал структурной химии. - 2006. - №47 (Приложение). - S.5-35.
10. Sciortino F., Fornili S. Hydrogen bond cooperativity in simulated water: Time dependence analysis of pair interactions // J. Chem. Phys. - 1989. - Vol.90, 5. -Р.2786-2792.
11. Guardia E., Marti J., Padro J.A. ets. Dynamic in hydrogen bonded liquids: water and alcohols // J. Mol. Liq. - 2002. - Vol.96-97. - Р.3-17.
12. Swiatla-Wojcik D. Evaluation of the criteria of hydrogen bonding in highly associated liquids // J. Chem. Phys. - 2007. - Vol.342. - Р.260-266.
THE INFLUENCE OF THE HYDROGEN BOND NETWORK IN LIQUID WATER ON ITS STRUCTURAL AND DYNAMIC PROPERTIES
R.R.Zaripov
In this work we present the hydrogen bond network in water by molecular dynamics for liquid water in the wide region of temperatures and pressures. The structural properties and hydrogen bond network are studied. We have detected that anomaly of viscosity in liquid water is connected with the deformation of hydrogen-bond network.
Key words: anomaly properties of water, TIP5P, hydrogen bonds, molecular dynamics.
Зарипов Равиль Радикович - аспирант кафедры теоретической физики Татарского государственного гуманитарно-педагогического университета.
E-mail: [email protected]
Поступила в редакцию 09.11.2010
