УДК 553.212, 552.111
1 1 Е.В. Коптев-Дворников , Д.А. Бычков
РАЗРАБОТКА ЛИКВИДУСНОГО ТЕРМОБАРОМЕТРА ДЛЯ МОДЕЛИРОВАНИЯ РАВНОВЕСИЯ ОЛИВИН-РАСПЛАВ
ФГБОУ ВО «Московский государственный университет имени М.В. Ломоносова», геологический факультет, 119991, Москва, ГСП-1, Ленинские Горы, 1
Lomonosov Moscow State University, Faculty of Geology, 119991, Moscow, GSP-1, Leninskiye Gory, 1
Получена система уравнений ликвидусного термобарометра оливин—силикатный расплав путем обработки методами многомерной статистики выборки из 772 экспериментальных равновесий кристаллов оливина с базитовыми расплавами. Термобарометр с низкой погрешностью воспроизводит экспериментальные составы в широком диапазоне состава базитов (от коматиитов до дацитов), температуры — от 1040 до 1500 °С и давления — от 1 бара до 30 кбар. Тестирование термобарометра продемонстрировало, что отклонения расчетной температуры ликвидуса от экспериментальной в большей части температурного интервала не превышают ±3 °С.
Ключевые слова: оливин, силикатный расплав, уравнение, термобарометр.
A system of equations of the liquidus thermobarometer of olivine — silicate melt was obtained by processing the sample of 772 experimental equilibria of olivines with basic melts using methods of multidimensional statistics. Equations reproduce with small error experimental data in a wide range of basite compositions (from komatiites to dacites), temperatures from 1040 to 1500 °С, pressures from 1 bar up to 30 kbar. Thermobarometer testing demonstrated that the deviations of the calculated liquidus temperature from the experimental one in most of the temperature range do not exceed ±3 °C.
Key words: olivine, silicate melt, equation, thermobarometer.
Введение. Крупные расслоенные мафит-уль-трамафитовые плутоны служат важными источниками меди, никеля, элементов платиновой группы, хрома, железа, титана, ванадия и др. Кроме придонных медно-никелевых руд, все остальные представлены выдержанными по простиранию горизонтами сплошных (хромитовые и титаномаг-нетитовые) или шлировидных (малосульфидные платинометалльные) руд, пространственно ассоциирующих с ритмически расслоенными пачками пород.
Для решения проблемы ритмической рассло-енности мафит-ультрамафитовых плутонов была предложена понятийная многослойно-суспензионная модель [Бычкова, Коптев-Дворников, 2004]. Разработанная ранее программа КриМинал (программа расчета равновесной кристаллизации (Кри) силикатных систем с использованием в расчетах суммы минеральных миналов (Минал)) [Бычков, Коптев-Дворников, 2005; ВуеЬкоу, Кор1еу-Dvornikov, 2014] — прообраз термодинамического блока для численной многослойно-суспензионной модели динамики внутрикамерной дифференциации, которая в свою очередь представляет дальнейшее развитие конвекционно-кумуляционной модели [Коптев-Дворников и др., 1979; Френкель, Ярошевский, 1978], реализованной затем
в динамическом блоке программы КОМАГМАТ [Арискин, Бармина, 2000; Френкель и др., 1988; Френкель, 1995].
Для термодинамического блока необходимо вывести уравнения для решения задачи термодинамического равновесия расплавов с породообразующими силикатными минералами и акцессорными в рассеянном состоянии фазами, скопления которых образуют экономически важные руды.
Обоснование принятого вида системы уравнений ликвидусного термобарометра. В работах М.Я. Френкеля [Френкель и др., 1988; Френкель, 1995] четко сформулировано существование двух подходов к построению алгоритмов решения задачи термодинамического равновесия: «Программы первого типа строятся как процедуры поиска минимума поверхности свободной энергии смеси на многограннике ограничений (условий баланса вещества). Программы второго типа реализуют численные схемы решения систем нелинейных алгебраических уравнений закона действующих масс, то есть уравнений, определяющих положение минимума поверхности термодинамического потенциала».
Первый подход использован в Melts. Несмотря на длительные усилия и значительные достижения коллектива исследователей [Asimow, Ghiorso,
1 Московский государственный университет имени М.В. Ломоносова, геологический факультет, кафедра геохимии, ст. науч. с., доцент, канд. геол.-минер. н.; e-mail: ekoptev@geol.msu.ru
2 Московский государственный университет имени М.В. Ломоносова, геологический факультет, кафедра геохимии, науч. с.; e-mail: krok@geol.msu.ru
1998], этот подход встречает трудности, связанные с отсутствием достаточно полных и надежных данных о соотношениях между значениями активности и концентрации компонентов силикатных расплавов.
Второй подход, реализованный, в частности, в программном комплексе КОМАГМАТ, имеет то преимущество, что позволяет использовать в качестве термодинамической основы не очень систематические и полные данные закалочных опытов по температурным зависимостям состава сосуществующих минералов и расплава. Решение уравнений закона действующих масс, связывающих между собой равновесные значения активности исходных веществ и продуктов реакции, происходит совместно с решением уравнений материального баланса. Таким образом, решение термодинамической задачи этим методом (часто называемым методом констант) получило распространение в моделях, используемых в магматической петрологии, например, в программах семейства MAGMOD-KOMAGMAT [Арискин, Бармина, 2000; Френкель и др., 1988; Френкель, 1995] и в программе Рйго^ ^апушИ^ку, Р1есИоу, 2011]. Такой метод позволяет использовать полуэмпирические уравнения равновесия, полученные статистической обработкой экспериментальных данных, но сохраняя черты процедуры минимизации свободной энергии. Разрабатываемый нами алгоритм поиска равновесия в силикатной системе использует второй подход.
История разработки основ термодинамического описания равновесий минерал—силикатный расплав достаточно полно изложена в монографии [Арискин, Бармина, 2000]. Первые уравнения, описывающие эти равновесия, получены в 1970-х гг. П.Л. Редером, Е. Эмсли, М.Дж. Дрейком, Р.Л. Нильсеном и другими (см. [Арискин, Бармина, 2000, с. 30]) путем статистической обработки закалочных экспериментов. Эти уравнения описывают зависимости констант равновесия реакций образования из расплава миналов как компонентов твердых растворов и имеют вид линейной зависимости от обратной температуры:
ькт=лц т+вт, (1)
где К ]т — константа равновесия реакции кристаллизации минала (т) твердого раствора минерала (/); Т — абсолютная температура; Ат и Вт/ — коэффициент и константа, получаемые в результате статистической обработки, соответственно. Следовательно, равновесие минерал—расплав определяется системой из т таких уравнений. Эти системы уравнений получили названия геотермометров минерал—расплав.
В то же время, для гетерофазной реакции закон действующих масс приобретает вид
кт=тп (а V, (2)
где aj — активность минала m в минерале j (обычно принимается равной Xj — мольной доле минала m в минерале j); ()''m — активность в расплаве компонента i (вычисляется согласно двухрешеточной модели расплава Р.Л. Нильсена и М.А. Дунгана [Nielsen, Dungan, 1983]) в степени стехиометрического коэффициента vi m в реакции образования минала m. Для членов экспериментальной выборки, в которой величины мольных долей минала m в минерале j поставлены в соответствие с температурой и составом сосуществующих расплавов, по уравнению (2) рассчитываются значения Kj, подставляются в уравнение (1) и методом линейной регрессии вычисляются параметры Aj и Bj в уравнении (1). Именно таков общепринятый сегодня метод определения этих параметров.
После подстановки уравнения (2) в левую часть уравнения (1) получаем выражение для логарифма содержания минала в минерале в зависимости от температуры и состава расплава:
ln XJm = Aj/T + BJm + Yv. ln ct.. (3)
m ml m ^^ ', m ' v /
Потенцирование этой величины дает значение концентрации минала в твердой фазе. Принимая во внимание термодинамическую природу зависимости фазовых равновесий от давления, А.А. Арискин с соавторами [Арискин, Бармина, 1990; Арискин, Цехоня, Френкель, 1991] учли влияние давления на константу равновесия реакции кристаллизации в следующем соотношении:
inK m = (Aj+ртл/т+вт, (4)
где ßm — коэффициент, подобранный для каждого минала графически (см. [Арискин, Бармина, 2000, с. 111]); P — давление, кбар. В настоящее время для систем этих уравнений используется, на наш взгляд, более точное название, а именно ликви-дусные термобарометры.
Уравнение для логарифма мольной доли ми-нала приобретает вид
inxm = (Am+ßmp)/T+Bm+^vm in c. (5)
Системы уравнений вида (5) составляют термодинамическую основу программного комплекса КОМАГМАТ, позволяющего рассчитывать составы твердых фаз, равновесных с расплавом.
Однако в таком виде термобарометры, как правило, удовлетворительно воспроизводят состав твердых фаз в ограниченном диапазоне состава расплавов, причем при обработке широких по составу выборок на графиках lnK—1/Т иногда возникают субпараллельные цепочки точек, что навело на мысль о зависимости константы в уравнениях (1)—(4) от состава расплава и летучести кислорода.
Для термодинамического блока новой модели динамики кристаллизационной дифференциации целесообразно вывести такие уравнения ликви-
дусных термометров, которые были бы пригодны для использования в широком диапазоне составов базитовых систем. Перебрав ряд подходов к учету состава системы при расчете константы равновесия, авторы статьи [Коптев-Дворников, Бычков, 2007] остановились на следующем виде уравнения ликвидусного термобарометра
lnK i =
Aj + VjP
T
+b¿+ci t +
+ Dm lg fO2 + El ln
'Al4 Si
+ FiW + £ Jl
(6)
Bj + CjT + Dj lg fO2 + Ej ln
'AP Si
+
нетита [Арьяева и др., 2018], ильменита [Черных, 2017] и жидкого моносульфида железа [Коптев-Дворников и др., 2012]. Помимо разработанных термобарометров, в результате этих работ были получены два методических вывода относительно способа нахождения коэффициентов и константы в уравнении (7) и оценки точности предсказания составов твердых фаз и значений температуры ликвидуса с помощью термобарометров.
Из уравнения (7) получаем выражение для вычисления содержания минала т в минерале j^.
xm=exp
где fO2 — летучесть кислорода, бар; W = ln((Na + K) Al/Si2); щ — мольная доля i-го компонента расплава; n — число учитываемых компонентов расплава. Параметры Al/Si и W (вычисляются с использованием атомных количеств элементов) предложены в работе [Арискин, Бар-мина, 2000] для уточнения оливиновых и плаги-оклазовых термобарометров соответственно; Cj, Dj, Ej, Fj, Jj ¡ — коэффициенты при соответствующих переменных. В основу нахождения этих коэффициентов положен принцип использования методов многомерной статистики. Из вида известных физико-химических уравнений следует термодинамический смысл этих параметров:
Aj ~ ДН1 /R Вj ~ - Д Vj /R
(A + ByP) /T + Bj + Dj lg fO2 +
y j ijj j j © J 2
Ej ln
Si
+ FJW + У JJ X¡ +tv . ln a'
j ^^ j, i . ^^ j,. .
i=1
(8)
+ + %31,Х, |~АБЦR,
где АН^, А Кт и АБ^ — энтальпийный, объемный и энтропийный эффекты реакции образования минала т минерала7 из расплава, R — универсальная газовая постоянная.
Из уравнения (6) следует выражение для логарифма мольной доли минала в минерале.
1П х]т = №+№ )Д + вт+стт+lgо +.
+ В1 1П(А1/ + % Г + % X + % V,,,, 1П а. (7)
!=1 ! = 1
Значения коэффициентов и константы в этих уравнениях могут быть найдены путем статистической обработки (оптимизации) соответствующей выборки методом многомерной линейной регрессии.
К настоящему времени нами опубликованы уравнения, позволяющие моделировать равновесия с силикатными расплавами в широком диапазоне состава, температуры, давления и летучести кислорода хромшпинелидов [Арьяева и др., 2016], маг-
Показано [Коптев-Дворников и др., 2012], что нахождение коэффициентов и константы в термобарометрах не традиционной оптимизацией линейного относительно логарифма константы уравнения вида (7), а экспоненциального уравнения зависимости содержания миналов от состава силикатного расплава, температуры, давления и летучести кислорода вида (8) приводит к существенному улучшению воспроизведения экспериментальных данных. Причина этого понятна — при оптимизации уравнения (8) минимизируются разности не между логарифмами содержания ми-налов, а разности между собственно содержанием миналов. Нахождение коэффициентов при переменных и константы в уравнении (8) выполняется путем минимизации суммы квадратов разностей между значениями расчетного и экспериментального содержания миналов с использованием надстройки «поиск решения» в программе Ехсе1. В статистике принято называть такие разности между наблюдаемыми и расчетными величинами «остатками», далее мы будем использовать этот термин.
Процесс оптимизации, как правило, не одноактный. Результаты оптимизации считаются удовлетворительными, если угловой коэффициент в уравнении линейной регрессии близок к 1, а свободный член — к 0, остатки распределены по нормальному закону. Те экспериментальные точки, которые нарушали эти требования, удалялись из выборки, и оптимизацию выполняли снова и снова до тех пор, пока эти условия не оказывались выполненными.
Нулевое значение свободного члена в уравнении регрессии и близость линейного члена к 1 свидетельствуют о том, что отсутствует существенная систематическая ошибка. Мера величины случайной ошибки (разброса значений) — критерий который показывает долю объясненной дисперсии наблюдаемой величины. Близость гистограммы остатков к нормальному распределению свиде-
i=1
тельствует о случайном характере распределения отклонений.
Анализ параметров использованных нами выборок показал, что ни одна из переменных в уравнениях вида (8) не распределена по нормальному или близкому к нему закону и ни с одной из переменных по отдельности остатки не коррелируют. Следовательно, различия между экспериментальными и расчетными величинами можно рассматривать как случайные отклонения, вызванные аналитическими погрешностями, погрешностями измерения экспериментальных параметров и т.п. Это дает нам право применять аппарат математической статистики для анализа и оценки достоверности полученного решения.
Мы предлагаем систему уравнений оливиново-го ликвидусного термобарометра, алгоритмически совместимую с термобарометрами для других фаз.
Вывод уравнений ликвидусного термобарометра оливин—расплав. Традиционно состав оливинов пересчитывают на два минала — форстеритовый MgSi0;5O2 ^о) и фаялитовый FeSi0;5O2 ^а). Однако как в природных образцах, так и в экспериментальных оливинах, микрозондовый анализ обнаруживает присутствие MnO, CaO и Сг^3 в количестве до нескольких процентов. В связи с этим мы пересчитывали составы оливинов на 5 миналов — форстеритовый MgSi0;5O2 ^о), фаялитовый FeSi0;5O2 (Fa), тефроитовый (Тер), ларнитовый CaSi0;5O2 ^гп), а также хром-оливиновый CгSi0;5O2 (CгO1), причем в последнем случае постулируется, что весь Сг в расплаве присутствует в двухвалентной форме.
Образование вышеперечисленных миналов из расплава происходит в результате следующих гетерофазовых реакций:
Mg01+0,5810'.,= =Mg8iо,502 (Ра),
БеО'+О^Ю, =Ре810,50, (Ра),
Мп0'+0,5810,= =Mn810,5O2 (Тер),
Са01+0,5810,= =Са810,50, (Ргп),
ХРо = еХР
+ ЕРо 1п
А1
( + Р РоР)/ т + ВРо + БРо 1н /0, +
п
+ X + 1п + 0,51п а
I
а8ю,
ХРа = еХР
(Л- + )/ т + ВРа + ЮРа 1н /0,
+ ЕРа 1п
81
+ X 3Ра,X + 1п а^0 + о, 51п а8 „
хтр = ехР
+ ЕТер 1п
I1
, (10)
Х1гп = еХР
+ ЕЫп 1п
ХСт01 = еХр
+ ЕСг01 1п
(( + кр)/ т + вТер + оТер 1н ¡0,
п
+ X 3'тер,,Iх'1 + 1п -МпО + 0,5511п а
1=1
( + Р шР)/ т + ВЬгп + Вш 1н /02 +
п
X 3ьгп,Х\ + 1п аСао + 0,51п а;
1=1
(ЛСг01 + $СЮ1Р)/ т + В Сг01 + Осг01 1н /0, + п
+ Х 3сгоих1 + 1п аСго15 + 0,51п а
А1
I1
I1
СЮ1+0,5810, (Сг01).
Из вида уравнения (8) и реакций (9) следуют выражения для расчета содержания оливиновых миналов (10).
Как упомянуто выше, коэффициенты при переменных и константы в уравнениях (10) находятся путем обработки выборки, включающей результаты закалочных опытов, методами многомерной статистики. Оптимизация заключается в минимизации суммы квадратов остатков миналов с использованием надстройки «поиск решения» в программе Ехсе1.
Для создания выборки мы использовали версию базы данных ИНФОРЭКС [Арискин и др., 1997], включающую результаты 14 473 закалочных экспериментов, извлеченных из 419 статей, опубликованных с 1977 по 2010 г. Из них оливин и закаленный расплав сосуществуют в результатах приблизительно 6100 опытов. Однако далеко не все из них пригодны для определения коэффициентов в уравнениях ликвидусных термобарометров. Условием включения экспериментов в выборку было выполнение их в «сухих» условиях, наличие данных о составе сосуществующих расплава и оливина, а также летучести кислорода, при которых проводился опыт (температура и давление определены для всех экспериментов в базе данных). Из опытов, выполненных при давлении 1 бар, отбирались те, в которых в качестве контейнера использовалась платиновая петля, а их длительность составляла не менее 48 ч. Среди высокобарных экспериментов крайне мало опытов, проведенных в условиях контролируемой летучести кислорода, поэтому для них критерий длительности был значительно смягчен — отобраны эксперименты продолжительностью более 12 ч. Кроме того, вид контейнера не служил критерием отбора. Даже при этих условиях с контролируемой летучестью кислорода при повышенных значениях давления оказалось всего 74 эксперимента. Мы предположили, что в том случае, если среди них попались
¡=1
¡=1
1=1
Рис. 1. Корреляция между экспериментальными и рассчитанными в КриМинал значениями температуры равновесия оливин-расплав (а) и гистограмма распределения остатков температуры (б), 568 экспериментальных точек. Сплошная линия — линейный тренд, штриховая линия — линия равных значений (практически совпадают). На гистограмме сплошная кривая — линия аппроксимации нормальным распределением. Заштрихованная область — доверительный интервал аппроксимации на 95%-ном
уровне надежности
слишком неравновесные опыты, они отсеются при дальнейшей статистической обработке. Поскольку разрабатываемая модель направлена на моделирование кристаллизации и динамики формирования интрузивов мафит-ультрамафитового ряда, в выборку не вошли эксперименты, в которых расплавы отвечали по составу щелочным породам, содержали ТЮ2 > 6,5 масс.% и (№20+К20) > 9 масс.%, т.е. находились вне диапазона состава от коматиитов до дацитов. Всего в выборке оказались результаты 651 эксперимента.
В результате последовательных оптимизаций и удаления тех опытов, для которых остатки ми-налов нарушали нормальный закон распределения этого параметра, размер выборки сократился до 568 экспериментов, из которых 53 выполнены при давлении от 1 до 19 кбар. Представление о качестве термобарометра можно получить на примере соответствия экспериментальных и рассчитанных значений температуры (рис. 1). На этом можно было бы остановиться, посчитав разработку термобарометра выполненной.
Однако обращает внимание небольшое число экспериментов в выборке, выполненных при повышенном давлении. В то же время оказалось, что значения коэффициентов при 1е/02 для форстери-тового и фаялитового миналов весьма близки по абсолютной величине и противоположны по знаку (^0= 0,02017; DFa= —0,02071). Поэтому возникла идея попытаться расширить число высокобарных экспериментов в выборке, восстановив в опытах с неконтролируемой летучестью ту «эффективную» летучесть кислорода, при которой составы распла-
вов окажутся в равновесии с составами оливинов. Принцип подбора летучести поясняет рис. 2.
Напомним, что принятый критерий равновесия в программе КриМинал — равенство единице (или 100%) суммы мольных долей миналов [Быч-
12
-8 Н-1-1-1
-12 -8 -4 О
№
Рис. 2. Зависимость разницы между экспериментальным и рассчитанным по уравнению (10) содержанием миналов от подбираемой летучести кислорода; Foэксп, Faэксп — содержание форстеритового и фаялитового миналов в экспериментальных оливинах соответственно; Foс, Faс — рассчитанные по термобарометру значения содержания миналов в оливинах; 100-Е — отклонение рассчитанной суммы миналов от 100%
Таблица 1
Список экспериментов в выборке для вывода системы уравнений термобарометра оливин—расплав
ц/ц Номер публикации, К* Номера опытов в публикации, п** Диапазон условий
Т, °С Р, кбар
1 9 5, 15, 18-21 1203-1252 0,001 13,4___ 11,59
2 20 1-13, 19-28, 33-41, 48, 50, 52, 53, 56-59 1040-1245 0,001 14,4... 12
3 26 17, 18 1106, 1133 0,001 -9,43, -9,09
4 27 1-9, 12-15, 17, 20, 22-24, 26-30, 32, 33, 40 1076-1228 0,001 —10,11_—8,15
5 28 1, 3-11, 13, 15, 16, 18-22, 24-26, 29, 31-38, 40, 44-47, 50 1137-1228 0,001 —9,14_—8,19
6 29 5, 6, 13-15, 21 1091-1125 0,001-8 —8,92_—7,2
7 30 3-5, 7, 17, 19, 21, 22 1153-1195 0,001 -9,13_-8,53
8 31 11-14, 21, 29, 40, 44-46, 53, 54, 58, 62, 65, 69-74, 77, 80, 83, 88, 90-93, 96-99, 101 1092-1235 0,001 -9,95_-7,92
9 39 40, 48 1275, 1300 10-15 -11,35, -10,2
10 49 11, 13 1230, 1325 10 -10,1, -6,3
11 54 1-3, 7, 10, 11, 13, 21-23, 28 1100-1250 0,001-8 —11,3_—4,75
12 63 2 1225 0,001 —9,80
13 67 1-11, 13-17, 21-29, 33-39, 45 1067-1200 0,001 —9,89_—5,93
14 69 1, 2, 6, 7, 8 1173-1187 0,001 -8,94_-8,41
15 70 7-11, 17, 18, 20, 24-26, 30, 31, 36, 43, 44 1074-1149 0,001 -9,3_-8,3
16 71 16-22, 28 1100-1149 0,001 -9,51_-6,05
17 75 10, 11 1160-1164 1 -12,77, -12,72
18 76 1, 8, 9, 14, 15, 29-33 1150-1350 0,001 —11,4_—6,77
19 77 2, 8, 9, 10 1215-1270 0,001 -8,3_-7,7
20 81 11-20, 27-34 1066-1193 0,001 -9,96_-8,38
21 82 54 1200 5 -6,95
22 83 1-7, 11, 12, 16-19, 23-26, 28-30, 33-35, 41, 42, 51 1090-1192 0,001 -12,87_-9,16
23 84 11-14 1126-1152 2 —8,1_—7,15
24 85 7, 16, 17, 29, 30 1260-1360 11-16 —9,2_—6,2
25 86 7, 31 1200, 1300 10, 15 —7,6, -6,5
26 87 5, 12, 18, 19, 24, 33, 54-56 1120-1200 0-15 -6,54_-2,6
27 88 23 1160 0,001 -8,04
28 90 1-19, 21-27 1160-1220 0,001 —11,3_—9,89
29 92 6-11 1108-1174 0,001 —9,7_—8,81
30 98 2-4 1170-1190 0,001 —8,94_—8,71
31 100 1-10, 16-18, 20, 21, 23-25, 29, 30, 32, 38, 43, 47, 48 1089-1248 0,001 —10,01_—7,76
32 158 1, 3, 4, 8, 14, 16, 18, 19 1270-1390 10 —8,65_—4,1
33 166 7-10, 12, 15, 18, 19 1170-1275 0,001 —13,32_—12,01
34 167 4, 6, 11 1200 0,1 —15,8_—14,2
35 169 5-12 1109-1137 0,001 —9,49_—8,23
36 176 14, 17, 19, 21-23, 26, 27, 29-33 1072 0,001 —11,78_—9,98
37 177 12-16, 18-20, 25-28, 32, 33 1065-1135 0,001-7 —13,52_—10,21
38 178 8-11, 18-20 1140-1203 0,001 —9,31_—8,38
39 189 4, 14-17 1116-1125 0,001 —9,31_—9,19
40 193 2, 3, 7-13, 20-27, 29-38 1050-1130 0,001 -12,44_-8,08
41 194 15-22 1149-1212 0,001 —8,90_—8,17
42 197 11-13 1149-1170 0,001 —8,90_—8,65
43 236 3-6, 8-10, 12-15, 17-19, 21-24, 29-32 1110-1225 0,001 —9,76_—8,44
44 247 13, 14, 16 1386-1416 19 -7,56_-7,38
45 255 1, 4, 5, 7-9, 16, 27, 28, 31, 32, 35-38 1100-1200 0,001 -9,51_-8,30
46 262 3, 11 1250, 1310 10 -7, -4,8
47 272 4-7, 9, 11-13, 15-17 1269-1388 15 —9,7_—5
48 273 1-3, 5, 8-10, 12, 13 1300-1500 15-30 —13,6_—5
49 279 9, 15, 16 1330-1370 12-16 —5,7_—5,2
50 280 11-14, 23-28 1111-1158 0,001 —9,58_—8,81
51 281 3-8, 14-17, 19 1103-1158 0,001 —9,58_—8,81
52 288 9 1200 6 —10,66
Окончание табл. 1
ц/ц Номер публикации, К* Номера опытов в публикации, п** Диапазон условий
Т, °С Р, кбар
53 291 29-35, 38, 39, 42, 45, 47 1325-1400 10-15 -5,5...-4
54 294 8, 10, 11 1237-1275 0,001 -7,90...-7,50
55 298 3-13, 15-26 1260-1390 10 -10,35.-3,6
56 299 7, 8, 10, 12, 15, 17, 27, 30, 32 1314-1400 10 -9,75.-4,4
57 300 16, 23 1301, 1319 10 -6, -5
58 302 1, 2 1360-1363 12 -7,95.-6,65
59 304 18 1350 10 -6,40
60 305 3, 5, 6, 18, 22, 23, 33, 34, 38, 39 1114-1184 0,001 -9,33.-7,75
61 309 2-7, 11-15, 17-22, 24-29 1270-1390 10 -11,5.-3,5
62 311 3, 10-13 1260-1295 10 -6,5.-5,2
63 312 14-18, 20, 28-30, 41, 43, 44, 46, 47, 49, 64, 69, 77 1270-1400 15 -14.-4,4
64 322 1-8, 10-12, 14, 16-19 1260-1405 7-17 -8,9.-3,2
65 323 1- 17 1270-1390 10 -9,95.-4,5
66 324 2, 7, 9, 13-17, 20-24 1240-1375 10 -15,55.-4
67 327 7 1350 15 -8,60
68 357 13, 21, 23, 31 1200-1400 11-20 -13,1.-5,95
69 367 5, 6, 11, 19, 21-24 1170-1260 1 -15,1.-9,5
70 379 1-29 1120-1280 2,8-9,3 -11,57.-9,40
71 419 1-3, 5-10 1210-1240 7,5 -7,73.-7,41
Примечания. * N — номера публикаций в файле bibl.txt базы данных ИНФОРЕКС; ** п — номера опытов из каждой публикации (К) в файлах базы данных ИНФОРЕКС.
Получить последнюю версию ИНФОРЭКС'а можно у авторов базы данных по запросу на электронный адрес агвкт@ geokhi.ru. Кроме того, полную информацию о выборке, включая составы сосуществующих фаз, экспериментальные значения температуры, летучести кислорода и давления могут быть предоставлены после запроса по адресу ekoptev@geol.msu.ru.
ков, Коптев-Дворников, 2005]. Оказалось, что не во всех, но во многих экспериментах уменьшение логарифма летучести кислорода приводит к уменьшению отклонения рассчитанной суммы миналов от 100% и к уменьшению остатков миналов.
Всего в базе данных ИНФОРЕКС нашлось 420 высокобарных экспериментов, выполненных в «сухих» условиях продолжительностью от 24 ч и более, в которых присутствовали расплав и оливин с известным составом. Для 221 из них удалось подобрать такие значения lg/O2, при которых отклонения экспериментальных сумм миналов от 100% не превышали ±3 мол.%. Полученные значения эффективных lg/O2 перекрывают весь диапазон используемых экспериментаторами буферов от СОС до ММО (90% значений lg/O2 лежит пределах от QFM —5 до QFM +5). После добавления высокобарных экспериментов с подобранной летучестью кислорода суммарная выборка включала 789 экспериментальных результатов.
В процессе оптимизации коэффициентов в уравнениях (10) размер выборки сократился до 772 экспериментов (табл. 1). Из выборки были удалены те эксперименты, для которых остатки миналов отклонялись от нормального распределения, а также эксперименты, в которых отклонения рассчитанных сумм миналов от 100% нарушали нормальное распределение этой величины.
Многогранник экспериментальных составов расплавов в координатах концентрации оксидов для окончательной выборки из 772 экспериментальных значений характеризуется следующими величинами (масс.%): 8Ю2 от 39 до 63; ТЮ2 от 0 до 6,5; А1203 от 3 до 21; FeO* от 2 до 32 ^О* — все железо, пересчитанное на Fe0); Mg0 от 1 до 19; Ca0 от 4 до 23; Ка^ от 0 до 7; К^ от 0 до 6. Значения активности компонентов в силикатном расплаве рассчитываются согласно двухрешеточ-ной модели силикатной жидкости Р.Л. Нильсена с соавторами (см. [Френкель и др., 1988]). Эта модель предполагает, что сумма содержания XN а0 + Х'ко меньше содержания Х^ 5. По этой причиине составы расплавов в выборке проверяли на положительное значение Х^ - Х^,а0 - Х^0 , где X1 — мольная доля г'-го компонента в р асплаве. Фигуративные точки составов расплавов выборки представлены на рис. 3. Характер их распределения исключает вероятность наличия ложных корреляций.
Диапазон интенсивных параметров выборки характеризуется температурой от 1040 до 1500 °С, давлением от 1 бара до 30 кбар, летучестью кислорода lg/o2 от —15,8 до —2,6 (от QFM —8,4 до QFM+5,2).
Полученные в результате оптимизации значения коэффициентов приведены в табл. 2. В ряде
случаев в результате статистической обработки линейные тренды на графиках корреляции расчетных и экспериментальных значений существенно отклоняются от линии равных значений, причем отсутствуют эксперименты, сильно отклоняющиеся от общего массива точек, при этом распределение остатков миналов носит нормальный характер. В этих случаях для улучшения согласования между расчетными и экспериментальными значениями вводится дополнительная поправка в виде линейного уравнения
Xm — aX m +
(11)
где X'm — содержание минала, рассчитанное по уравнению вида (8), а и b — коэффициенты в уравнении линейного тренда зависимости между экспериментальными и расчетными содержаниями минала.
Таблица 2
Значения коэффициентов и константы для уравнений (10), найденные c применением надстройки «поиск решения» в программе Excel для миналов оливина, а также поправочных коэффициентов a и b в уравнении (11)
Ко-
эф-фи-циент
A
D
J
JFe3+
JMg
JNa
MgSio 5O2
2472,42
1,17835
-2,42839
0
0,020173
0,209025
2,47999
4,57343
2,63147
0,80597
2,92489
0
3,67513
1,56403
FeSi0 5O2
2359,15
21,3298
-1,91080
0
-0,020714
0,161143
1,09697
0,08051
-0,26096
6,89001
0
-2,47828
1,06578
-0,548423
MnSi0 5O2
4623,86
6,45791
-2,26951
0
0,035951
0,400065
0
-0,66754
-2,44440
-9,67833
2,37490
-2,46664
1,88472
0
5,24420 0,753846 3,53496 25,9473
CaSi0 5O2
5842,07
22,4830
-22,6142
0
0,026191
0,241826
14,1649
13,1513
11,6460
13,4825
19,8625
15,3031
21,5803
18,8068
CrSi0,5O2
0
-27,6742
-34,7209
0
-0,027630
0,008427
32,1139
29,7251
37,2774
44,7893
36,5313
35,2576
34,4670
36,8602
39,4796
Поправочные коэффициенты
b
0,97935
0,014413
0,98080
0,0063191
0,90146
0,0003997
0,98016
0,0001420
0,97443
0,0001048
Тестирование оливинового термобарометра в программе КриМинал. Корреляция экспериментальных составов оливинов с рассчитанными в программе КриМинал [Бычков, Коптев-Дворни-ков, 2005; ВусЬкоу, Кор1еу^уогшкоу, 2014] по уравнениям (10) с коэффициентами из табл. 2 показаны на рис. 4 (для наглядности и облегчения
Рис. 3. Содержание оксидов в экспериментальных стеклах, равновесных с оливином (772 эксперимента в выборке)
20-
о о
<Й
\ 16-PQ О tt К о
g 12-
<и К К
й *
§ и
4-
35
25-
40 45 50 55 60 65 Si02, масс.%
• • • MgO
40 45 50 55 60 65 Si09, масс.%
40 45 50 55 60 65 Si(X масс.%
R
C
E
т;
Л1
"CVO-I-
Га
V
п
Содержание компонента, масс. %, расчет
Содержание компонента, масс. %, расчет
Содержание компонента, масс. %, расчет
Содержание компонента, масс. %, расчет
Содержание компонента, масс. %, расчет
Содержание компонента, масс. %, расчет
Рис. 4. График соответствия расчетных и экспериментальных значений содержания оксидов в оливинах, равновесных с расплавами (772 эксперимента в выборке), и гистограммы распределения соответствующих остатков. Черные кружки — значения, для которых «восстановлена» летучесть кислорода, белые кружки — остальные эксперименты. Остальные условные обозначения см. на рис. 1
Рис. 5. Корреляция между экспериментальными и рассчитанными в КриМинал значениями температуры равновесия оливин-расплав (а) и гистограмма распределения остатков температуры (б) для выборки из 772 экспериментов. Условные обозначения
см. на рис. 4
сравнения с химическими анализами экспериментальных оливинов расчетные мольные содержания миналов пересчитаны на масс.% оксидов).
Верификация воспроизведения экспериментальных значений температуры ликвидуса выполнена также в программе КриМинал (рис. 5).
Судя по гистограммам на рис. 4 и 5, распределение остатков близко к нормальному. Это подтверждается и критерием согласия. Нормальность распределения остатков позволяет использовать для оценки качества термобарометров доверительные интервалы. Преимущество доверительных интервалов — их нацеленность, в отличие от дисперсии, не на оценку качества единичного измерения, а на определение при заданной вероятности границ отклонения расчетной величины от истинного значения. Полезные свойства доверительных интервалов заключаются в возможности определять их для зависимостей (в случае нормального распределения остатков) и возможности их сужения путем увеличения числа измерений.
Средние значения остатков содержания для 8Ю2 0,06±0,03 масс.%, для MgO -0,05+0,07 масс.%, для FeO -0,01+0,09 масс.%, а для МпО, СаО и Сг2О3 они практически равны нулю при ненулевых размерах доверительных интервалов. Это говорит о том, что из всех отклонений значимо только завышение концентрации 8Ю2 в оливине. Но значение завышения, равное 0,06%, много меньше аналитической погрешности измерения этой величины. Средний остаток температуры равен 0,6+1,5 °С, что тоже подтверждает несмещенность оценок с помощью предлагаемого термобарометра.
Для более надежной оценки достоверности описания состава и температуры кристаллиза-
ции оливинов с помощью предложенного нами термобарометра мы построили графики корреляции между расчетными и экспериментальными составами оливинов и значениями температуры равновесия (рис. 4 и 5). В случае отсутствия закономерных отклонений точек на этих графиках от линии равных значений можно утверждать, что предлагаемый набор уравнений адекватно описывает равновесие между расплавом и оливином во всем диапазоне состава оливина и температуры равновесия. Для оценки значимости отклонения от линии равных значений мы использовали доверительный интервал линейной регрессии и доверительные интервалы коэффициентов регрессионного уравнения.
Выяснилось, что для всех компонентов оливина, кроме 8Ю2, нет значимого отклонения линии регрессии от линии равных значений, т.е. линии равных значений лежат внутри узких доверительных коридоров линейных регрессий. Для 8Ю2 угловой коэффициент составляет 0,989+0,009, а свободный член — 0,357+0,333, что приводит к занижению содержания 8Ю2 на 0,03+0,14 масс.% при содержании 8Ю2 в оливине 30% и к завышению на 0,10+0,07 масс.% при 42 масс.%. Таким образом, максимальные погрешности этого термобарометра с учетом величин доверительных интервалов и отклонения линий регрессии от линий равных значений составляют (масс.%) < +0,17 для 8Ю2, +0,38 для FeO, +0,31 для MgO, +0,01 для MnO, +0,06 для CaO, +0,02 для Сг^. Очевидно, что эти погрешности не превышают аналитические.
Для температуры угловой коэффициент, равный 0,981+0,008, а свободный член — 23,3+9,5,
44
с
о М
л
0х О
о св
40
св н ж 1) я
0
8 36^ §
и
1
св
К
<а
и
• = 0,936х+2,564
-2,4 -1,6 -0.8 0 0,8 1,6 2,4
¥еО
^ = 0,953дс+1,38
Я2 = 0,982
32 36 40 44
Содержание компонента, масс. %, расчет
0 20 40 60
Содержание компонента, масс. %, расчет
60 -1
с
о И
л
ч® о4
о о я
Рис. 6. График соответствия рассчитанных по термобарометру из КОМАГМАТ и экспериментальных значений содержания оксидов в оливинах, равновесных с расплавами, и гистограммы распределения соответствующих остатков. Условные обозначения см. на рис. 1
„ 40"
св Н К <и Я о с
2 §
ё20Н
Я
са
Ы о.
<а §
и
у = 0,957х+1,25 Я2 = 0,982
-6 -5-4-3-2-101234
20 40
Содержание компонента, масс.
60
/о, расчет
приводит к занижению температуры на 3,6±2,8 °С при температуре равновесия 1050 °С и к завышению на 5,2±4,9 °С при 1500 °С, а в диапазоне от 1150 до 1450 °С лежит в пределах доверительного интервала. Таким образом, максимальные отклонения расчетной температуры ликвидуса лежат в диапазоне от —7,4 до 10,1 °С, а в большей части температурного интервала не превышают ±3 °С.
Подчеркнем, что фигуративные точки 221 высокобарного эксперимента с «восстановленными» значениями log/O2 естественным образом попадают в общие корреляции для всех компонентов оливина и температуры.
К настоящему времени различными исследователями предложено много способов расчета
равновесия оливин—расплав, по крайней мере 16 из них предложены к использованию в программе Ре1хо1о§ 3 РапушЬеУ8ку, Р1есЬоу, 2011]. Сравнение со всеми вариантами в рамках одной публикации потребует слишком много места и не отвечает нашим целям. Как отмечено во введении, разрабатываемые нами ликвидусные термобарометры предназначены для программы КриМинал, являющейся термодинамическим блоком будущей программы для моделирования динамики формирования мафит-ультрамафитовых интрузивов. Единственное на сегодня программное средство, моделирующее динамику процессов тепломассо-переноса, пусть и в упрощенном виде, — программные комплексы семейства МАГМОД—КО-
Т, °С, КриМинал Т, °С, эксперимент - Т, °С, КриМинал
Рис. 7. Корреляция между экспериментальными и рассчитанными по термобарометру из КОМАГМАТ значениями температуры равновесия оливин—расплав (а) и гистограмма распределения остатков температур (б). Условные обозначения см. на рис. 1
МАГМАТ 3.73. В связи с этим здесь мы сравнили результаты расчета равновесия оливин—расплав с использованием полученного нами оливинового ликвидусного термобарометра с возможностями подхода, реализованного в термодинамическом блоке программы КОМАГМАТ 3.73 [Шрв:// comagmat.web.ru/apps-comagmat.html]. На рис. 6 и 7 приведены результаты воспроизведения состава оливинов и значения их температуры ликвидуса из нашей выборки.
Поскольку в КОМАГМАТ'е версии 3.73 состав оливина рассчитывается на 2 минала — форстерит и фаялит, то расчет воспроизводит содержание только трех оксидов.
В целом КОМАГМАТ неплохо воспроизводит экспериментальные составы и значения температуры. Однако судя по гистограммам на рис. 6 и 7, распределение остатков отклоняется от нормального. Это подтверждается и критерием согласия. Следовательно, доверительные интервалы для оценки качества использованного в КОМАГМАТЕ 3.73 термобарометра использовать не корректно.
Среднее значение остатков содержаний для 8Ю2 равно 0,07 масс.%, для MgO — 0,48 масс.%, для FeO — 0,42 масс.%. Эти значения невелики и сопоставимы с аналитической погрешностью измерения их концентрации. Средний остаток температуры составляет 12,9 °С, что говорит об относительно небольшом, но систематическом занижении расчетной температуры равновесия.
При анализе регрессионных уравнений на графиках корреляции между расчетными и экспериментальными составами оливинов и значениями температуры равновесия (рис. 6 и 7) выяснилось,
что для всех компонентов коэффициенты линии регрессии значимы, но остатки концентраций во всем диапазоне составов оливинов не превосходят типичных аналитических погрешностей. Для температуры угловой коэффициент составляет 1,03+0,012, а свободный член — 23,1+14,2, что приводит к занижению расчетной температуры от 9 °С при температуре равновесия 1050 °С до 22 °С при 1500 °С.
Эти обстоятельства следует учитывать при использовании КОМАГМАТ'а для моделирования динамики формирования расслоенных интрузивов.
Заключение. Для 221 высокобарного эксперимента «восстановлены» значения летучести кислорода, что существенно увеличило выборку экспериментов, пригодных для вывода термобарометров.
Разработан оливиновый ликвидусный термобарометр, который воспроизводит экспериментальные составы без систематических сдвигов. Расчетные оценки состава отличаются от истинных значений на величины, не превышающие аналитическую погрешность. Отклонения рассчитанных значений температуры ликвидуса от истинных не превышают —3,6+2,8 °С при температуре 1050 °С и +5,2+4,9 °С при 1500 °С, а в диапазоне от 1150 до 1450 °С не превышают +3 °С.
Благодарности. Авторы благодарны С.П. Крашенинникову (ГЕОХИ РАН) за полезные советы и замечания. Выражаем признательность А.А. Арискину (МГУ имени М.В. Ломоносова), Г.С. Барминой и Г.С. Николаеву (ГЕОХИ РАН) за создание базы данных ИНФОРЕКС, без которой работа была бы сильно осложнена.
СПИСОК ЛИТЕРАТУРЫ
Арискин А.А., Бармина Г.С. Термометрия равновесий плагиоклазов с расплавами базальтов и андезитов // Геохимия. 1990. № 3. С. 441-447.
Арискин А.А., Бармина Г.С. Моделирование фазовых равновесий при кристаллизации базальтовых магм. М.: МАИК «Наука/Интерпериодика», 2000.
Арискин А.А., Мешалкин С.С., Альмеев Р.Р. и др. Информационно-поисковая система ИНФОРЕКС: анализ и обработка экспериментальных данных по фазовым равновесиям изверженных пород // Петрология. 1997. Т. 5, № 1. С. 32-41.
Арискин А.А., Цехоня Т.И., Френкель М.Я. ЭВМ-барометрия и генетическая интерпретация базальтовых стекол Центральной Атлантики // Геохимия. 1991. № 7. С. 1038-1047.
Арьяева Н.С., Коптев-Дворников Е.В., Бычков Д.А. Ликвидусный термобарометр для моделирования равновесия хромшпинелиды-расплав: метод вывода и верификация // Вестн. Моск. ун-та. Сер. 4. Геология. 2016. № 4. С. 30-39.
Арьяева Н.С., Коптев-Дворников Е.В., Бычков Д.А. Ликвидусный термобарометр для моделирования равновесия магнетит-расплав // Вест. Моск. ун-та. Сер. 4. Геология. 2018. № 1. С. 70-79.
Бычков Д.А., Коптев-Дворников Е.В. Программа КриМинал для моделирования равновесия расплав— твердые фазы при заданном валовом составе системы // Мат-лы междунар. конф. Улан-Удэ: Изд-во БурНЦ СО РАН, 2005. С. 122-123.
Бычкова Я.В., Коптев-Дворников Е.В. Ритмическая расслоенность киваккского типа: геология, петрография, гипотеза формирования // Петрология. 2004. Т. 12, № 3. С. 281-302.
Коптев-Дворников Е.В., Арьяева Н.С., Бычков Д.А. Уравнение термобарометра для описания сульфид-силикатной ликвации в базитовых системах // Петрология. 2012. Т. 20, № 5. С. 1-18.
Коптев-Дворников Е.В., Бычков Д.А. Геотермометры для широкого диапазона составов базитов // Мат-лы междунар. конф. «Ультрабазит-базитовые комплексы складчатых областей». Иркутск: Изд-во СО РАН, 2007. С. 178-181.
Коптев-Дворников Е.В., Ярошевский А.А., Френкель М.Я. Кристаллизационная дифференциация интрузивного магматического расплава. Оценка реальности седиментационной модели // Геохимия. 1979. № 4. С. 488-508.
Френкель М.Я. Тепловая и химическая динамика дифференциации базитовых магм. М.: Наука, 1995. 239 с.
Френкель М.Я., Ярошевский А.А. Кристаллизационная дифференциация интрузивного магматического расплава // Геохимия. 1978. № 5. С. 643-668.
Френкель М.Я., Ярошевский А.А., Арискин А.А. и др. Динамика внутрикамерной дифференциации базитовых магм. М.: Наука, 1988. 215 с.
Черных Н.С. Влияние физико-химических параметров на отделение рудных фаз от базитовых магм (по результатам математического моделирования): Автореф. канд. дисс. М., 2017.
Asimow P.D., Ghiorso M.S. Algorithmic modifications extending MELTS to calculate subsolidus phase relations // Amer. Mineralogist. 1998. Vol. 83, N 9-10. P. 1127-1132.
Bychkov D.A., Koptev-Dvornikov E.V. The Software for Simulation of Equilibrium Crystallization // Gold-schmidt2014 Abstr. URL: https://whiteiron.org/uploads/ conferences/24/abstracts/A-Z.pdf: 2014. P. 319-319.
Danyushevsky L.V., Plechov P. Petrolog3: Integrated software for modeling crystallization processes: PETRO-LOG3 // Geochemistry, Geophysics, Geosystems. 2011. Vol. 12, N 7. P. 1-32.
Nielsen R.L., Dungan M.A. Low pressure mineral-melt equilibria in natural anhydrous mafic systems // Contrib. Mineral. Petrol. 1983. Vol. 84, N 4. P. 310-326.
Поступила в редакцию 18.02.2019
Поступила с доработки 26.04.2019
Принята к публикации 26.04.2019