CC BY
82
100
Химия
Вестник Нижегородского университета им. Н.И. Лобачевского, 2013, № 2 (1), с. 100-106
УДК 544.272+544.034.1
ДИНАМИКА И ДИФФУЗИЯ АТОМОВ И РАДИКАЛОВ, ОБРАЗУЮЩИХСЯ ПРИ УФ-ОБЛУЧЕНИИ ВОДНОГО ЛЬДА.
МОЛЕКУЛЯРНО-ДИНАМИЧЕСКОЕ ИССЛЕДОВАНИЕ
© 2013 г. Н.В. Свинков1, С.К. Игнатов\ М.Ю. Куликов2, А.Г. Разуваев1
1 Нижегородский госуниверситет им. Н.И. Лобачевского 2 Институт прикладной физики РАН, Н. Новгород
n.v. svinkov@gmail. com
Поступила в редакцию 10.10.2012
Методом классической молекулярной динамики проведено исследование поведения молекул воды, атомов водорода и гидроксил-радикалов внутри кристаллического водного льда в интервале температур от 50 до 150 K. Определен характер движения частиц, рассчитаны коэффициенты диффузии и са-модиффузии (для молекул воды) в указанном интервале температур.
Ключевые слова: водный лед, моделирование, классическая молекулярная динамика, процессы фо-
толиза, радикалы, диффузия.
Введение
Водный лед в кристаллическом и аморфном состояниях присутствует в больших количествах в верхних слоях атмосферы нашей планеты [1], а также межзвездном пространстве [2, 3]. Изучение воздействия ультрафиолетового излучения на водный лед имеет большое значение для химии атмосферного и межзвездного льда. Под воздействием жесткого космического УФ-излучения во льду интенсивно протекают процессы фотолиза, первичными продуктами которого являются атомы водорода и радикалы ОН. Скорость и механизм диффузии образовавшихся атомов и радикалов внутри льда и на его поверхности определяют вероятность процессов их дальнейшей рекомбинации, а также протекания реакций с другими атомами и молекулами. Установлено [4-7], что главным продуктом фотолиза, накапливающимся внутри льда, является пероксид водорода, однако его выход и скорость накопления зависят от температуры и состава внешней среды. При поглощении больших доз УФ-фотонов а-серии Лаймана в интервале температур 35-100 К обнаружена фотодесорбция молекул воды [6, 7] и других продуктов, таких как Н2 и О2 [6]. Десорбция радикалов может приводить к «испарению» льда под воздействием УФ-излучения и будет влиять на химический состав льда [8]. В целом механизмы химических и физических процессов, следующие за стадией поглощения УФ-излучения, изучены слабо. В частности, скорость и механизм диффузии атомов водорода и радикалов ОН после первичной фотодиссоциации молекулы воды, их динамика внутри твердого льда и кон-
станты скорости элементарных реакций с окружающими молекулами остаются неизвестными.
Процессам фотодиссоциации молекул воды посвящен ряд теоретических работ. Процессы десорбции и рекомбинации продуктов фотолиза кристаллического и аморфного водного льда исследованы Андерсоном с сотрудниками методом классической молекулярной динамики (МД) [810]. В работах представлены данные в интервале температур от 10 до 90 К о величинах свободного пробега атомов Н и радикалов ОН, вероятности рекомбинации фотофрагментов и фотодесорбции радикалов и молекул воды с поверхности ледяных частиц. Авторы указывают на большую подвижность и склонность к фотодесорбции атомов Н, чем ОН. В работе [11] представлены результаты МД исследования диффузии атома водорода внутри гексагонального льда при температуре 8 К, рассчитаны коэффициенты диффузии Н при данной температуре.
С практической точки зрения исследование процессов образования и перемещения продуктов фотолиза имеет важное значение для изучения баланса перекисных соединений в атмосфере, влияния на климат и атмосферной химии в целом. Однако большинство работ, связанных с моделированием диффузии фотофрагментов, посвящено низкотемпературной области (~10 К), что соответствует условиям космического пространства. Существует необходимость детального теоретического исследования динамики и диффузии фотофрагментов в области температур 10-250 К, при которых ледяные частицы существуют в различных слоях атмосферы Земли. Данная работа посвящена изучению особенностей диффузии продуктов фотолиза
водного льда методом молекулярной динамики в интервале температур 50-150 К.
Методика моделирования
Твердая фаза воды, или водный лед, может существовать в различных кристаллических и аморфных разновидностях. При обычных (низких) давлениях стабильную кристаллическую разновидность льда называют льдом I. Он представляет собой совокупность двух кристаллических фаз - фазы Ih (обычный, или гексагональный, лед), обладающей гексагональной кристаллической решеткой, и фазы Ic (кубический лед), обладающей кубической решеткой. Лед Ih, или обычный лед, является наиболее распространенной формой льда на Земле. Параметры решетки вблизи точки плавления составляют a = 0.4532 нм, с = 0.7367 нм. Отношение c/a (1.628) очень близко к идеальному отношению гексагональной упаковки металлов и практически не зависит от температуры. Таким образом, решетка напоминает искаженную структуру графита.
Особенностью кристалла льда является подвижность структуры его поверхности. Поскольку водородные связи, соединяющие молекулы воды, имеют энергию около 20 кДж-моль"1, температурные колебания приводят к высоким амплитудам колебаний молекул в кристалле и, особенно, в приповерхностных слоях. Повышение температуры до 260 К приводит к практически полному деструктурированию поверхности и образованию так называемого квазижид-кого слоя, в котором движение молекул уже не может описываться как колебательное движение вблизи кристаллических узлов, а скорее напоминает движение молекул в жидкости. При температурах ниже 220 К (по некоторым данным - ниже 250 К) квазижидкий слой перестает играть существенную роль, и поверхность льда может быть с хорошей степенью точности описана как поверхность молекулярного кристалла.
Для моделирования методом молекулярной динамики наиболее подходящими являются два типа моделей - модель изолированной частицы (молекулярного кластера) и модель с трехмерными периодическими условиями. Первый тип модели в большей степени отражает свойства наночастиц, в то время как частицы большего размера лучше описываются вторым типом моделей [12]. В связи с тем, что в наших исследованиях предполагается участие частиц большого размера (около 1 мкм и выше), для моделирования был выбран второй тип модели.
Моделирование атмосферных ледяных частиц проводилось с использованием системы на основе кристаллического гексагонального льда
Ш. Эта модификация характеризуется способностью к протонному разупорядочиванию, в результате на основе одной и той же кристаллической решетки атомов водорода возникает большое число ориентационных изомеров, различающихся ориентацией водородных связей. Для решения вопроса о выборе подходящего ориентационного изомера используются специальные «упорядоченные» модели, в которых при одной и той же структуре кристаллической подрешетки атомов кислорода выбираются одна или несколько упорядоченных (т.е. обладающих дальним порядком) подрешеток атомов водорода. Одной из таких систем является так называемый лед Пизани [13]. Более современные модели, основанные на оптимизации энергии среди различных ориентационных изомеров, были разработаны в работе Хирша [14]. В настоящее время эти структуры часто используются в теоретических исследованиях, обеспечивая простоту моделирования, структурную устойчивость и низкий поверхностный дипольный момент (свойство, важное для воспроизведения электрических свойств льда). В связи с этим, в данной работе в качестве модели была выбрана ромбическая низкотемпературная равновесная форма с упорядоченным расположением протонов, так называемый лед XI Хирша. По данным
[14], лед XI обладает минимальной энергией среди других ориентационных изомеров (оценка методом DFT), а также имеет чередование терминальных ОН-групп и обеспечивает минимальный дипольный момент кристаллических слоев, что гарантирует отсутствие электрической поляризации всего кристалла.
В данной работе была использована модель кристалла, состоящая из 10 бислоев (1800 молекул Н2О). Эта структура была получена транслированием элементарной ячейки льда XI и представляла собой параллелепипед размерами 40x39x36 А в х, у и г направлениях соответственно. Параметры элементарной ячейки и координаты атомов были взяты из работы [14]. Ось г соответствует направлению нормали к бислою. Начало координат располагалось в центре системы, для создания пространственной структуры льда на модель были наложены трехмерные периодические условия. Плотность кристалла, соответствующего данной модели, составляет
0.96 г/см3.
Молекулы воды рассматривались как жесткие. Для описания взаимодействия молекул Н2О использовалась трехточечная модель воды SPC
[15]. Размерные параметры молекул воды были скорректированы в соответствии с выбранной моделью. Для модели SPC длина связи О-Н составляет 1.0 А, валентный угол - 109.47°. За-
ряды на атомах H и O равны 0.41e и -0.82e соответственно [16]. Взаимодействие молекул H2O в рамках модели SPC описывалось O-O Леннард-Джонсовским взаимодействием.
Для исследования динамики диффузии частиц фотофрагментов на основе описанной модели льда были рассмотрены две системы. Первая модель была получена введением в модель льда 20 атомов H, вторая - 20 радикалов OH. Начальные координаты частиц внутри ячейки выбирались случайным образом. Ориентация OH-радикалов выбиралась случайно в направлении координатных осей. Длина связи O-H гидроксил-радикалов составляла 0.97 Â и была фиксирована в ходе расчета, заряды равны -0.4e и 0.4e для атомов O и H соответственно. Описание взаимодействий Н и ОН с молекулами воды проводилось при помощи потенциальных функций, использованных ранее в работе [8].
Расчеты выполнялись с использованием программного пакета для молекулярно-динамических расчетов DL_POLY [17]. Система представляла собой NVT-ансамбль с наложенными трехмерными периодическими условиями. Интегрирование уравнений движения проводилось при помощи алгоритма Leapfrog Verlet [12]. Поддержание постоянства длин связей осуществлялось при помощи алгоритма SHAKE [18]. Для расчета электростатических взаимодействий использовался метод суммирования по Эвальду [19]. Радиус действия неэлектростатических потенциалов составлял 10 Â. Для последующего анализа результатов динамики проводилась запись текущих координат и смещения от первоначального положения для всех частиц с периодичностью 50 шагов. Молекулярная система предварительно уравновешивалась, как правило, на протяжении 40000 шагов. Для поддержания постоянства температуры использовался термостат Берендсена [20] с периодом релаксации 0.1 пс. При анализе результатов динамики период уравновешивания не учитывался. Коэффициенты самодиффузии рассчитывались с помощью обобщенного уравнения Эйнштейна
2
< r >
D = lim t — ’ (1)
6t
где D - коэффициент самодиффузии, <r2> - средний квадрат смещения (MSD) частицы к моменту времени t. При расчете D значение смещения усреднялось по всем частицам данного типа.
Результаты и обсуждение
Для проверки адекватности используемой модели водного льда был проведен расчет исходной системы (H2O)18oo без радикальных частиц. Система уравновешивалась в течение 10 пс, расчет МД проводился в течение 0.5 нс с шагом 0.25 фс. Расчеты проводились в температурном интервале 50-150 К. Полученные при наивысшей температуре (150 К) зависимости среднего квадрата смещения и коэффициента самодиффузии молекул воды от времени представлены на рис. 1.
Из представленных зависимостей видно, что при 150 К молекулы воды внутри ледяной структуры совершают колебательные движения около положения равновесия, коэффициент самодиффузии H2O стремится к нулю, плавления льда не происходит. Аналогичные результаты были получены при более низких температурах. Этот факт согласуется с результатами предыдущих исследований [21], в которых было показано, что при использовании модели SPC плавление льда происходит при температуре 190 К. Рассчитанные в интервале температур 50-150 К коэффициенты самодиффузии воды представлены в табл. 1.
Таблица 1
Коэффициенты самодиффузии молекул воды при
различных температурах
Температура, К D •109, mV1
50 4.410-4
75 5.310-4
100 6.2-10-4
125 8.210-4
150 1.010-3
В
0>
S
U
0.4-,
0.2
0.1-
0.0
1.4-і
1.2-
1.0-
'о 0.8-
°о 0.6-
Q 0.4-
0.2-
0.0-
100
400
500
0
100
200 300
Время, пс
400
500
200 300
Время, пс
Рис. 1. Зависимость от времени среднего квадрата смещения (а) и коэффициента самодиффузии (б) молекул воды при 150 К
Полученные значения коэффициентов самодиффузии находятся в интервале (10-3^10-4)х10-9 м •с- . В разработанной модели кристаллическая структура льда будет сохраняться, а диффузия водных молекул будет обусловлена их колебаниями около положения равновесия. Определенные экспериментально коэффициенты само-диффузии для молекул воды, приводимые в литературе, как правило, относятся к более высоким температурам. В [22] был экспериментально определен коэффициент самодиффузии воды при 273 К. Полученная величина (1.1 -10-9 м2-с-1) значительно выше, чем полученные нами значения при 50-150 К, что также подтверждает адекватность выбранной нами модели.
Для системы, содержащей 20 атомов H, молекулярно-динамический расчет проводился в интервале 100 пс после уравновешивания в течение 10 пс. Для оценки воспроизводимости расчеты проводились с шагом 0.5, 0.25 и 0.1 фс. Средний квадрат смещения атомов водорода за время расчета составлял 12.0, 7.1 и 7.0 Â соответственно. Поскольку уменьшение шага интегрирования ниже 0.25 фс не приводит к существенному изменения среднего смещения, дальнейшие расчеты проводились при этом значении данного параметра.
Молекулярно-динамическое исследование диффузии атомов H и радикалов OH проводилось после уравновешивания в течение 10 пс, время расчета после уравновешивания составляло 0.5 нс. Расчеты проводили в температурном интервале 50-150 К с шагом 25 К. Зависимости смещения частиц радикалов от времени при 150 K представлены на рис. 2 а и 2б для Н и OH соответственно. В процессе диффузии в системе, содержащей радикалы OH, происходило небольшое разупорядочение кристаллического окружения, при моделировании диффузии атомов водорода разупорядочения практически не происходило.
Анализ полученных данных показывает, что движение частиц в кристалле складывается из низкоамплитудных осцилляций между бислоями водного льда и резких смещений (скачков), которые соответствуют периодическим перескокам через границу бислоя. Перескоки происходят в прямом и обратном направлениях, их разность дает число «эффективных» перескоков, обеспечивающих перемещение частиц внутри кристалла. Таким образом, диффузия Н и ОН представляет собой совокупность частых низкоамплитудных колебаний (медленного дрейфа) между бислоями и значительно более редких перескоков через бислой. Перескоки атомов водорода более выражены, чем радикалов ОН. В случае атомов Н основной вклад в диффузию вносят перескоки, для радикалов ОН влияние этих двух составляющих примерно одинаково. Температурная зависимость количества и частоты перескоков, а также частота эффективных перескоков представлены в табл. 2 и на рис. 3а и 3б.
С увеличением температуры частота перескоков возрастает. Характер температурной зависимости близок к экспоненциальному: для атомов Н в интервале 75-150 К график этой зависимости в полулогарифмических координатах представляет прямую с коэффициентом корреляции 0.942. Для ОН эта зависимость является линейной во всем интервале температур (с отклонением от прямой точки при 75 К), коэффициент корреляции 0.931. При температурах ниже 100 К перескоки водорода на протяжении всего расчета не были зарегистрированы. Для ОН перескоки регистрируются уже при 75 К, что, по-видимому, объясняется большей кинетической энергией ОН как более массивной частицы. Число эффективных перескоков в целом пропорционально общему числу перескоков и увеличивается при возрастании температуры, что особенно заметно для ОН.
Параметры диффузии продуктов фотолиза при различных температурах
Таблица 2
Температура, К Общее число перескоков Число обратных перескоков Число эффективных перескоков Частота перескоков, пс-1 Частота эффективных перескоков, пс-1 D-109, м2-с-1
Атомы H
50 0 0 0 0 0 1.110-3
75 0 0 0 0 0 1.310-3
100 7 2 5 7.0^ 10-4 5.010-4 8.810-3
125 18 3 15 1.810-3 1.510-3 4.8^ 10-2
150 98 30 б8 9.810-3 б.810-3 3.1-10-1
Радикалы OH
50 0 0 0 0 0 5>10-4
75 72 б7 5 7.2-10-3 5.010-4 3.810-2
100 141 133 8 1.410-2 8.010-4 3>10-2
125 38 28 10 3.810-3 1.010-3 б.910-2
150 254 234 20 2.5-10-2 2.010-3 1.010-1
Время, пс Время, ПС
Рис. 2. Характерное перемещение атомов водорода (а) и радикалов ОН (б) внутри ледяного кристалла при 150 К
Всех Ш'рбСКйЕШ
Г, К
3.0-10"
= 2.5-103-§ 2.04 О Ч
¡Л
I1 1.5-10=
| 1.0-10*5
5.0-10 О.П-1
50
■ Вчхк ГШриСКУКИВ — ф— Эффект и внык гкрсскйков
75
100
ГК
] 25
150
Рис. 3. Зависимость от температуры общего числа и числа «эффективных» перескоков атомов водорода (а) и
радикалов ОН (б) внутри кристаллического льда
Т, К
Коэффициенты диффузии, рассчитанные по уравнению (1) в интервале 50-150 К, представлены в табл. 2, их температурная зависимость показана на рис. 4. Отметим, что для расчета коэффициентов диффузии использовались финальные участки МД-траекторий, соответствующие конечным 10% времени расчета, что являлось дополнительной гарантией полного уравновешивания системы.
Зависимость коэффициента диффузии атома водорода от температуры имеет ярко выражен-
ный экспоненциальный характер, то есть подчиняется уравнению Аррениуса. Аналогичная зависимость для гидроксил-радикала близка к линейной. Подобное различие, по-видимому, объясняется большей подвижностью атома водорода внутри ледяного кристалла вследствие его меньшего размера.
Абсолютные значения коэффициентов диффузии атомов водорода DH в интервале температур 75-125 К несколько меньше, чем для ОН, однако, вследствие экспоненциальной зависи-
мости от температуры, быстро возрастает и при 150 К почти в три раза превышает DOH. Подобное различие, возможно, связано с тем, что в указанном интервале при повышении температуры влияние перескоков на величину коэффициента диффузии снижается, а основной вклад в диффузию вносит их перемещение вдоль бислоев. Возможно, медленная диффузия атомов Н при низких температурах объясняется их более сильной склонностью к захвату со стороны атома кислорода окружающих молекул Н2О. Разрушение «мгновенных» комплексов Н...О2Н требует большего запаса кинетической энергии и происходит гораздо легче при повышении температуры. В поддержку этой точки зрения свидетельствует тот факт, что при перемещении гидроксил-радикалов в направлении, перпендикулярном бислоям, их перемещение медленнее, чем атомов водорода, так как частота перескоков ОН ниже. В случае перемещения ОН вдоль бислоя скорость диффузии ОН выше, о чем свидетельствуют коэффициенты диффузии, полученные в интервале температур 75-125 К. Также следует учитывать, что при температурах ниже 100 К движение атомов Н уже нельзя рассматривать как чисто классическое, и, таким образом, заметный вклад в коэффициент диффузии должен давать туннельный эффект [11].
Заключение
В настоящем исследовании проведено молекулярно-динамическое моделирование динамики молекул воды, атомов водорода и радикалов ОН в кристаллической частице льда при температурах от 50 до 150 К, что соответствует условиям предполагаемых гетерогенных фотореакций генерации пероксида водорода, гидроксил-радикалов и других окислителей в верхних слоях атмосферы Земли (в мезосфере) на высотах порядка 80 км. На основании проведенного моделирования показано, что характер движения атомов водорода и гидроксил-радикалов в ледяной частице значительно различается. Движение атомов водорода представляет собой совокупность низкоамплитудного дрейфа (частых движений) и редких высокоамплитудных скачков, причем основной вклад в коэффициент диффузии вносит второй тип движения. В отличие от предшествующих работ, методом МД определены частоты перескоков, установлена их зависимость от температуры для всего интервала температур 50-150 К и оценены коэффициенты диффузии атомов к поверхности ледяной частицы. В случае радикалов ОН два типа движений дают примерно одинаковые вклады в коэффициент диффузии и, таким образом,
движение можно рассматривать как обычную диффузию, характеризуемую общим коэффициентом. Коэффициенты самодиффузии молекул воды, определенные в данной работе, согласуются с известными экспериментальными данными и подтверждают правильность проведенного моделирования.
Работа выполнена при поддержке Российского фонда фундаментальных исследований (проекты 10-0501112, 11-03-00085).
Список литературы
1. Klán P., Holoubek I. // Chemosphere. 2002. V. 46. Р. 1201-1210.
2. Boogert A., Ehrenfreund P. Astrophysics of Dust: ASP Conference Series. San Francisco: ASP, 2004. V. 309. 547 p.
3. Ehrenfreund P., Fraser H. Solid State Astrochemistry. Dordrecht: Kluwer, 2003. 317 p.
4. Ghormley J.A., Hochanadel C.J. // J. Phys. Chem. 1971. V. 75. P. 40-44.
5. Gerakines P.A., Schutte W.A., Ehrenfreund P. // Astron. Astrophys. 1996. V. 312. P. 289-305.
6. Westley M.S., Baragiola R.A., Johnson R.E., Baratta G.A. // Planet. Space Sci. 1995. V. 43. P. 1311-1315.
7. Westley M.S., Baragiola R.A., Johnson R.E., Baratta G.A. // Nature (London). 1995. V. 373. P. 405-407.
8. Andersson S., Al-Halabi A., Kroes G.-J., Dishoeck E.F. // J. Chem. Phys. 2006. V. 124. P. 6471501-6471514.
9. Andersson S., Kroes G.-J., Dishoeck E.F. // Chem. Phys. Lett. 2005. V. 408. P. 415-421.
10. Arasa C., Andersson S., Cuppen H.M., Dishoeck E.F., Kroes G.-J. // J. Chem. Phys. 2010. V. 132. P. 1845101-1845112.
11. Markland T.E., Habershon S., Manolopoulos D.E. // J. Chem. Phys. 2008. V. 128. P. 19450601-19450611.
12. Allen M.P., Tildesley D.J. Computer Simulation of Liquids. Oxford: Clarendon Press, 1989. 385 p.
13. Pisani C., Casassa S., Ugliengo P. // Chem. Phys. Lett. 1996. V. 253. P. 201-208.
14. Hirsch T.K., Ojamae L. // J. Phys. Chem. B. 2004. V. 108. P. 15856-15864.
15. Berendsen H.J.C., Postma J.P.M., Gunsteren W.F., Hermans J. Intermolecular Forces. Dordrecht: Reidel, 1981. 331 p.
16. Mark P., Nilsson L. // J. Phys. Chem. A. 2001.
V. 105. P. 9954-9960.
17. Todorov I.T., Smith W., Trachenko K., Dove M.T. // J. Mater. Chem. 2006. V. 16. P. 1911-1918.
18. Forester T.R., Smith W. // J. Comput. Chem. 1998. V. 19. P. 102-111.
19. Fennell C.J., Gezelter D.J. // J. Chem. Phys. 2006. V. 124. P. 23410401-23410412.
20. Berendsen H.J.C., Postma J.P.M., van Gunsteren W.F., DiNola A., Haak J.R. // J. Chem. Phys. 1984. V. 81. P. 3684-3690.
21. Vega C., Sanz E., Abascal J.L.F. // J. Chem. Phys. 2005. V. 122. P. 1145071-1145079.
22. Chenevert T.L. // Journal of Magnetic Resonance Imaging. 2011. V. 34. P. 983-987.
DYNAMICS AND DIFFUSION OF RADICALS FORMED DURING UV IRRADIATION OF WATER ICE.
MOLECULAR DYNAMICS STUDY
N. V. Svinkov, S.K Ignatov, M.Yu. Kulikov, A.G. Razuvaev
A study of the behaviour of water molecules, hydrogen atoms and hydroxyl radicals in the crystalline water ice in the temperature range from 50 to 150 K has been carried out by the method of classical molecular dynamics. The particle motion character has been determined and diffusion and self-diffusion (for water molecules) coefficients in the given temperature range have been obtained.
Keywords: water ice, simulation, classical molecular dynamics, photolysis processes, radicals, diffusion.
