Письма в ЖЭТФ, том 116, вып. 5, с. 313 - 317

Модифицированная молекулярная модель для описания

конденсированных фаз тяжелой воды

В.Р.Белослудов, К.В.Гец¹⁾, Р.К.Жданов, Ю.Ю.Божко, О.С.Субботин

Институт неорганической химии им. А. В. Николаева Сибирского отделения РАН, 630090 Новосибирск, Россия

Поступила в редакцию 15 июля 2022 г.

После переработки 27 июля 2022 г.

Принята к публикации 27 июля 2022 г.

В рамках метода решеточной динамики предложена простая трехточечная модель тяжелой воды,

полученная путем оптимизации теплового расширения и температуры плавления дейтерированого льда

гексагонаьной структуры I_h к экспериментальным значениям через параметры леннард-джонсовского и

кулоновского потенциалов между узлами с фиксированными зарядами, описывающих межмолекуляр-

ные взаимодействия молекул воды. Метод решеточной динамики учитывает энергию нулевых колебаний,

что позволяет более точно производить построение термодинамических функций. В работе приведены

результаты расчета плотности колебательных состояний, температурной зависимости линейных разме-

ров элементарной ячейки и температуры плавления льда I_h при различной концентрации молекул D2O.

DOI: 10.31857/S1234567822170086, EDN: jtaopa

Введение. Тяжелая вода - это вода, в которой

тролиза, вплоть до настоящего времени тяжелая во-

атомы водорода (протий,¹H) замещены изотопами

да находится в центре внимания большого количе-

водорода²H или, другими словами, дейтерием (D)

ства физико-химических исследований [9], так как

с ядром, состоящим из одного протона и одного ней-

тяжелая вода играет значительную роль во многих

трона. В результате молекулярная масса тяжелой во-

природных и технологических процессах. Одним из

ды выше, чем у обычной воды примерно на 11 %.

важных свойств тяжелой воды является то, что она

Несмотря на почти идентичные структуры, образу-

практически не поглощает нейтроны, поэтому может

емые молекулами H₂O и D₂O, тяжелая и обычная

использоваться в качестве замедлителя нейтронов и

вода имеют различие в некоторых физических свой-

одновременно теплоносителем в ядерных реакторах.

ствах. При этом различия свойств тяжелой воды по

Для проведения численных экспериментов требу-

сравнению с обычной водой более значительные, чем

ется иметь точную модель молекулы воды, которая

просто изменение массы молекулы, в рамках класси-

воспроизводит требуемые физико-химические свой-

ческой концепции механической теории жидкой воды

ства, например, плотность, электрическую воспри-

[1]. Дейтерий, из-за своей массы имеет более низкую

имчивость, энергетические параметры, коэффици-

энергию нулевых колебаний, что делает “дейтерие-

ент диффузии и пр. За последние несколько десяти-

вую связь” более сильной, чем “водородная связь”

летий было разработано множество различных моде-

[2]. Эта разница может быть причиной более проч-

лей H₂O [10]. Например, трехточечная модель воды

ной структуры тяжелой воды [3, 4]. При относитель-

SPC (Single Point Charge), предложенная Берендсе-

но высоких температурах вязкость и теплоемкость

ном и соавторами в 1981 г., была одной из первых мо-

D₂O более значительны, чем аналогичные парамет-

делей, но до сих пор широко используется [11]. Позд-

ры в H₂O [5]. Кроме этого, D₂O имеет повышенную

нее данная модель была доработана, что вылилось в

температуру плавления по сравнению с H₂O, а мак-

создание модели SPC/E [12].

симальная плотность D₂O достигается при темпера-

Создание моделей тяжелой воды, необходимое

туре 11.2^◦C вместо 4^◦C [5]. Изотопический эффект

для воспроизведения экспериментальных свойств

проявляется не только в изменении макроскопиче-

D₂O, проходило не так интенсивно и стратегия

ских свойств, но может быть обнаружен и на кван-

исследований заключалась не в создании новых

товом уровне [6,7].

моделей с нуля, а в перепараметризации уже

Начиная с открытия дейтерия [8] в 1932 г. и его

существующих моделей H₂O, сохраняя их специфи-

выделения в чистом виде в 1933 г. с помощью элек-

ческие характеристики. Так, модель для тяжелой

воды на основе модели SPC/E, обозначенная как

¹⁾e-mail: gets@niic.nsc.ru

SPC/HW [3], была создана путем изменения заряда

Письма в ЖЭТФ том 116 вып. 5 - 6

2022

313

314

В.Р.Белослудов, К.В.Гец, Р.К.Жданов и др.

Таблица 1. Сравнительные значения зарядов на атомах кислорода q_O и водорода q_H, пространственной σ и энергетической

ε констант потенциала Леннард-Джонса, а также внутримолекулярных угла водород-кислород-водород a_HOH и расстояния

кислород-водород d_OH для потенциалов SPC/E_mod-HW, SPC/E_mod, SPC/HW и SPC/E

SPC/E_mod-HW

SPC/E_mod

SPC/HW

SPC/E

q_H, |e|

+0.435

+0.4238

+0.435

+0.4238

q_O, |e|

-0.87

-0.8476

-0.87

-0.8476

σ,Å

3.1781

3.1556

3.166

ε, кДж/моль

0.65063

0.650

aHOH,^◦

109.47

d_OH,Å

на атоме кислорода q_O с -0.8476 |e|, принятом в

SPC/E [12] и другими часто используемыми потенци-

модели SPC/E, до -0.87 |e|, оставляя все остальные

алами. Величины заряда предлагаемого потенциала

параметры за исключением массы водорода без из-

SPC/E_mod-HW были взяты из потенциала SPC/HW,

менений. Моделирование при помощи молекулярной

разработанного для описания тяжелой воды [3], в то

динамики, выполненное с потенциалом SPC/HW

время как остальные силовые константы остались

показывает хорошее согласие с экспериментальными

равными константам из потенциала SPC/E_mod. Для

значениями по основным свойствам тяжелой воды.

описания молекул H₂O в данной работе использует-

Кроме этого, имеются адаптации под тяжелую воду

ся потенциал SPC/E_mod. Параметры разработанного

других известных потенциалов обычной воды [3].

и исходных потенциалов приведены в табл. 1. Моле-

Несмотря на то, что методы рентгеноструктурно-

кулы воды считались жесткими. Параметры потен-

го анализа [13] значительно продвинулись, класси-

циала Леннард-Джонса между разными типами ато-

ческие рентгеноструктурные методы до сих пор ак-

мов определялись при помощи комбинаторных пра-

тивно используются. При изучении водных кристал-

вил Лоренца-Берселота [19, 20].

лических структур, образование исследуемого веще-

Модель льда I_h состояла из 128 молекул H₂O.

ства может производиться с использованием D₂O

Протоны были расставлены случайным образом, но в

вместо H₂O, так как такая замена позволяет опре-

соответствии с правилом Бернала-Фаулера [21] с ми-

делить положение не только атомов кислорода, но

нимальным значением суммарного дипольного мо-

и дейтерия, что связано с большей площадью сече-

мента решетки. Составы с различным содержанием

ния рассеяния дейтерия [14], а также она позволя-

молекул тяжелой воды получались случайной заме-

ет исследовать изотопический эффект. Так, напри-

ной определенного числа молекул H₂O на молекулы

мер, представлены результаты исследования при по-

D₂O.

мощи рентгеновских и нейтронных методов структур

Восстановление термодинамических потенциалов

льдов [15] и гидратов [16], полученных из тяжелой

льдов производилось с использованием методов ре-

воды.

шеточной динамики в квазигармоническом прибли-

В данной работе ставится задача построения

жении [22], которые позволяют учитывать квантовые

потенциала дейтерированной воды для описания

эффекты, например, нулевые колебания. Для этого

структурных и термодинамических свойств льдов, а

был получен оптимизированный по потенциальной

также исследовать изотопический эффект, основы-

энергии набор структур с различным набором объе-

ваясь на экспериментальных данных (температура

мов элементарной ячейки для получения уравнения

плавления и тепловое расширение).

состояния. Использование метода решеточной дина-

Детали вычислений. Для описания межмоле-

мики позволяет непосредственно вычислять вектора

кулярного взаимодействия молекул тяжелой воды в

собственных колебаний частиц, для применения кри-

рамках метода решеточной динамики в квазигармо-

терия Линдеманна [23,24] при нахождении темпера-

ническом приближении был разработан трехточеч-

туры плавления льдов, а также вычислять энергию

ный потенциал SPC/E_mod-HW на основе модифи-

Гельмгольца на основе колебательных данных.

цированного потенциала SPC/E, который здесь на-

Результаты и обсуждение. Собственные час-

зывается SPC/E_mod. Потенциал SPC/E_mod [17] по-

тоты колебаний, как и собственные вектора коле-

казал наибольшую точность при описании термоди-

баний, являются критически важными для восста-

намических свойств льда и газовых гидратов [18]

новления термодинамических параметров системы

по сравнению с немодифицированным потенциалом

и расчета температуры плавления кристаллической

Письма в ЖЭТФ том 116 вып. 5 - 6

2022

Модифицированная молекулярная модель для описания конденсированных фаз тяжелой воды

315

Рис. 1. (Цветной онлайн) Расчетная плотность колебательных состояний σ (сплошные линии) льда I_h, состоящего

только из молекул H2O (черные линии), только из молекул D2O (синие линии), и смеси 50 % H2O + 50 % D2O (оран-

жевые линии). Пунктирными линиями показаны результаты экспериментальных исследований (а) для обычного и

“тяжелого” льда I_h [25]. Изменение спектра вращательных колебаний (b) при повышении концентрации D2O. Данные

приведены для температуры 270 К и давления 1 бар

решетки, используя критерий Линдеманна. На ри-

сунке 1 приведено сравнение плотности колебатель-

ных состояний обычного льда и льда, образованного

тяжелой водой в сравнении с экспериментальными

данными по аналогичным системам (рис. 1a). Мож-

но видеть, что расчетные данные достаточно точно

[25] описывают положение и ширину щели между

трансляционными и вращательными колебаниями.

Более того, положение основных расчетных пиков

в области до 300 см^-1 совпадает с эксперименталь-

ными данными. Таким образом, можно говорить о

точности разработанного теоретического описания.

Можно видеть, что рассчитанная трансляционная

область спектра от 0 до ∼ 300 см^-1 не претерпева-

Рис. 2. (Цветной онлайн) Температурная зависимость

ет значительных изменений при замене H₂O на D₂O.

объема элементарной ячейки Vu.c. D2O льда Ih при

Основное различие заключается в изменении обла-

атмосферном давлении. Черная линия соответствует

сти вращательных колебаний (рис.1b): при повы-

расчетным данным. Оранжевой линией с квадратны-

шении концентрации D₂O появляются низкочастот-

ми маркерами и синей линией с круглыми маркера-

ми показаны результаты экспериментальных исследо-

ные вращательные колебания в области от ∼ 380

ваний [15] и [26] соответственно

до ∼ 520 cм^-1 из-за значительно возросшего момен-

та инерции молекул. В чистом D₂O льду, где отсут-

ствуют молекулы H₂O, верхняя граница области вы-

сокочастотных вращательных колебаний снижается

рение системы, так и температуру плавления D₂O

с ∼1000 до ∼750см^-1. Таким образом, спектр вра-

льда I_h. На рисунке 2 представлены расчетные объ-

щательных колебаний льда в общих чертах является

емы элементарных ячеек тяжелого льда в сравнении

суперпозицией спектров предельных вариантов кон-

с данными эксперимента при атмосферном давле-

центрации.

нии. Видно, что отличие расчетного значения объе-

Расчет свободной энергии позволяет получить

ма элементарной ячейки от экспериментальных дан-

уравнение состояния системы, из которого можно

ных не превосходит -0.2Å3 при низкой температу-

определить зависимость объема элементарной ячей-

ре и +0.8Å3 при высокой температуре, что являет-

ки от температуры при заданном давлении. В данной

ся достаточно точным для описания термодинамиче-

работе фокус исследования был направлен на мо-

ских свойств льда на основе тяжелой воды, так как

дификацию потенциала SPC/E таким образом, что-

для SPC/E_mod расхождение с экспериментальными

бы он хорошо описывал как температурное расши-

данными больше [18]. При этом потенциал SPC/HW

Письма в ЖЭТФ том 116 вып. 5 - 6

2022

316

В.Р.Белослудов, К.В.Гец, Р.К.Жданов и др.

недооценивает объем элементарной ячейки более чем

термодинамические свойства конденсированных фаз

на 2Å3.

воды с высокой точностью. Показано соответствие

Рассчитанный максимум плотности льда I_h до-

объема элементарной ячейки в зависимости от тем-

стигается при температуре ∼ 60 K, что показано на

пературы с экспериментальными данными. Помимо

рис. 2 и совпадает с экспериментальными значения-

этого, разработанный потенциал позволяет находить

ми 59-60 К [26], что так же указывает на высокую

температуру плавления льда в зависимости от кон-

точность новой модели.

центрации тяжелой воды в нем.

На рисунке 1 было показано, что в H₂O и D₂O

Предложенный подход использует описание ди-

льдах VDOS трансляционных колебаний совпадают.

намических и термодинамических свойств твердого

Мы предполагаем, что в рассматриваемой области

тела в отличие от других подходов, в которых ис-

температуры именно эти колебания имеют наиболь-

пользуется описание свойств жидкой фазы. Поэтому

шие вектора смещения, увеличение которых приво-

представленный потенциал в большей степени пред-

дит к разрушению структуры. Несмотря на то, что

назначен для описания твердых фаз тяжелой воды.

молекулы H₂O и D₂O описываются разными потен-

Работа поддержана Министерством науки и выс-

циалами, величины межмолекулярных расстояний

шего образования Российской Федерации, проект

остаются близкими, поэтому можно предположить,

#121031700321-3.

что критерий Линдеманна H₂O и D₂O льдов дол-

жен совпадать в точке плавления. Взяв за репер-

W. E. Thiessen and A. H. Narten, J. Chem. Phys. 77(5),

ное значение критерия для льда H₂O при темпера-

2656 (1982).

туре 273.15 К, были найдены значения температуры,

S. Herrig, M. Thol, A. H. Harvey, and E. W. Lemmon,

при которых величины критерия Линдеманна льдов

J. Phys. Chem. Ref. Data 47(4), 043102 (2018).

различного состава равны реперному значению. На

J. R. Grigera, J. Chem. Phys. 114(18), 8064 (2001).

рисунке 3 представлена температура плавления льда

J. L. Abascal and C. Vega, J. Chem. Phys. 123(23),

Ih с ее линейной аппроксимацией в зависимости от

234505 (2005).

концентрации D₂O. Можно видеть, что эта зависи-

D. R. Lide, CRC Handbook of Chemistry and Physics: A

мость имеет схожий с линейным характер, учитывая,

Ready-Reference Book of Chemical and Physical Data,

что экспериментальные образцы льда D₂O содержат

CRC, Boca Raton FL, London (2003), v. 84.

небольшую долю H₂O (∼ 0.1 %).

P. V. Enkovich, V. V. Brazhkin, and S. G. Lyapin, JETP

Lett. 110(10), 687 (2019).

P. V. Enkovich, V. V. Brazhkin, S. G. Lyapin, and

S. M. Stishov, JETP 128(2), 207 (2019).

H. C. Urey, F. G. Brickwedde, and G. M. Murphy, Phys.

Rev. 39(1), 164 (1932).

V. C. Chamorro, C. Tempra, and P. Jungwirth, J. Phys.

Chem. B 125(17), 4514 (2021).

10.

J. B. Linse and J. S. Hub, J. Chem. Phys. 154(19),

194501 (2021).

11.

H. J. C.

Berendsen,

J. P. M.

Postma,

W. F. van Gunsteren, and J. Hermans, Interaction

Models for Water in Relation to Protein Hydration,

in Intermolecular forces, Springer, Dordrecht (1981),

p. 331.

Рис. 3. (Цветной онлайн) Зависимость температуры

12.

H. J. C. Berendsen, J. R. Grigera, and T. P. Straatsma,

плавления TM льда Ih от процентной доли молекул

J. Phys. Chem. 91(24), 6269 (1987).

D2O в его структуре при атмосферном давлении. Чер-

13.

M. K. Eseev, V. I. Matveev, and D. N. Makarov, JETP

ная линия аппроксимирует результаты расчета (плю-

Lett. 114(7), 387 (2021).

сы). Круг показывает температуру плавления H2O

льда I_h, звезда показывает температуру плавления

14.

E. O. Wollan, W. L. Davidson, and C. G. Shull, Phys.

Rev. 75(9), 1348 (1949).

D2O льда Ih

15.

K. Röttger, A. Endriss, J. Ihringer, S. Doyle, and

Выводы. В работе предложена трехточечная мо-

W. F. Kuhs, Acta Crystallogr. B: Struct. Sci. 50(6), 644

дель, основанная на модели SPC/E, для описания

(1994).

молекулы D₂O, которая позволяет в рамках мето-

16.

M. Kida, M. Watanabe, Y. Jin, and J. Nagao, J. Chem.

да решеточной динамики описывать динамические и

Eng. Data 60(6), 1939 (2015).

Письма в ЖЭТФ том 116 вып. 5 - 6

2022