CC BY
49
УДК 541.8
МОЛЕКУЛЯРНО-ДИНАМИЧЕСКОЕ МОДЕЛИРОВАНИЕ ЖИДКОЙ ВОДЫ С ЭМПИРИЧЕСКИМ ПОТЕНЦИАЛОМ ВОДОРОДНЫХ СВЯЗЕЙ *
В.Е. Петренко2, Н.А. Абакумова1, М.Л. Антипова2, Ю.М. Кесслер2,
А.В. Боровков2, О.В. Ведь2
Кафедра химии, ТГТУ (1);
Институт химии растворов РАН, г. Иваново (2)
Представлена членом редколлегии профессором В.И. Коноваловым
Ключевые слова и фразы: водородные связи; потенциалы парного взаимодействия; функции радиального распределения.
Аннотация: Рассмотрены эмпирические парные потенциалы для воды, учитывающие неэлектростатическую составляющую водородной связи. Проанализирована чувствительность полученных в компьютерном эксперименте функций радиального распределения к этой составляющей. Получены модельные потенциалы, для которых положение и высота первого и второго пиков функции радиального распределения кислород-кислород совпадают с экспериментальными.
Введение
Наиболее широкое применение при расчете свойств жидкой воды нашли сравнительно простые эмпирические потенциалы парного взаимодействия молекул [1 - 4], параметризованные по какому-либо свойству реальной системы. Часто в роли таких свойств выступают экспериментальные функции радиального распределения (ФРР) [2] или внутренняя энергия ивн [4] при Т = 298 К. Тем не менее, ни один модельный потенциал не может полностью воспроизвести экспериментальную ФРР кислород-кислород gOO(r) [5]. Возникает противоречие [6]: либо модель корректно описывает высоту и положение первого пика gOO(r), но при этом отсутствует второй пик в области rOO~4,5 Â, а значение внутренней энергии менее отрицательно, чем ивн(эксп.) = -41,417 кДж/моль [7] (например, TIPS [2], TIP3P [3]), либо хорошо воспроизводит второй пик gOO(r) и ивн, но дает завышенный первый пик gOO(r) (например, SPC [1], TIP4P [3], SPC2 [4]). Как выяснилось, возможность преодолеть это противоречие существует только для моделей с четырехцентровой геометрией молекулы [8]. О том, что четырехточечные модели более перспективны в плане дальнейшего совершенствования, писали и авторы работ [9, 10]. Они пришли к выводу, что трехточечные модели типа SPC и TIP3P из-за своей геометрии не могут быть параметризованы по свойствам льда Ih даже при учете поляризуемости, но такая параметризация осуществима для четырехточечных моделей типа TIP4P [9]. Позднее [10] этот вывод был полностью подтвержден при более детальном исследовании возможности параметризации таких же моделей по свойствам жидкой воды при Т = 298 К.
*
Работа выполнена при поддержке Российского фонда фундаментальных исследований (код проекта № 02-03-32287) и Фонда содействия отечественной науке.
Ю.М. Кесслер2
Модели и особенности моделирования
Ранее [7, 8, 11, 12] мы исследовали потенциальные функции типа
U ( r ) = 4 6OO
CTOO
rOO
CTOO
л6
rOO
-т
т}
и пришли к выводу, что невозможно достичь лучшего согласия с экспериментом, оставаясь в рамках модели, которая при описании водородной связи учитывает только электростатическое взаимодействие. В плане дальнейшего совершенствования простых модельных потенциалов для жидкой воды нам представляется перспективным учет вклада в водородную связь неэлектростатического взаимодействия ОН.
За основу была взята эмпирическая жесткая модель Т1Р4Р [3]. Геометрия модельной молекулы и параметры потенциала (1) остались без изменений, как в классическом варианте, предложенном Йоргенсеном [3]: ZHOH = 104,52°, ЯОН = = 0,9572А, ЛОм = 0,15А, еОО = 0,155 ккал/моль, стОО = 3,1536 А. Но величина эффективного точечного заряда была уменьшена до дН = 0,4415|е|, чтобы дипольный момент модельной молекулы был равен 1,85 Б (величина дипольного момента молекулы воды в газовой фазе [7]). В потенциальную функцию (1) было введено дополнительное слагаемое, учитывающее неэлектростатическое взаимодействие ОН. Мы рассмотрели различные варианты формы записи этого слагаемого: степенную, экспоненциальную, в виде потенциала Бакингема [13] и потенциала Морзе [14], и остановились на степенной функции 12 - 10
U
OH
non-el
( rOH) = 4
U
OH
( OH Y2 ( OH Y0
ra ra
rOH
rOH
Значения параметров (см. табл. 1) менялись в диапазоне
U
OH
(2)
= 1,00... 2,50
ккал/моль и = 1,20________1,90 А. Полученный набор модельных потенциалов
можно разделить на группы, внутри которых изменяется лишь один параметр, что позволяет оценить степень его влияния на результаты моделирования.
Моделирование проводилось в МУТ-ансамбле с использованием программного пакета МОБУ8. Ячейка содержала 216 молекул воды при температуре 298 К и плотности р = 0,997 г/см3, что соответствует экспериментальному значению для жидкой воды при нормальных условиях. Время выхода системы в равновесие составляло 80 пс, время моделирования составляло 100 пс с шагом 2 фс. Интегрирование ньютоновских уравнений движения проводилось с помощью алгоритма Верле [15], с использованием периодических граничных условий. Электростатические взаимодействия учитывались с помощью суммирования по Эвальду. Радиус обрезания (сиМ^радиус) составил 9,3 А.
Обсуждение результатов
Наблюдается следующая закономерность: на модельных ФРР О... О могут существовать три максимума (рис. 1). Один соответствует первому пику на экспериментальной goo(r), два других расположены при rOO ~ 4,5 Â и при rOO > 6 Â (будем называть их соответственно вторым и третьим пиком модельной ФРР
О. О). Для потенциалов, не учитывающих неэлектростатическое взаимодействие ОН, как было упомянуто выше, существуют два варианта. Первый вариант: на gOO(r) присутствуют все три пика, но высота первого максимума завышена по сравнению с экспериментом (модели SPC, TIP4P и т. п., см. кривую 1 на рис. 1).
Рис. 1 ФРР кислород-кислород для моделей:
1 - Т1Р4Р стандартная [3]; 2 - Т1Р4Р с измененным дН = 0,4415|е|; 3 - Т1Р4Р с дН = 0,4415|е|, дополненная учетом неэлектростатического взаимодействия ОН (2),
1^1 = 9,2048 кДж/моль, г°Н = 1,4 А
Второй вариант: модель корректно воспроизводит высоту первого максимума, но второй пик отсутствует (модели TIPS, TIP3P, а также TIP4P, у которой точечный заряд на атомах водорода равен 0,4415|е|, но потенциальная функция (1) не дополнена неэлектростатическим взаимодействием ОН, см. кривую 2 на рис. 1). При вводе в потенциальную функцию (1) дополнительного слагаемого (2) появляется третий, промежуточный, вариант. Условие (3) может соблюдаться одновременно и для первого, и для второго пиков goo(r).
Нашей целью было получить модельный потенциал, ФРР которого удовлетворяли бы следующим условиям:
Imax < 3 0 r 2max
gOO <3,0, rOO
4,5a , gOmæ > gOmax > «Omax.
Imax
p = < і
1 g 2max
gOH
Imax
P2 = gHH- > 1
^2max
gHH
справедливость которых подтверждена экспериментально [5, 16].
Были рассмотрены различные комбинации параметров
U
OH
OH
(3)
(см.
модели А0 - Б4, табл. 1). Полученный набор модельных потенциалов можно разделить на группы, внутри которых изменяется лишь один параметр, что позволяет оценить степень его влияния на результаты моделирования. Результаты расчетов приведены в табл. 2. Как можно видеть, для некоторых из рассматриваемых потенциалов (модели В2-2, В2-11, Б4) условие (3) выполняется практически полностью, за исключением соотношения первого и второго пиков ФРР О_О.
Убедившись, что дальнейшее изменение параметров
U
OH
OH
и ra при со-
хранении той же геометрии модели не приводит к желаемому результату, мы начали варьировать расстояние между центром атома кислорода и смещенным по
1030 ISSN 0136-5835. Вестник ТГТУ. 2004. Том 10. № 4А. Transactions TSTU.
Параметры моделей
Модель АО Al A2 A3 A4 A4-1 A5 A5-1 A6 Al A8 A9 A10 Ail A12
|J.,D 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85
Чя - Є 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415
Гоя, А 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572
г ом. А 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15
/нон, ° 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52
sqq , ккал/моль 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155
Ооо. А 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536
иоя mm , ккал/молі 2,20 2,20 2,20 2,20 2,20 2,20 2,20 2,20 2,20 2,20 2,20 2,20 2,20 2,20 2,20
r°H. 1,600 1,400 1,350 1,300 1,325 1,500 1,275 1,800 1,200 1,450 1,360 1,380 1,430 1,410 1,370
Модель А13 A14 Bl-1 Bl-2 B2-1 B2-2 B2-3 B2-4 B2-5 B2-6 B2-7 B2-8 B2-9 B2-10 B2-11
M-.D 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85
Чя - е 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415
гон. А 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572
f ом. А 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15
/НОН, ° 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52
Sqq , ККал/мОЛЬ 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155
Ооо. А 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536
иоя пип , ккал/молі 2,20 2,20 2,00 2,50 2,00 2,50 2,30 2,10 2,15 2,25 2,40 2,18 2,12 2,16 2,45
r°H-Â 1,700 1,750 1,300 1,300 1,400 1,400 1,400 1,400 1,400 1,400 1,400 1,400 1,400 1,400 1,400
ISSN 0136-5835. Вестник ТГТУ. 2004. Том 10. № 4А. Transactions TSTU. 1031
Продолжение табл. 1
Модель ВЗ-1 ВЗ-2 В4 В4-1 В4-2 В5 С1 С2 сз D1 D2 D3 D4 El Е2
M-,D 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85
е 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4415 0,4519 0,4519
Пэн. Â 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572
г ОМ, Â 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,15 0,16 0,16
/нон, ° 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52
е00, ккал/моль 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155
Ооо, к 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536
t/OH mm , ккал/мол! 2,40 2,45 2,093 2,3572 2,3024 2,591 1,50 1,20 1,00 2,50 2,50 2,50 2,50 2,20 2,50
Г™, к 1,450 1,450 1,4142 1,400 1,4142 1,415 1,800 1,800 1,900 1,360 1,420 1,500 1,410 1,400 1,400
Модель ЕЗ F1 F2 F3 G1 G2 Н1 Н2
M-.D 1,85 1,85 1,85 1,85 1,85 1,85 1,85 1,85
9н , е 0,4519 0,4627 0,4627 0,4627 0,3961 0,3961 0,4987 0,4987
Гон. à 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572 0,9572
г ОМ. к 0,16 0,17 0,17 0,17 0,10 0,10 0,20 0,20
/нон, ° 104,52 104,52 104,52 104,52 104,52 104,52 104,52 104,52
е00, ккал/моль 0,155 0,155 0,155 0,155 0,155 0,155 0,155 0,155
Ооо. к 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536 3,1536
иоя mm , ккал/мол! 2,50 2,20 2,50 2,40 2,20 2,50 2,20 2,50
r°H. 1,350 1,400 1,400 1,400 1,400 1,400 1,400 1,400
1032 ISSN 0136-5835. Вестник ТГТУ. 2004. Том 10. № 4А. Transactions TSTU.
Характеристики функций радиального распределения и значения внутренней энергии, полученные в МД эксперименте
Модель А0 Al А2 A3 A4 A4-1 A5 A5-1 A6 A7 A8 A9 A10 All A12
(7ВН, кДж/моль -41,1 -30,9 -28,8 -27,4 -28,1 -36,1 -26,8 -48,6 -25,5 -33,5 -29,2 -30,0 -32,2 -31,5 -29,5
lmax д 'ОО , л 2,744 2,744 2,767 2,813 2,790 2,720 2,813 2,860 2,860 2,720 2,767 2,767 2,744 2,744 2,767
lmax goo 4,318 2,632 2,416 2,313 2,350 3,458 2,277 5,400 2,249 3,002 2,463 2,539 2,828 2,727 2,503
2 max і гОО , А 4,580 4,325 4,511 4,394 4,394 4,441 4,883 5,557 4,604 4,418 4,255 4,301 4,301 4,185 4,487
2 max §00 1,082 1,017 0,974 0,940 0,958 1,080 0,930 1,299 0,895 1,058 0,990 1,006 1,031 1,023 0,986
Зтах д 'ОО -л 6,626 6,510 6,464 6,231 6,301 6,650 6,440 7,556 6,231 6,487 6,347 6,510 6,510 6,650 6,510
Зтах Soo 1,053 1,037 1,039 1,040 1,047 1,053 1,048 1,048 1,060 1,045 1,043 1,047 1,044 1,045 1,046
Рх 1,413 0,880 0,781 0,716 0,746 1,163 0,693 1,398 0,643 1,016 0,800 0,839 0,960 0,910 0,814
Р2 1,083 0,991 0,980 0,985 0,977 1,037 0,980 1,199 0,979 1,019 0,975 0,982 0,998 0,997 0,981
Модель А13 А14 Bl-1 Bl-2 B2-1 B2-2 B2-3 B2-4 B2-5 B2-6 B2-7 B2-8 B2-9 B2-10 B2-11
UBH, кДж/моль -43,9 -45,3 -27,0 -29,7 -29,9 -32,2 -31,3 -30,6 -30,6 -31,1 -31,7 -30,8 -30,6 -30,8 -32,1
Ітах Д гОО ' л 2,813 2,837 2,813 2,790 2,767 2,720 2,744 2,744 2,744 2,744 2,720 2,744 2,744 2,744 2,720
Ітах 8 оо 4,694 4,954 2,289 2,335 2,535 2,783 2,697 2,606 2,609 2,651 2,723 2,610 2,594 2,634 2,771
2 max д 'ОО , А 5,487 5,487 4,720 4,743 4,348 4,301 4,418 4,325 4,255 4,371 4,325 4,441 4,325 4,325 4,278
2 max Soo 1,051 1,158 0,925 0,945 0,999 1,039 1,028 1,007 1,008 1,018 1,031 1,004 1,009 1,025 1,037
Зтах X >00 - А 6,789 8,672 6,231 6,208 6,394 6,510 6,650 6,533 6,510 6,673 6,650 6,580 6,533 6,417 6,580
Зтах Soo 1,019 1,010 1,052 1,039 1,043 1,055 1,047 1,046 1,045 1,046 1,045 1,047 1,042 1,045 1,053
Рх 1,469 1,414 0,700 0,739 0,842 0,944 0,908 0,865 0,867 0,887 0,920 0,874 0,867 0,880 0,994
Р2 1,094 1,267 0,980 0,976 0,980 1,002 0,994 0,987 0,984 0,992 0,997 0,987 0,981 1,027 1,000
ISSN 0136-5835. Вестник ТГТУ. 2004. Том 10. № 4А. Transactions TSTU. 1033
Модель ВЗ-1 ВЗ-2 В4 В4-1 В4-2 В5 С1 С2 СЗ D1 D2 D3 D4 El Е2
(7ВН, кДж/моль -34,8 -35,0 -31,1 -31,6 -32,0 -33,7 -38,1 -34,9 -33,5 -30,5 -33,4 -38,6 -32,8 -31,6 -33,1
lmax Д '00 - л 2,720 2,720 2,744 2,744 2,744 2,720 2,883 2,883 2,953 2,744 2,720 2,697 2,720 2,744 2,720
lmax ёоо 3,177 3,206 2,693 2,715 2,770 2,988 4,115 3,716 3,705 2,587 2,995 3,802 2,865 2,689 2,881
2 max Д гоо ,л 4,394 4,371 4,232 4,255 4,301 4,464 5,580 5,696 5,743 4,278 4,278 4,394 4,371 4,348 4,348
2 max ёоо 1,082 1,092 1,021 1,042 1,035 1,071 1,145 1,091 1,081 1,024 1,061 1,122 1,053 1,017 1,055
Зтах Д гОО ' л 6,673 6,650 6,580 6,487 6,533 6,650 8,579 8,905 8,603 6,394 6,487 6,719 6,533 6,533 6,650
Зтах ёоо 1,056 1,056 1,042 1,048 1,050 1,052 1,032 1,028 1,042 1,048 1,050 1,064 1,050 1,045 1,050
Pi 1,084 1,101 0,904 0,925 0,941 1,020 1,187 1,078 1,015 0,864 1,010 1,291 0,976 0,904 0,978
Р2 1,030 1,039 0,991 0,994 1,000 1,020 1,124 1,086 1,146 0,987 1,011 1,079 1,006 0,987 1,007
Модель ЕЗ F1 F2 F3 G1 G2 Н1 Н2
(7ВН, кДж/моль -30,7 -32,7 -34,1 -33,7 -27,6 -29,0 -36,4 -38,2
lmax Д гоо ,л 2,744 2,720 2,720 2,720 2,767 2,744 2,720 2,697
lmax ё00 2,537 2,777 2,977 2,907 2,358 2,525 3,129 3,323
2 max д гОО , А 4,348 4,418 4,371 4,418 4,650 4,604 4,325 4,301
2 max ёоо 1,015 1,044 1,074 1,073 0,973 0,993 1,109 1,142
Зтах Д гОО ' л 6,487 6,510 6,487 6,487 6,278 6,278 6,626 6,394
Зтах 600 1,045 1,049 1,051 1,051 1,050 1,031 1,065 1,076
Pi 0,848 0,941 1,012 0,991 0,770 0,848 1,070 1,138
Р2 0,981 0,998 1,016 1,010 0,982 0,990 1,029 1,074
биссектрисе угла HOH отрицательным точечным зарядом Лом. Одновременно менялась и величина эффективного точечного заряда qH таким образом, чтобы соблюдалось условие ц = 1,85 D. Таким образом, были получены модели E1 - H2 (табл. 1). Как показали результаты расчетов (см. табл. 2), условию (3) полностью соответствуют два парных потенциала: E2 и F3.
Полученные данные позволили обнаружить ряд закономерностей и тем самым существенно упростили процесс параметризации.
Были получены 43 модели, отличающиеся между собой только значениями
параметров
U
OH
он
(модели A1 - D4). Среди них можно выделить три
группы. Для моделей группы I значение
U
OH
изменялось в интервале от 2,0
до 2,5 ккал/моль при постоянном rc = 1,4 Á (модели B2-1, B2-4, B2-9, B2-5,
B2-8, B2-6, B2-3, B2-7, B2-11, B2-2, табл. 1, 2). В группах II и III величина
OH
U
он
оставалась постоянной,
изменялась. Для группы II -
U
min
OH
= 2,2 ккал/моль r°H = 1,2...1,8 Â (модели A6, A5, A3, A4, A2, A8, A12, A9, A1,
A11, A10, A7, A4-1, A0, A13, A14, A5-2, табл. 1, 2); III -
U
OH
= 2,5 ккал/моль,
OH
г^" = 1,3... 1,5 А (модели В1-2, Б1, Б4, Б2, Б3, табл. 1, 2). Таким образом, в
пределах каждой группы изменялся лишь один параметр, что позволило оценить его влияние на результаты численного эксперимента.
На рис. 2 - 5 приведена графическая иллюстрация влияния изменения параметров функции (2) на результаты компьютерного моделирования (,
4Г, Р и Р2).
Как можно видеть, модели группы I ведут себя вполне прогнозируемо:
гон
gOT и Р\ монотонно увеличиваются с ростом
на всем рассматриваемом
интервале. На зависимости gO™aX от
U
он
наблюдается небольшое отклонение,
но мы пока затрудняемся объяснить, с чем она связана. Что же касается зависимого, вообще говоря, трудно ожидать четкой корреляции количе-
сти Р2 от
U
он
ственной характеристики ФРР Н.Н с параметром потенциала ОН.
На примере моделей группы II можно проследить, как меняется поведение
ОН
первого и второго пиков на gOO (г) с изменением параметра гс . Как известно, все предложенные к настоящему времени в литературе модельные потенциалы для воды дают два варианта поведения ФРР О...О: либо эксперименту соответствует gO1Oax , но при этом полностью отсутствует второй пик в области 4,5 А, либо
в согласии с экспериментом находятся гОтах и g^ax , но при этом gOoax > 3,0. Введение в потенциальную функцию (1) слагаемого (2) привело к тому, что при любом наборе параметров и, соответственно, при любой gOoaX на gOO (г) присутствует второй пик при г ~ 4,5 А, хотя во многих случаях его высота и не превышает 1,0.
Рассмотрим рис. 2, б. При г(°)Н < 1,38 А (точки 11 - 17) gOoaX < 1,0 . Такие
потенциальные функции недооценивают взаимодействия частиц во второй сольватной оболочке и описывают систему, аналогичную дипольной жидкости без
^ 1 max
S oo 2,85 -, 2,80 2,75 2,70 -2,65 -2,60 -2,55 -2,50
©10
©8
©9
©7
©6
’-Є&Р 05
1,9
©1 —I—
2,0
I
2,1
g
1 т ах
I
2,2
а)
I
2,3
U
2,4
он
і 1 і
2,5 2,6
ккал/моль
g
1 max
б)
Г 0 Н !
а ’ А
в)
„ОН о
а Л
Рис. 2 Зависимости gOoaX :
ІГ7-ОНІ ~ т л ОН
а - от величины р тт в ряду моделей группы I; б - от величины га в ряду моделей
группы II; в - от величины г(°Н в ряду моделей группы III. Пунктиром обозначена граница, определенная условиями (3). Номера точек соответствуют моделям: 1 - В2-1, 2 - В2-4, 3 - В2-9, 4 - В2-5, 5 - В2-8, 6 - В2-6, 7 - В2-3, 8 - В2-7, 9 - В2-11, 10 - В2-2, 11 - А6, 12 - А5, 13 - А3, 14 - А4, 15 - А2, 16 - А8, 17 - А12, 18 - А9, 19 - А1, 20 - А11, 21 - А10, 22 - А7, 23 - А4-1, 24 - А0, 25 - А13, 26 - А14, 27 - А5-2, 28 - В1-2, 29 - Б1, 30 - Б4, 31 - Б2, 32 - Б3
~2 max
£00 1,04—
1,031,021,011,00-0,99
1,9
g
08
07
Об
>gP 04
о,-
05
2,0
I
2,1
2,2 a)
—і—'—і—'—і—■—і
2,3 2,4 2,5 2,6
U°m"n > ккал/моль
^ r> 2max
Рис. 3 Зависимости g0Q :
а - от величины
и
OH
Т , ОН
в ряду моделей группы I; б - от величины га в ряду моделей
ОН
группы II; в - от величины в ряду моделей группы III. Пунктиром обозначена граница, определенная условиями (3). Номера точек соответствуют моделям: 1 - В2-1, 2 - В2-4, 3 - В2-9, 4 - В2-5, 5 - В2-8, 6 - В2-6, 7 - В2-3, 8 - В2-7, 9 - В2-11, 10 - В2-2, 11 - А6, 12 - А5, 13 - А3, 14 - А4, 15 - А2, 16 - А8, 17 - А12, 18 - А9, 19 - А1, 20 - А11, 21 - А10, 22 - А7, 23 - А4-1, 24 - А0, 25 - А13, 26 - А14, 27 - А5-2, 28 - В1-2, 29 - Б1, 30 - Б4, 31 - Б2, 32 - Б3
Р\ 0,975 -і 0,950 -0,925 -0,900 -0,875 -0,850 -0,825
Ol -<--------1—
1,9 2,0
iC&oP5
07
06
—I---------1--------1—
2,1 2,2
і
2,3
09
Oio
08
—I----'---1---»---1
2,4 2,5 2,6
t/“" , ккал/моль
r » » а
a ’ ^
/QH
■^min
Рис. 4 Зависимости Р
rQH
'a
группы II; в - от величины га в ряду моделей группы III. Пунктиром обозначена грани-
ца, определенная условиями (3). Номера точек соответствуют моделям: 1 - В2-1, 2 - В2-4, 3 - В2-9, 4 - В2-5, 5 - В2-8, 6 - В2-6, 7 - В2-3, 8 - В2-7, 9 - В2-11, 10 - В2-2, 11 - А6, 12 - А5, 13 - А3, 14 - А4, 15 - А2, 16 - А8, 17 - А12, 18 - А9, 19 - А1, 20 - А11, 21 - А10, 22 - А7, 23 - А4-1, 24 - А0, 25 - А13, 26 - А14, 27 - А5-2, 28 - В1-2, 29 - Б1, 30 - Б4, 31 - Б2, 32 - Б3
а
г°н А
a ? ^
.0 H 9
a Д
Рис. 5 Зависимости Р2:
ОН ОН
от величины \итт в ряду моделей группы I; б - от величины га в ряду моделей
группы II; в - от величины га в ряду моделей группы III. Пунктиром обозначена грани-
ца, определенная условиями (3). Номера точек соответствуют моделям: 1 - В2-1, 2 - В2-4, 3 - В2-9, 4 - В2-5, 5 - В2-8, 6 - В2-6, 7 - В2-3, 8 - В2-7, 9 - В2-11, 10 - В2-2, 11 - А6, 12 - А5, 13 - А3, 14 - А4, 15 - А2, 16 - А8, 17 - А12, 18 - А9, 19 - А1, 20 - А11, 21 - А10, 22 - А7, 23 - А4-1, 24 - А0, 25 - А13, 26 - А14, 27 - А5-2, 28 - В1-2, 29 - Б1, 30 - Б4, 31 - Б2, 32 - Б3
а
водородных связей [7]. У моделей, соответствующих точкам 18 - 24, (1,38 А < г°Н < 1,60 А) gOSax > 1,0 , г(2,тах ~ 4,5 А, то есть поведение модельной жидкости во второй сольватной оболочке близко к поведению реальной системы, условие (3) выполняется. Эти потенциалы представляют для нас наибольший интерес. И, наконец, для моделей, у которых Г°Н > 1,60 А (точки 25 - 27) можно видеть, что второй пик на gOO (г) сливается с третьим, обычно расположенным при гОО > 6,0 А (значения г^хх, gOoaX приведены в табл. 2),
gOOaX > 4,0, гОО1- > 5,0 А. То есть, мы имеем дело с жидкостью, в которой настолько сильны межмолекулярные взаимодействия, что первая сольватная оболочка гораздо более насыщена, чем в реальной воде, а вторая и третья слились между собой.
На графиках зависимостей gOoax , gooaX , Р и Р2 от параметра г°Н в группе II (рис. 2, б - 5, б) можно видеть следующую особенность: на интервале г<°Н = 1,3...1,5 А все они носят линейный или почти линейный характер. Такая
7 ОН
же закономерность наблюдается и при
U
min
= 2,5 ккал/моль для моделей груп-
пы III (см. рис. 2, в - 5, в).
Итак, полученные модельные потенциалы достаточно корректно воспроизводят структуру жидкой воды. Чтобы оценить их возможности при описании динамических и энергетических свойств, мы провели для них расчет внутренней энергии ивн и коэффициента самодиффузии Б.
Как известно, экспериментальное значение внутренней энергии составляет -9,92 ккал/моль [7]. Из табл. 2 можно видеть, что для наших моделей расхождение с экспериментом составляет как минимум 2 ккал/моль. Это еще раз подтверждает невозможность одновременной параметризации потенциальных функций для воды по ФРР и внутренней энергии (ранее [6] мы достаточно подробно исследовали этот факт на примере модельных потенциалов типа (1)).
Список литературы
1. Berendsen H. J. C. Interaction models for water in relation to protein hydration / H.J.C. Berendsen, J.P.M. Postma, W.F. van Gunsteren, J. Hermans // in Intermolecular Forces (ed. Pullman B.). - Reidel: Dordrecht, 1981. - Pp. 331-338.
2. Jorgensen W.L. Transferable intermolecular potential functions for water, alcohols and ethers. Application to liquid wate / W. L.Jorgensen // J. Amer. Chem. Soc. -1981. - Vol. 103, No. 2. - Pp. 335-340.
3. Jorgensen W. L. Comparison of simple potential functions for simulating of liquid water / W. L. Jorgensen, J. Chandrasekhar, J.D. Madura // J. Chem. Phys. - 1983. -Vol. 79, No. 2. - Pp. 926-935.
4. Berendsen H.J.C., Grigera J. R., Straatsma T.P. The missing term in effective pair potential H.J.C. Berendsen, J.R. Grigera, T.P. Straatsma // J. Phys. Chem. - 1987. -Vol. 91, No. 24 - Pp. 6269-6271.
5. Soper A. K. Orientational correlation functions molecular liquids: The case of liquid water / A.K. Soper // J. Chem. Phys. - 1994. - Vol. 101, No. 8. - Pp. 6888-6901.
6. Петренко В.Е. Вода в компьютерном эксперименте: противоречие в параметризации потенциалов / В.Е. Петренко, М.Л. Дубова, Ю.М. Кесслер, М.Ю. Перова // Журн. физ. химии. - 2000. - Т. 74, № 11. - С. 1957-1961.
7. Эйзенберг Д. Структура и свойства воды / Д. Эйзенберг, В. Кауцман - Л.: Гидрометеоиздат, 1975. - 280с.
8. Антипова М. Л. Параметры моделей и условные коэффициенты чувствительности функций радиального распределения воды по данным молекулярной динамики / М.Л. Антипова, В.Е. Петренко, Ю. М. Кесслер // Известия РАН, Серия химическая. - 2003. - №2. - С. 316-321.
9. Brodholt J. Parametrizing polarizable intermolecular potentials for water with the ice Ih phase / J. Brodholt, M. Sampoly, R. Vallauri // Mol. Phys. - 1995 - Vol. 85, No. 1. - Pp. 81-90.
10. Brodholt J. Parametrizing a polarizable intermolecular potentials for water / J. Brodholt, M. Sampoly, R. Vallauri // Mol. Phys. - 1995 - Vol. 86, No. 1. - Pp. 149-158.
11. Кесслер Ю. М. Трехточечные и четырехточечные эмпирические жесткие модели воды: новая информация о свойствах / Ю.М. Кесслер, В.Е. Петренко, М. Л. Антипова // Журн. физ. химии. - 2002. - Т. 76, № 2. - С. 280-286.
12. Кесслер Ю. М. Соотношение производных в двух зарактерных точках функции потенциальной энергии парного взаимодействия молекул U(r) / Ю.М. Кесслер, В.Е. Петренко, М.Л. Антипова // Доклады Академии наук. - 2002. Т. 382, № 5. - С. 653-656.
13. Смирнова Н.А. Молекулярные теории растворов / Н.А. Смирнова. - Л.: Химия, 1987. - 336 с.
14. Robinson W.G. Water in Biology, Chemistry and Physics / W. G. Robinson,
S.-B. Zhu, S. Singh, M.W. Evans. - World Scientific, Singapore, 1996. - 509 p.
15. Allen M. P. Computer simulation of liquids / M.P. Allen, D.J. Tildesley. -Clarendon Press, London, 1987. - 387 p.
16. Soper A. K. A new determination of the structure of water at 25 °C / A.K. Soper, M.G. Phillips // Chem. Phys. - 1986. - Vol. 107, No. 1. - Pp. 47-60.
17. Ben-Naim A. Aspects of the statistical-mechanical theory of water / A. Ben-Naim, F.H. Stillinger // in Water and Aqueous Solutions: Structure, Thermodynamics, and Transport Processes, Wiley-Intersc., N.-Y., 1972.
Molecular-Dynamic Modeling of Liquid Water with Empiric Potential of Hydrogen Bonds
V.E. Petrenko2, N.A. Abakumova1, M.L. Antipova2, Yu.M. Kessler2,
A.V. Borovkov2, O.V. Ved2
Department of Chemistry, TSTU (1);
Institute of Solutions Chemistry RAN, Ivanovo (2)
Key words and phrases: potentials of double influence; hydrogen bonds; radial distribution functions.
Abstract: Empiric double water potentials accounting for non-electrostatic constituent of hydrogen bond are studied. The sensitivity of obtained through computer experiment functions of radial distribution to this constituent is analyzed. Model poten-
tials in which location and height of the first and second peaks of radial distribution oxygen-hydrogen function coincide with experimental ones are obtained.
Molekulardynamische Modellierung des flüssigen Wassers mit dem empirischen Potential der Wasserstoffverbindungen
Zusammenfassung: Es sind die die unelektrostatische Komponente der Wasserstoffverbindungen berücksichtigenden empirischen Doppelpotentiale für das Wasser betrachtet. Es ist die Empfindlichkeit der im Computerexperiment erhaltenen Funktionen der radialen Einteilung zu dieser Komponente analysiert. Es sind die Modellpotentiale, für die die Lage und die Höhe der ersten und zweiten Spitzenwerte der Funktion der radialen Einteilung “Sauerstoff - Sauerstoff’ mit den experimentellen Potentialen übereinstimmen, erhalten.
Modélage moléculaire et dynamique de l’eau liquide avec un potentiel empirique des liaisons d’hydrogène
Résumé: Sont envisagés les potentiels formant la paire pour l’eau tenant compte de la composante non électrostatique de la liaison d’hydrogène. Est analysée la sensibilité des fonctions de la répartition radiale envers cette composante qui sont reçues au cours d’un expériment d’ordinateur. Sont obtenus les potentiels modélés pour lesquels la situation et la hauteur de la première et de la deuxième pointes de la fonction de la répartition radiale hydrogène - hydrogène coïncident avec celles de l’expériment.
