Физика Земли, 2022, № 1, стр. 130-139

ВВЕДЕНИЕ

Мусковит является геофизическим материалом, участвующим во многих петрологических процессах, происходящих в глубинных слоях Земли. В связи с этим экспериментальному исследованию его теплофизических свойств и уравнений состояния при высоких динамических и статических давлениях посвящен ряд работ (см. [Sekine, 1991; Faust, Knittle, 1994; Curetti, 2006], их цитирование и ссылки в них). Так, в работе [Sekine, 1991] определена ударная адиабата мусковита в диапазоне давлений 20–140 ГПа. При этом авторы [Sekine, 1991] сделали вывод, что слюда устойчива вдоль ударной адиабаты до давлений 80 ГПа. В работе [Faust, Knittle, 1994] выполнены исследования сжатия кристаллического мусковита вдоль его комнатной изотермы до 20 ГПа. Авторы [Faust, Knittle, 1994] установили, что при комнатной температуре слюда теряет дальний кристаллический порядок при P_tr = 18 ГПа и аморфизируется вблизи 27 ГПа. Изменение кристаллографических характеристик мусковита при давлении 18 ГПа также установлено в работе [Curetti, 2006]. Авторы работы [Faust, Knittle, 1994] предположили также, что ударная адиабата в работе [Sekine, 1991] относится к аморфному состоянию мусковита при давлениях, превышающих ≈20 ГПа.

Таким образом, согласно процитированным работам в диапазоне давлений 0–80 ГПа область 0–18 ГПа принадлежит кристаллическому мусковиту, а в области 18–80 ГПа располагается смесь фаз и аморфный мусковит.

Построение уравнений состояния (в смысле отыскание функциональной связи между давлением, объемом и температурой) фаз мусковита носят преимущественно теоретический характер. Так, в работе [Соколова, 2014] предложен полуэмпирический подход к построению уравнений состояния кристаллического мусковита, учитывающий сильную анизотропию этого материала. Первопринципные расчеты уравнений состояния кристаллического мусковита представлены в работах [Teich-McGoldrick, 2012; Ulian, Valdrè, 2015]. Что же касается уравнений состояния аморфного мусковита, то в доступной литературе они не обнаружены. В связи с этим цель данной статьи заключалась в конструировании полуэмпирических уравнений состояния аморфного мусковита, согласованных с уравнениями состояния его кристаллической фазы.

УРАВНЕНИЯ СОСТОЯНИЯ КРИСТАЛЛИЧЕСКОГО МУСКОВИТА

Индивидуальные уравнения состояния мусковита строились на основе полуэмпирического выражения для свободной энергии всесторонне сжатого твердого тела F = F(V, T) в виде функции удельного объема V и температуры T. Для кристаллического мусковита использовалось выражение:

где: R – удельная газовая постоянная; Θ = Θ(V) – характеристическая температура Эйнштейна; E_x = = E_x(V) – потенциальная энергия.

Вид функций E_x = E_x(V) и Θ = Θ(V) определяется формулами из работы [Молодец, 2017а]. Так, характеристическая температура Θ определяется как

В (2) параметр υ₀, имеющий смысл характеристического объема, определяется соотношением:

где γ₀ = γ₀(V₀, T₀) термодинамический параметр Грюнайзена при начальном удельном объеме V₀, начальной комнатной температуре T₀ = 298.15 К и начальном атмосферном давлении P₀. Соответственно, параметр Θ₀ представляет собой начальное значение характеристической температуры Θ₀= Θ(V₀).

Выражение для E_x= E_x(V) имеет вид:

Выражение (1) содержат два типа опорных параметров. Во-первых, это начальные справочные теплофизические свойства материала в (2): γ₀, Θ₀, V₀. В табл. 1 представлены используемые значения теплофизических свойств кристаллического мусковита KAl₂(AlSi₃)O₁₀(OH)₂ вместе с его формульной молярной массой M.

Вторую группу составляют параметры для потенциальной составляющей свободной энергии – это подгоночный параметр υ_x, первое приближение которого есть υ₀ и константы C₁, C₂, C₃, выражающиеся через υ_x и справочные теплофизические свойства материала (см. [Молодец, 2017а]).

Согласно термодинамическим тождествам, частные производные свободной энергии позволяют построить калорическое (зависимость внутренней энергии E = E(V, T) от объема и температуры) и термическое (зависимости давления P = = P(V, T) от объема и температуры) уравнения состояния. Рассмотрим в этом разделе построение термического уравнения состояния P = –(∂F/∂V) кристаллического мусковита с использованием (1) и параметров из табл. 1.

Дифференцирование (1) по объему дает

где P_x = –dE_x/dV есть потенциальное давление, которое в соответствии с (4) имеет вид:

В (7) Θ = Θ(V) – характеристическая температура (2), γ = γ(V) – объемная зависимость коэффициента Грюнайзена:

Таким образом, в полуэмпирическом выражении термического уравнения состояния (7) для кристаллического мусковита оказываются определенными все опорные константы: параметры γ₀, Θ₀, V₀, K₀ заданы в табл. 1, расчет υ₀ осуществляется по формуле (3), расчет C₁, C₂ осуществляется согласно работе [Молодец, 2017а], для υ_x принимается первое приближение υ_x = υ_0.

Заключительный этап построения термического уравнения состояния (7) заключается в уточнении подгоночного параметра υ_x. Для кристаллического мусковита подгонка осуществлялась следующим образом. На экспериментальной комнатной изотерме [Curetti, 2006] выбиралась средняя точка с координатами {P = 10.54 ГПа; V = 0.317 см³/г} (см. звезду 1 на рис. 1). Затем “вручную” задавались различные значения υ_x, при каждом из которых рассчитывалась комнатная изотерма (7). Из полученного семейства расчетных изотерм выбиралась та, которая проходила через выбранную среднюю точку. В рассматриваемом случае такой изотермой оказалась кривая 2, для которой значение υ_x составило υ_x = 228.112 см³/моль. Найденная величина υ_x и соответствующие константы C₁, C₂ вместе с υ₀ представлены в табл. 2. Константа C₃, необходимая для полного определения функции свободной энергии (1) находится (см. [Молодец, 2017а]) из условия, согласно которому минимальное значение потенциальной энергии E_x(4) должно быть равно нулю. Значение C₃ также приведено в табл. 2.

На рис. 1 показана оптимальная расчетная комнатная изотерма кристаллического мусковита на фоне экспериментальных точек [Curetti, 2006; Faust, Knittle, 1994]. Как видно, расчетная изотерма удовлетворительно согласуется с экспериментами [Curetti, 2006; Faust, Knittle, 1994] до 18 ГПа и заметно отличается от экспериментальных данных в работе [Curetti, 2006] при более высоких давлениях. На этом основании можно обозначить верхнее давление области применимости разработанного уравнения состояния кристаллического мусковита величиной ≅18 ГПа, где расчет согласуется с экспериментом обеих работ.

Вместе с этим отметим, что в области давлений ≅5 ГПа, где наряду с экспериментальными результатами [Curetti, 2006] имеются данные и других авторов – это экспериментальная [Comodi, 2002] и теоретическая [Ulian, Valdrè, 2015] комнатные изотермы. В этой области давлений расчетная изотерма 2 практически совпадает с данными из работ [Curetti, 2006; Ulian, Valdrè, 2015], но на 30% отличается от данных в работе [Comodi, 2002]. Согласно [Ulian, Valdrè, 2015], различие экспериментальных изотерм разных авторов частично связано с химическим составом образцов, а близость экспериментальной и теоретической изотермы обязано степенью кристалличности образов. В этой связи разумно полагать, что во всем диапазоне давлений 0–18 ГПа расчетная изотерма 2 на рис. 1 оптимальна для кристаллического мусковита из работ [Curetti, 2006; Ulian, Valdrè, 2015].

ОБСУЖДЕНИЕ ПРОГНОСТИЧЕСКИХ ВОЗМОЖНОСТЕЙ УРАВНЕНИЯ СОСТОЯНИЯ КРИСТАЛЛИЧЕСКОГО МУСКОВИТА

Отметим, что уравнение состояния (7) в комплекте с (2) и (9) позволяет провести численные расчеты объемно-температурных зависимостей нужного свойства материала, например, изохорной теплоемкости C_V = C_V(V, T):

изотермического модуля объемного сжатия B = = B(V, T) коэффициента объемного теплового расширения β = β(V, T) изобарной теплоемкости C_P = C_P(V, T)

Пользуясь этой возможностью и параметрами из табл. 1 и табл. 2, обсудим применимость разработанного полуэмпирического описания B, β, и C_P кристаллического мусковита при высоких давлениях и температурах путем сравнения расчетных и соответствующих литературных данных.

Заметим предварительно, что (7) представляет собой явную функцию объема и температуры, позволяющая вычислить давление и его частные производные по объему вдоль изотермы (P = P(V, T_i), где T_i = const). Однако аналитическое вычисление, например, объема вдоль изобары (V = V(P_i, T) где P_i = const) с помощью (7) затруднительно, так как преобразовать (7) в аналитическую формулу V = V(P, T) не удается. Тем не менее в таких случаях задача нахождения изобары V = V(P_i, T) легко решается для (7) численно при каждом значении T.

На рис. 2 представлен расчет 1 и 2 изотермического модуля объемного сжатия мусковита по формуле (11) вдоль двух изотерм (7). На этом же рисунке приведены соответствующие результаты первопринципных расчетов 3 и 4 из работы [Teich-McGoldrick, 2012]. Как видно, 1 и 2 согласуются с 3 и 4 качественно – оба результата демонстрируют нелинейное возрастание модуля объемного сжатия кристаллического мусковита с увеличением давления и его уменьшение c возрастанием температуры. В количественном отношении оба результата совпадают между собой с точностью 10–15%.

На рис. 3 цифрой 1 обозначен график (12) на фоне экспериментальных 2 и теоретических 3 литературных данных для коэффициента объемного теплового расширения мусковита. Как видно, 1 согласуется с 2 и 3 в качественном и в количественном отношении – наибольшее различие не превышает 15%.

Для расчета температурной зависимости изобарной теплоемкости по (13) привлекалась расчетная изобара V = V(P_i, T) при P_i = 1.0 атмосфер. График формулы (13) представлен кривой 1 на рис. 4. Здесь же приведены табличные данные из работ [Robie, Hemingway, 1995; Robie, 1976] и теоретические из работы [Ulian, Valdrè, 2015] результаты температурной зависимости изобарной теплоемкости кристаллического мусковита 2, 3, 4. Как видно, график 1 дает заниженные значения теплоемкости в области низких температур, начиная с ≈300 K. Это занижение обусловлено неточностью приближения Эйнштейна в области низких (T < ≈0.4Θ) температур. В диапазоне температура 300–1000 К расчетные значения 1 отличаются от результатов из работ [Robie, Hemingway, 1995; Robie, 1976; Ulian, Valdrè, 2015] не более, чем на 5%.

Таким образом, построенное уравнение состояния кристаллического мусковита (7) отображает объемно-температурные зависимости теплофизических свойств этого материала с точностью 5–15% в при давлениях до 18 ГПа и в температурном диапазоне 300–1000 К. Рассмотрим теперь построение полуэмпирического уравнения состояния аморфного мусковита, считая основным состоянием мусковита его кристаллическую фазу.

УРАВНЕНИЯ СОСТОЯНИЯ АМОРФНОГО МУСКОВИТА

Как хорошо известно, при построении полуэмпирических уравнений состояния аморфного твердого тела определенным препятствием является трудность получения его экспериментальной изотермы высокого давления, которая для кристаллических материалов в настоящее время сравнительно легко определяется в алмазных наковальнях, оснащенных сильным рентгеновским источником. Однако в измерении ударных адиабат аморфного материала трудностей не возникает. В этой связи отметим, что формализм построения полуэмпирических уравнений состояния твердых тел из предыдущего раздела может быть применен к аморфной фазе с опорой на экспериментальную ударную адиабату. Рассмотрим эту возможность на примере аморфного мусковита и построим его уравнения состояния, используя приемы [Молодец, 2017б] и данные ударно-волнового эксперимента [Sekine, 1991].

Будем считать, что выражение для свободной энергии аморфного материала можно записать в виде:

где первое слагаемое F(V, T) задается прежними соотношениями (1)–(6) с тем лишь различием, что набор справочных свойств, подгоночный параметр υ_x и константы C₁, C₂, C₃ принадлежат аморфному материалу. Второе слагаемое E_m представляет разность потенциальных энергий кристаллического и аморфного материала. Третье слагаемое представляет аналог конфигурационной части свободной энергии жидкости, в котором R – газовая постоянная, Т – температура. Таким образом, для аморфного состояния, подобно данным в работе [Молодец, 2017б], помимо υ_x подгоночными считаются также и коэффициенты a_S и E_m.

Как и выше, термическое уравнение состояния P = P(V, T) определяется частной производной свободной энергии по объему. Для аморфного мусковита использование (14) дает P = –(∂F_a/∂V) = = –(∂F/∂V), а его функциональный вид не отличается от (7). В свою очередь частная производная свободной энергии по температуре позволяет записать калорическое уравнение состояния E = = E(V, T) как E = F_a − T(∂F_a/∂T) = F − T(∂F/∂T) + E_m, что с использованием (1) и (10) дает:

Заметим далее, что при работе с данными ударно-волнового эксперимента предпочтительно уравнение состояния в форме Ми-Грюнайзена.

которое получается путем комбинации соотношений (7) и (15).

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

и подставить энергию E = E_h из (17) в (16), то получим давление P = P_h(V) вдоль ударной адиабаты в виде функции только объема:

Соотношение (18) представляет собой стандартную формулу физики ударного сжатия конденсированных сред с индивидуальными свойствами материала, характеризуемого значением V₀ и функциями E_x(V) и γ(V). Например, для кристаллического мусковита параметры этих функций приведены в табл. 1 и табл. 2, а на рис. 5 представлены их графики.

Как уже отмечалось, в данной работе принято, что функции E_x(V) и γ(V) имеют одинаковый вид для обеих фаз мусковита, но отличаются значениями справочных свойств, подгоночного параметра υ_x, а также имеют дополнительные подгоночные параметры a_S и E_m. Значения подгоночных коэффициентов в (14) для аморфного мусковита находились в следующем порядке. Во-первых, предполагалось, что энтропия разупорядоченного мусковита S_0d = S_d(V₀, T₀) из работы [Robie, Hemingway, 1995] равна энтропии аморфного мусковита. Затем величина параметра a_S подбиралась так, чтобы энтропия S₀ = S(V₀, T₀) аморфного мусковита при нормальных начальных условиях равнялась величине энтропии разупорядоченного мусковита S_0d = 306.4 Дж/(мольК) из работы [Robie, Hemingway, 1995]. В качестве первого приближения для E_m была взята величина разности энтальпий 15.2 кДж/моль разупорядоченного и упорядоченного мусковита из работы [Robie, Hemingway, 1995].

Величина подгоночного параметра υ_x аморфного мусковита подбиралась с опорой на экспериментальную ударную адиабату [Sekine, 1991]. Для этого вначале определялась средняя точка [Sekine, 1991] с координатами {P = 75.75 ГПа; V = = 0.2233 см³/г} (см. звезду 5 на рис. 1). Затем “вручную” задавались различные значения υ_x, при каждом из которых рассчитывалась ударная адиабата (18). Из полученного семейства расчетных ударных адиабат выбиралась та, которая проходила через выбранную среднюю точку. В рассматриваемом случае такой ударной адиабатой оказалась кривая 6, для которой значение υ_x составило υ_x = 334.229 см³/моль. Величины C₁, C₂, C₃ для аморфного мусковита находились с использованием прежних формул из работы [Молодец, 2017а]. Наконец, величина E_m уточнялась так, чтобы химические потенциалы ${{\mu }} = F\left( {V,T} \right) + PV$ кристаллического ${{{{\mu }}}_{с}} = {{{{\mu }}}_{с}}\left( {{{P}_{{tr}}},{{T}_{0}}} \right)$ и аморфного ${{{{\mu }}}_{a}} = {{{{\mu }}}_{a}}\left( {{{P}_{{tr}}},{{T}_{0}}} \right)$ мусковита равнялись между собой при давлении P_tr = 18 ГПа и температуре T₀ = = 298.15 K. Полученное значение E_m = 0.3105 кДж/г вместе с a_S, υ_x и набором констант C₁, C₂, C₃ для аморфного мусковита также представлены в табл. 2. На рис. 5 представлены графики ключевых функций E_x(V) и γ(V) уравнения состояния Ми-Грюнайзена для аморфного мусковита.

Таким образом, соотношения (1) и (14) вместе с коэффициентами табл. 1 и табл. 2 представляют собой взаимосогласованное описание термодинамических свойств не только кристаллической, но и аморфной фазы мусковита в области высоких статических и динамических давлений.

ИНТЕРПРЕТАЦИЯ ТЕПЛОФИЗИЧЕСКИХ СВОЙСТВ АМОРФНОГО МУСКОВИТА НА ОСНОВЕ ПОCТРОЕННЫХ УРАВНЕНИЙ СОСТОЯНИЯ

Интерпретация ряда свойств кристаллического мусковита с позиций разработанного полуэмпирического описания была дана в предыдущем разделе. Обсудим в этой логике описание некоторых теплофизических свойств аморфного мусковита при высоких давлениях и температурах.

Расчетная комнатная изотерма (7) аморфного мусковита представлена на рис. 1 графиком 7. Как видно, комнатная изотерма аморфного мусковита проходит через максимальную по давлению точку из работы [Faust, Knittle, 1994]. При более высоких давлениях данные по изотермическому сжатию аморфного мусковита в доступной литературе не обнаружены. Однако в области высоких давлений, как отмечалось выше, имеются данные ударно-волновых экспериментов. Рассмотрим результаты расчетов для аморфного мусковита, которые можно сопоставить с экспериментом в работе [Sekine, 1991] в координатах скорость ударной волны(D)–массовая скорость(u).

Как хорошо известно, решение (18) и законы сохранения импульса и массы на ударном скачке позволяют обычным способом рассчитать скорость ударной волны D:

и массовую скорость u:

Отметим, что расчет ударных адиабат как аморфного, так и кристаллического мусковита осуществляется по общим формулам (18)–(20) с тем лишь уточнением, что в качестве начальных значений объема V₀ и энергии E₀ в этих формулах в обоих случаях берутся величины для кристаллического мусковита V₀ = 0.3535 см³/г и E₀ = 0.58173 кДж/г. Кроме этого, для аморфного мусковита расчет начинается с давлений, превышающих область двухволнового ударного сжатия ~0–30 ГПа.

Рассчитанная с помощью (19), (20) ударная адиабата аморфного мусковита (с соответствующими параметрами из табл. 1 и табл. 2) показана графиком 1 на рис. 6. Здесь же представлен первичный эксперимент 2 и его линейная аппроксимация 3 из работы [Sekine, 1991]. Как видно, расчетная ударная адиабата аморфного мусковита за исключением точки {2.44 км/с; 8.01 км/с} совпадает с экспериментом 2 в пределах его погрешности и практически не отличается от эмпирической аппроксимации 3 из работы [Sekine, 1991]. Таким образом, точность разработанного термодинамического описания применительно к ударному сжатию аморфного мусковита находится на уровне точности 3% имеющихся экспериментальных данных.

На рис. 6 показана также аналогично рассчитанная по тем же формулам (19), (20) ударная адиабата 4 кристаллического мусковита. Отметим, что хотя экспериментальных ударно-волновых данных по сжатию кристаллического мусковита в этой области не обнаружено, тем не менее можно предполагать, что точность расчета ударной адиабаты кристаллического мусковита также составляет ~3%.

Заметим далее, что соотношения (7) и (18) позволяют рассчитать температуру ударного сжатия. Действительно, если в (7) положить, что объем и температура принадлежат ударной адиабате, то величина давления P в (7) будет равна значению P_h (18), то есть P = P_h. Решая это уравнение относительно температуры T_h найдем T_h = T_h(V, P_h) в виде функции объема и давления вдоль ударной адиабаты:

где

В (21), (22) функции Θ(V), P_x(V), γ(V), P_h(V) задаются формулами (2), (8), (9), (18) соответственно.

На рис. 7 представлен расчет температуры ударного сжатия (21) кристаллического (график 1) и аморфного (график 2) мусковита в зависимости от давления ударного сжатия. На этом же рисунке штриховой линией представлен участок линии разложения мусковита в области высоких давлений и температур из работы [Sekine, 1991]. Как видно, линия 2 пересекает линию 3 в точке {P = 72 ГПа; T = 1770 K}. На этом основании можно обозначить область фазовой диаграммы мусковита, где применимо предлагаемое описание аморфного мусковита как {25 ГПа < P < 70 ГПа; 300 К < T < 1770 K}.

В заключение оценим местоположение границы начала аморфизации кристаллического мусковита в координатах давление–температура по аналогии с нахождением линии равновесия между фазами. Для этого заметим, что формулы (1) и (14) в комплекте с параметрами из табл. 1 и табл. 2 позволяют вычислять химические потенциалы кристаллического μ_c = F(V_c, T) + PV_c и аморфного μ_a = F_a(V_a, T) + PV_a мусковита. При этом термическое уравнение состояния (7) для каждой фазы позволяет исключить ее объем (кристаллической фазы V_c и аморфной фазы V_a) и представить химические потенциалы кристаллической μ_c и аморфной μ_a фаз в виде функции своих переменных – давления и температуры. Поэтому граница T_b = = T_b(P) между кристаллической и аморфной фазой может быть рассчитана как линия пересечения поверхностей химических потенциалов кристаллической и аморфной фаз в пространстве {μ, P, T}, задаваемая уравнением:

Определенная таким способом граница начала аморфизации кристаллической фазы мусковита T_b = T_b(P) показана на рис. 7 линией 4.

Отметим, что в выполненных расчетах линия 4 имеет нижнюю границу по температуре, обусловленную используемой моделью Эйнштейна, которая при низких (по сравнению с характеристической температурой Θ) температурах отличается от более точной теории Дебая. Но, как известно, теории Эйнштейна и Дебая дают существенно различные результаты лишь при температурах ниже T_min ≈ 0.4Θ и практически совпадающие результаты при высоких температурах T_min > ≈0.4Θ, что для мусковита составляет T_min > 300 K. Таким образом, в диапазоне температур T_min < T < 450 К граница T_b между кристаллической и аморфной фазами мусковита может быть представлена прямой с отрицательным наклоном dT_b/dP ≈ –47.8 К/ГПа.

ЗАКЛЮЧЕНИЕ

Реконструированы полуэмпирические уравнения состояния (термическое, калорическое, Ми-Грюнайзена) кристаллического и аморфного мусковита, предназначенные для расчета теплофизических свойств мусковита в диапазоне давлений 0–70 ГПа и температур 300–1700 К. Найдены количественные значения опорных параметров таких ключевых функций, как объемные зависимости коэффициентов Грюнайзена и потенциальных энергий обеих фаз мусковита. Проведено сравнение между разнообразными литературными данными (изотерма высокого давления, ударная адиабата, модуль объемного сжатия, теплоемкость, коэффициент теплового расширения) и собственными прогностическими расчетами объемно-температурных теплофизических свойств при статическом и ударном сжатии мусковита. В целом результаты выполненных расчетов позволяют заключить, что в рассмотренном диапазоне давлений предложенное термодинамическое описание обеих фаз мусковита соответствует эксперименту в пределах его погрешности. Представлена оценка местоположения границы начала аморфизации кристаллической фазы мусковита в виде прямой линии с отрицательным наклоном dT_b/dP ≈ –47.8 К/ГПа.