Геохимия, 2021, T. 66, № 11, стр. 1050-1066

ВВЕДЕНИЕ

Исследования изотопных эффектов кислорода занимают одно из центральных мест в изотопной геохимии. Подавляющее большинство изотопных геотермометров являются кислородными. Данная статья посвящена совместному расчету химических и кислородных изотопных равновесий на основе программного комплекса GEOCHEQ_Isotope. Основные подходы и методы расчетов кислородных изотопных эффектов в данной работе в целом аналогичны методам и подходам, изложенным нами ранее для изотопных эффектов углерода. Впервые совместное моделирование химических и изотопных равновесий путем минимизации свободной энергии, получившее название: метод изотопно-химической системы (идея принадлежит Б.Н. Рыженко), было выполнено для изотопов серы в работах Д.В. Гричука с соавторами (Банникова и др., 1987; Гричук, 1987, 1988, 2000; Гричук, Лейн, 1991). Этот метод позволяет проводить расчет химических и изотопных равновесий в один этап в отличие от обычных подходов (Ohmoto, 1972), когда на первом этапе проводится расчет химических равновесий, а на втором проводят расчет изотопных равновесий.

Как отмечалось ранее (Mironenko et al., 2018), одной из основных причин, тормозящих распространение совместных расчетов химических и изотопных равновесий, является отсутствие внутренне-согласованных баз термодинамических данных по редким изотопологам. Данное исследование направлено на устранение этого недостатка и создание базы данных для моделирования геохимических процессов с учетом кислородных изотопных эффектов на основе критического анализа имеющейся информации по изотопным равновесиям.

ПОДХОДЫ И МЕТОДЫ

Также, как и в случае углеродных изотопных эффектов (Mironenko et al., 2018), основу для расчетов составляет программный комплекс (Мироненко и др., 2000), использующий базу данных SUPCRT92 (Johnson et al., 1992) c более поздними дополнениями и исправлениями. Минимизация энергии Гиббса в программном комплексе G-EOCHEQ проводится по алгоритму выпуклого симплекса (de Capitani, Brown, 1987). Для учета изотопных эффектов в расширении программного комплекса GEOCHEQ – GEOCHEQ_Isotope вместо химического элемента, для которого рассчитываются изотопные эффекты, введены его изотопы в качестве независимых компонент, a в список веществ добавлены соответствующие изотопологи (Mironenko et al., 2018). В случае кислорода – это изотопы ¹⁶O и ¹⁸O и их изотопологи.

Все расчеты изотопных эффектов в GE-OCHEQ_Isotope проводятся в приближении идеальности смеси изотопов, в котором замещения ядер изотопов в разных позициях не влияют друг на друга, и изотопные эффекты не зависят от кратности изотопного замещения. Для молекул, содержащих несколько атомов, это приближение обычно формулируют в виде правила геометрического среднего (Галимов, 1973, 1974). Для всех элементов за исключением водорода приближение идеальности смеси изотопов выполняется при температурах, превышающих 100 K (Поляков, 1993). Это приближение по умолчанию считается выполненным в подавляющем большинстве геохимических расчетов, если не проводится специальный расчет перекрестных (или “clumped”) изотопных эффектов (Гричук, 2000; Hill et al., 2014, 2020; Schauble, Young, 2020). Как и в случае изотопных эффектов углерода расчеты кислородных изотопных эффектов опираются на формализм β-фактора. Понятие β-фактора было введено Я.М. Варшавским и С.Э. Вайсбергом (Варшавский, Вайсберг, 1957) на основе предложенного ранее (Bigeleisen, Mayer, 1947) редуцированного по числам симметрии отношения статистических сумм изотопологов (the reduced isotopic partition function ratio). Однако, в практику геохимических исследований понятие β-фактора было внедрено Э.М. Галимовым (Галимов, 1973, 1974, 1982; Galimov, 2006). β-фактор вещества, содержащего атомы рассматриваемого химического элемента в эквивалентных положениях удобно выразить в логарифмической форме в терминах энергии Гиббса

где β – β-фактор, G – энергия Гиббса, P – давление, T – температура, R – универсальная газовая постоянная, z – кратность изотопного замещения, m – масса изотопа, “*” помечены величины, относящиеся к тяжелому изотопу. В соответствии с этим равновесный коэффициент разделения изотопов (α_A/B) между двумя веществами A и B равен отношению их β-факторов: или в логарифмической форме для изотопного сдвига (Δ_A/B):

В последнем равенстве (2б) учтено что |β – 1| $ \ll $ 1.

В том случае, когда атомы элемента, подвергающегося изотопному замещению, находятся в неэквивалентных положениях формулы (1) и (2а) требуют уточнения. Поскольку в природных процессах в изотопном обмене участвуют все изотопологи, то необходимо понимание того каким образом β-фактор соединения в целом связан с β-факторами атомов в определенных неэквивалентных положениях. Развитая Э.М. Галимовым (Галимов, 1973, 1974, 1981, 1982; Galimov, 2006) на примере изотопов углерода теория внутримолекулярных изотопных эффектов носит общий характер и имеет непосредственное применение к фракционированию изотопов кислорода, поскольку в геохимических исследованиях часто приходится иметь дело с веществами, которые содержат несколько атомов кислорода, причем не все они находятся в эквивалентных положениях¹¹. В частности Э.М. Галимовым было показано что для замещения редким изотопом β-фактор соединения в целом связан с β-факторами однозамещенных изотопологов (β_i) следующим соотношением (Галимов, 1973, 1981, 1982):

Соотношение (3а) было названо Э.М. Галимовым первым законом аддитивности (Галимов, 1981; 1982)²². Впоследствии было показано, что соотношение (3а) выполняется не только в случае замещения редкими изотопами, но и в том случае, когда β_i не сильно отличаются от единицы (Поляков, 1987; Polyakov, Horita, 2021). Кроме того, логарифмическая форма оказывается незначительно точнее (Polyakov, Horita, 2021):

Поскольку в представленной версии GE-OCHEQ_Isotope внутримолекулярные эффекты кислорода не рассчитываются, использование β-фактора соединения в целом в соответствии с формулой (3б) позволяет избежать излишнего “разбухания” базы данных за счет введения данных по изопологам разной кратности и локализации. Таким образом, при расчете коэффициентов равновесного фракционирования изотопов кислорода в комплексе GEOCHEQ_Isotope достаточно наряду с наиболее распространенным изотопологом, содержащим ¹⁶O, рассмотреть ¹⁸O изотополог, на один атом кислорода которого приходится изменение энергии Гиббса в соответствии с (1), аналогично тому, как это было сделано ранее для фракционирования изотопов углерода ¹³C/¹²C (Mironenko et al., 2018).

Кислородные β-факторы обнаруживают слабую зависимость от давления (Polyakov, Kharlаshina, 1994), поэтому в GEOCHEQ_Isotope учитывается только их температурная зависимость. Температурные зависимости β¹⁸O-факторов представлены в виде полинома по обратным температурам, как это сделано для углеродных β-факторов (Mironenko et al., 2018):

в котором $x = {{{{{{{10}}^{6}}} \mathord{\left/ {\vphantom {{{{{10}}^{6}}} T}} \right. \kern-0em} T}}^{2}}$, a A_i – коэффициенты полинома (4), представленные в табл. 1.

Еще одним добавлением к расчету химических и изотопных равновесий в системах, включающих изотоп ¹⁸O, является изменение формулы для моляльности растворенных в воде компонентов в связи с присутствием системе изотополога H₂¹⁸O:

где m_i – моляльность i-го компонента, n_i – количество молей компонента 55.51 – количество молей в 1 кг воды.

β¹⁸O-ФАКТОРЫ В ТЕРМОДИНАМИЧЕСКОЙ БАЗЕ ДАННЫХ GEOCHEQ_ISOTOPE

Также, как и в случае углеродных изотопных факторов (Mironenko et al., 2018), точность и обоснованность расчетов по кислородным изотопным равновесиям определяется точностью данных по β¹⁸O-факторам, которые используются в программном комплексе GEOCHEQ_Isotope. Для оценки кислородных равновесных изотопных факторов фракционирования применяют как экспериментальные, так и теоретические подходы. Для экспериментального определения кислородных равновесных изотопных факторов фракционирования наряду с методом частичного изотопного обмена Нортропа-Клейтона (Northrop, Clayton, 1966) используют также трехизотопный метод (Matsuhisa et al., 1978). Для калибровки кислородных изотопных геотермометров используют и метод эмпирической калибровки (Clayton et al., 1968; Устинов и др., 1978; Macey, Harris, 2006). Обзор этих методов с оценкой их преимуществ и недостатков дан в (Chacko et al., 2001). Теоретические расчеты проводят как с помощью модельных колебательных спектров (Urey, 1947; Bigeleisen, Mayer, 1947; Kieffer, 1982) и межатомных потенциалов (Patel et al., 1991), так и из “первых принципов” на основе метода функционала плотности (Schauble, 2004; Blanchard et al., 2009; Hill et al., 2014; Крылов и др., 2017; Schauble, Young, 2021). В последние годы метод расчета из “первых принципов” стал наиболее употребительным. Кислородные β-факторы многих минералов рассчитаны с помощью полуэмпирического метода инкрементов (Schütze H., 1980; Zheng, 1991; 1999), который представляет собой распространение метода изотопических чисел связей Э.М. Галимова (Галимов, 1973, 1974) на кристаллические тела. Однако к этим результатам надо относится с осторожностью, поскольку они часто противоречат экспериментальным данным и теоретическим расчетам (Horita, Clayton, 2007).

СОВМЕСТНЫЙ РАСЧЕТ РАВНОВЕСНЫХ КИСЛОРОДНЫХ И УГЛЕРОДНЫХ ИЗОТОПНЫХ ЭФФЕКТОВ И ХИМИЧЕСКОГО СОСТАВА В ГИДРОТЕРМАЛЬНОЙ СИСТЕМЕ

В настоящее время GEOCHEQ_Isotope позволяет моделировать равновесные химические процессы с учетом кислородных и углеродных изотопных эффектов. В качестве примеров подобного моделирования мы выбрали систему, аналогичную рассмотренной ранее (Mironenko et al., 2018). В целом методика расчета химических равновесий с учетом как углеродных, так и кислородных изотопных равновесий совпадает с той, которая применялась к моделям, включающим оценку только углеродных равновесных изотопных эффектов (Mironenko et al., 2018; Мироненко и др., 2021; Сидкина и др., 2019). Естественным изменением явилось добавление изотопа ¹⁸O в качестве независимого компонента и ¹⁸O-изотопологов в список веществ, для которых проводится минимизация энергии Гиббса. Кроме того, моляльность компонентов водного раствора рассчитывали по формуле (5).

В качестве примера использования G-EOCHEQ_Isotope мы рассчитали химические и изотопные равновесия в углекислой гидротермальной системе валового состава “1 моль СаО + + 1 кг Н₂О + 1.4 моля СО₂” для температур 15, 25, 40°С и давления 40 бар в зависимости от водородного показателя pH кислотности. Соотношения ¹⁸О/¹⁶О и ¹³С/¹²С в системе отвечали природным распространенностям этих изотопов. Рассматривались ¹⁶O-, ¹⁸O-, ¹²C- и ¹³C-изотопологи кальцита, воды, СО₂(aq), ${\text{HCO}}_{3}^{ - }$, ${\text{CO}}_{3}^{{2 - }}$, CaCO₃(aq), Ca${\text{HCO}}_{3}^{ - }$, CO₂(газ), СН₄(газ), H₂O(газ). Для подобной системы имеется экспериментально-теоретическое исследование зависимости состава растворенного неорганического (карбонатного) углерода и равновесного фракционирования изотопов кислорода между растворенными карбонатными частицами и водой от pH (Beck et al., 2005). Между этим экспериментальным исследованием и карбонатной системой, которая рассмотрена в настоящей статье, есть очень незначительные отличия. При моделировании с использованием GEOCHEQ_Isotope, мы рассматривали возможность появления растворенных частиц CaCO₃(aq), Ca${\text{HCO}}_{3}^{ - }$, в то время, как в работе (Beck et al., 2005) эти частицы не принимали во внимание. С другой стороны, при теоретическом моделировании в работе (Beck et al., 2005) учитывали присутствие в растворе частиц H₂CO₃. Следует подчеркнуть, что концентрация H₂CO₃ по результатам моделирования составляла ~0.2–0.3% в общем количестве растворенного карбонатного углерода при pH < 5 и резко уменьшалась при понижении кислотности. Результаты моделирования состава карбонатных частиц средствами GEOCHEQ_Isotope представлены на рис. 7. Наблюдается хорошее согласие между результатами нашего моделирования с предшествующими результатами (Beck et al., 2005).

Доминирование той или иной формы растворенного карбонатного углерода в зависимости от pH с учетом разницы в значениях β¹⁸O- и β¹³C-факторов растворенных карбонатных частиц (табл. 1, рис. 1) определяет зависимость равновесных коэффициента фракционирования изотопов кислорода и углерода между общим растворенным неорганическим углеродом и сосуществующими фазами. Примеры таких зависимостей по результатам расчетов с использованием GEOCHEQ_Isotope показаны на рис. 8. Для изотопов кислорода имеются надежные экспериментальные данные (Beck et al., 2005) по равновесному фракционированию растворенного карбонатного углерода и воды в интервале температур 15–40°С. Результаты нашего моделирования обнаруживают хорошее согласие с этими экспериментальными данными (рис. 8). Уменьшение обогащения растворенных карбонатов тяжелым изотопом кислорода относительно воды с ростом pH (рис. 8) объясняется последовательной сменой доминирующей карбонатной формы в ряду CO₂(aq) → ${\text{HCO}}_{3}^{ - }$ → ${\text{CO}}_{3}^{{2 - }}$ наряду с тем, что равновесные факторы фракционирования этих карбонатных форм с водой (а значит и их β¹⁸O-факторы) уменьшаются в той же последовательности (рис. 1). Для изотопов углерода имеет место обратная ситуация: β¹³С-факторы указанных карбонатных форм имеют обратный порядок (Mironenko et al., 2021) и, поэтому, обеднение растворенного карбонатного углерода тяжелым изотопом ¹³С уменьшается с ростом pH. Причем углеродный изотопный сдвиг между растворенныч карбонатным углеродом и кальцитом при высоких pH претерпевает инверсию, поскольку значения β¹³С-фактора для доминирующей карбонатной формы ${\text{CO}}_{3}^{{2 - }}$ больше, чем для кальцита (рис. 8).

ВЫВОДЫ

Представлено расширение базы данных и программного комплекса GEOCHEQ_Isotope, предназначенного для одновременного моделирования химических и изотопных равновесий методом минимизации свободной энергии Гиббса, на расчет изотопных равновесий кислорода. Также, как и в случае изотопов углерода, представленная версия основана на формализме β-фактора в приближении. идеальной смеси изотопов.

Представленная версия GEOCHEQ_Isotope использует базу данных GEOCHEQ (Мироненко и др., 2000) для описания неизотопных термодинамических свойств веществ. Энергия Гиббса для изотопологов рассчитывается по известным значениям β‑факторов в соответствии с уравнением (1). При этом учитывается аддитивная связь между lnβ соединения в целом и lnβ_i однозамещенных изотопологов.

Проведена адаптация расчетных алгоритмов на случай использования изотопов кислорода и повышения точности расчетов. В частности, в уравнении для расчета моляльности учтено присутствие изотополога H₂¹⁸O.

Проведен критический анализ, оптимизация и согласование имеющейся информации по равновесным факторам фракционирования изотопов кислорода. Температурные зависимости углеродных β-факторов унифицированы в виде полиномов (4) по степеням обратной температуры (табл. 1).

Применение комплекса GEOCHEQ_Isotope к совместному расчету химических кислородных и углеродных изотопных равновесий протестировано на моделировании зависимости состава растворенных карбонатных частиц и коэффициента равновесного фракционирования ¹⁸O/¹⁶O между растворенным карбонатным углеродом и водой от pH раствора. Получено хорошее согласие с экспериментальными результатами (Beck et al., 2005). С увеличением pH обогащение растворенного карбонатного углерода тяжелым изотопом ¹⁸O уменьшается за счет последовательной смены доминирующей карбонатной формы в ряду CO₂(aq) → ${\text{HCO}}_{3}^{ - }$ → → ${\text{CO}}_{3}^{{2 - }}$, характеризующейся меньшим значением β¹⁸O-фактора (${{\beta }^{{{\text{18}}}}}{{{\text{O}}}_{{{\text{C}}{{{\text{O}}}_{2}}}}} > {{\beta }^{{18}}}{{{\text{O}}}_{{{\text{HCO}}_{3}^{ - }}}} > ~~{{\beta }^{{18}}}{{{\text{O}}}_{{{\text{CO}}_{{\text{3}}}^{{{\text{2 - }}}}}}}$). Для углерода рост pH приводит к увеличению коэффициента фракционирования растворенный карбонатный углерод–кальцит, поскольку ${{\beta }^{{{\text{13}}}}}{{{\text{C}}}_{{{\text{CO}}}}}_{{_{{\text{2}}}}} < {{\beta }^{{13}}}{{{\text{C}}}_{{{\text{HCO}}_{3}^{ - }}}} < ~~{{\beta }^{{13}}}{{{\text{C}}}_{{{\text{CO}}_{{\text{3}}}^{{{\text{2 - }}}}}}}$. При этом при высоких значениях рН (~11) происходит инверсия равновесного изотопного сдвига и растворенный неорганический углерод обогащен изотопом ¹³С относительно кальцита.

Исследование выполнено в рамках госзадания ГЕОХИ РАН (№ 0137-2019-0015) и частично поддержано РФФИ грант № 19-05-00865a.