close

Вход

Забыли?

вход по аккаунту

?

bd000101575

код для вставкиСкачать
Российская академия н а у к
И н с т и т у т физики Земли и м . О . Ю . Ш м и д т а
На правах рукописи
УДК550.311, 552.11,551.21
С и м а к й н Александр Геннадьевич
МОДЕЛИРОВАНИЕ ФАЗОВЫХ ПЕРЕХОДОВ
И МАССОПЕРЕНОСА
В МАГМАТИЧЕСКИХ КАМЕРАХ
Специальность 25.00.10 геофизика, геофизические методы поисков полезных ископаемых
Автореферат диссертации на соискание ученой степени
доктора физико-математических наук
Москва 2005
Работа вьтолнена в Институте экспериментальной минералогаи Р А Н
Официальные оппо,невты:
Доктор геолого-минералогических наук Шмулович К.И. ( И Э М РАН)
Доктор физико-математических наук Каракин А.В. (ВНИИ Геосистем)
Доктор физико-математических наук Геншафт Ю,С. (ИФЗ РАН)
Ведущая организация: Институт геохимии и аналитической химии
им. В.И.Вернадского РАН
Защита состоится "23" ноября в 11 часов на заседании диссертационного
совета Д.002.001.01 в Институте физики Земли им. О.Ю.Шмидга РАН
(123810 Москва, Б.Грузинская, 10).
С диссертацией можно знакомиться в библиотеке Института физики Земли
РАН.
Автореферат разослан " 7 7 "
октября
2005 г.
Ученый секретарь Диссертационног^р совета
кандидат физ.-мат. наук
А.П.Трубицын
oni^o
101ЪО
бр^'^^с^/^ 7
ВВЕДЕНИЕ
Необходимость и актуальность исследований. Наша живая планета
находится в постоянном движении. Мантийной конвекцией сдвигаются
континенты, продолжаются процессы плавления и перемещения
расплавленной тверди, периодически магма выплескивается наружу,
выступая наглядным доказательством подземного жара. Логика познания
глобальных геодинамических процессов, закономерностей эволюции
магматических камер, будь то субвулканические магматические очаги,
лавовые озера, палеоинтрузии, несущие руду, или линзы частично
расплавленных
пород,
подстилающих
активные
современные
геотермальные системы - ведет к рассмотрению конвекции в гетерофазных
средах с фазовым переходом. Особенное практическое значение в
последние время признается за необходимостью предсказания причины,
даты и места извержения супервулканов. В исторически обозримом
прошлом (75 тьгс. лет назад) извержение такого супервулкана в Индонезии,
которое превосходит на порядки зарегистрированные современные
извержения, послужило причиной многолетней зимы и вымирания более
90% жившего тогда человека разумного. Есть кандидаты на подобные
извержения и в настоящее время. Перед учеными, изучающими
геотермальные системы, также есть возможность внести свой вклад в
решение
энергетической
проблемы
будущего
и
преодолеть
углеводородный голод, толкающий целые нации на путь разбоя и войн.
Тепловая и композиционная конвекция в магматических и связанных с
ними гидротермальных системах, как правило, сопряжена с фазовыми
переходами, с гетерогенизацией растворов и расплавов и процессами
гравитационной сепарации фаз. Эти явления в последнее время стали очень
модным объектом для исследований. Если предыдущие декады можно
было бы назвать эпохой «двойной диффузионной» конвекции, то в
настоящее время большое внимание придается пузырению, кристаллизации
природных магм, процессам гетерогенизации солевого флюида. В нашей
стране застывание интрузий первоначально рассматривалось B.C.
Голубевым и В.Н. Шараповым по аналогии с направленной
кристаллизацией
металлических
слитков
с
диффузионным
перераспределением примесей без учета конвекции и осаждения
кристаллов. На полуколичественном уровне моделирование становления
магматических
интрузий
с
учетом
гравитационной
сегрегации
кристаллических фаз было предпринято группой М.Я.Френкеля в 70-90ые
годы 20 века. Но конвективные процессы при этом рассматривались не на
строго физическом уровне. В настоящее время сопряженные процессы
пузырения и кристаллизащи в ходе течения магмы от очага к поверхности
рассматриваются О. Мельником и А. Барминым. При теоретическом
описании классическая термическая кондяутшя и готярпфячныр течения
РОС НАЦИОНАЛЬНА
6ИБЛИОТЕКА
Cnetcntyf
О» ЩРш>^^у'^
противопоставляются. Последние рассматриваются как одна из форм
течений Рэлея-Тейлора в многослойной системе, состоящей из чистых
жидкостей, расположенных в пространстве в последовательности, ведущей
к гравитационной неустойчивости (легкая жидкость расположена под
тяжелой). Тот факт, что в гетерофазной системе гравитационная
неустойчивость
достигается
распределением
переменной
доли
диспергированной фазы (газообразной, твердой) в объеме жидкости,
игнорируется. А диспергированная фаза может перемещаться в поле силы
тяжести, что придает этому фактору плотности сходство с переносом тепла.
Основная иель настоящего исследования состояла в изучении
динамики процессов фазового перехода в расплавах на локальном уровне, в
погранслое течения, процессов конвекции в гетерогенных средах с
перемещающимися фазами и приложении теоретических результатов к
геотермальному и вулканическому процессам. В соответствии с
поставленной целью решались следующие конкретные задачи:
1. Изучение процесса пузырения магматического расплава, измерение
распределения флюидных пузырей по размерам, определение фрактальной
размерности
распределений
пузырей
в
природных
лавах
и
экспериментальных образцах, анализ физических механизмов нуклеации и
взаимодействия пузырей при их росте.
2. Изучение кристаллизации, сопряженной с декомпрессией и
дегазахщей магматического расплава на примере гавайитовых магм
характерных для вулкана Этна, крупнейшего действующего вулкана
Европы.
Содержательный
анализ
механизмов
формирования
экспоненциального распределения по размерам кристаллов магматических
минералов.
3. Проведение исчерпывающего анализа решений задачи о
пограничном слое застывающего пластообразного магматического тела в
равновесном приближении бинарного расплава с диаграммой плавкости
эвтектического типа.
4. Приложение модели межфазного гравитационного течения в среде с
фазовыми переходами к расчету тепло и массопереноса в гидротермальной
двухфазной флюидной системе соль-вода.
5. Количественный анализ конвекции в гетерогенной жидкости при
заданном потоке дисперсной фазы на границах. Нахождение критических
условий возникновения конвекции в терминах седиментационного числа
Рэлея по аналогии с термической конвекцией.
6. Разработка алгоритма и написание программы для численного
моделирования процессов гетерогенной конвекции в магматической
системе. Проверка работы программы на тестовых примерах.
7. Использование прохраммы для расчета процесса конвекции,
вызванной пузырением и сопряженной кристаллизацией в магматических
камерах под вулканами.
Научная новизна.
1. Впервые проанализированы различные микрорежимы пузырения
магматического расплава. Равномерная гетерогенная нуклеация при
пониженном пороге активационного барьера фазового перехода приводит к
заполнению пространства пузырями разных генераций с приближением
распределения по размеру к фрактальному, с фрактальной размерностью
2.5-2.8. Режим фронтального пузырения возможен из-за существенных
вязких напряжений у растущих пузырей малых, близких к микронным,
размеров.
2. Впервые комплексно проанализирована диаграмма плавкости и
кинетика роста главных кристаллических фаз гавайитов при низких
давлениях водяного флюида. Рассмотрение скоростей роста и особенностей
состава кристаллов плагиоклаза и клинопироксена по микропримесям дало
возможность связать экспоненциальный участок функции распределения
этих минералов в гавайитовых лавах вулкана Этна с процессом дегазации и
сопряженной кристаллизации лавы в жерле вулкана Этна.
3. Впервые получено точное решение задачи Стефана, сопряженной с
объемной кристаллизацией оседающей фазы в двухкомпонентной
двухфазной системе. В зависимости от скорости оседания кристаллов и
продвижения границы полного застывания возможны различные варианты
распределения оседающих кристаллов, включая разрывные близкие по
природе к ударным волнам, о существовании которых было впервые
заявлено М.Я.Френкелем.
4. Впервые математически описано стационарное одномерное течение
в двухфазном флюиде на примере системы NaCl-H20. Оценены потоки
тепла и массы через колонну флюида, ограниченную снизу границей
кипения солевого флюида, а сверху конденсации газообразной фазы.
5. Впервые получены решения задачи о кристаллизации в пограничном
слое с оседающими кристаллами на уровне распределений кристаллов по
размерам.
6. Впервые показано, что помимо кратковременного течения типа
Рэлея-Тейлора в гетерофазных жидкостях с относительным перемещением
фаз возможна перманентная конвекция при постоянном потоке оседающей
(всплывающей) фазы. Найдено минимальное значение седиментационного
числа Рэлея, при которым возможна стационарная конвекция, и которое
оказалось близко к бпределенному ранее В,П.Трубицынам и Е.Харыбиным
(1989) методом анализа устойчивости бескнечно малых возмущений.
7. Разработан метод решения задачи моделирования конвекции при
существенно адвективном переносе диспергированной фазы с внутренним
источником. Впервые в геофизических приложениях маркерная методика
применена для решения на характеристиках задач о перемещении границы
двухфазных областей, ограниченных фронтом растворения.
8. Впервые получены численные решения, описывающие процессы
конвекции, сопряженной с пузырением и кристаллизацией в магматических
очагах под вулканами.
Практическое значение. Экспериментальные данные и модели дегазации
магматических расплавов могут быть практически применены в
вулканологии. Они дают представление о том, насколько неравновесным
может быть процесс дегазации, и как кинетика пузырения может привести
к катастрофическим извержениям, позволяют интерпретировать динамику
современных и палеоизвержений по распределению по размерам пузырей и
кристаллов в лавах. Модель одномерного гравитационного течения в
колонне гетерофазного реагирующего флюида была применена для оценки
скорости тепло- и массопереноса в корневых частях геотермальных систем.
Теория седиментационной конвекции применима ко многим конвективным
системам как технического, так и геологического значения. Программу для
расчета конвекции, вызванной дегазацией и сопряженной кристаллизацией
магмы, можно использовать для расчета флюидного потока на некоторых
вулканах андезитового состава Курильской и Японской островной дуг. Как
показано автором, этот код также можно использовать для численного
моделирования
солевьге
экструзий
в
условиях
перманентного
поверхностного растворения.
Исходные материалы и методы исследований. Экспериментальная часть
работы базируется на более 150 экспериментах с гранитным и
гавайитовым водосодержащим расплавом. Экспериментальные образцы
охарактеризованы сотнями снимков в отраженных электронах, десятками
ИК
спектров, многими сотнями анализов на рентгеновском
микроанализаторе.
Расчетные
результаты
получены
только
с
использованием программ автора.
Построение диссертации
Работа состоит из введения, 8 глав и заключения, изложенных на 183 стр.,
включает 75 рис., 6 табл., список литературы из 226 наименований.
Аппробаиия работы
Результаты работы докладывались на международных конференциях
в Москве (1991), Пизе (1994), в Страсбурге (1995), в Гейдельберге (1996),
в Ницце (20(Ю), в Стэнфорде (2002). На всесоюзных совещаниях в
Черноголовке (1998), в Москве (1999-2004), в Миассе (1991). Материалы
работы использовались в лекциях для студентов университетов г. Пиза
(2002), Флоренция (2002), Мадрид (1999).
По теме диссертации опубликовано более 40 работ с тезисами,
включая 6 полноценных публикаций в рецензируемых международных
журналах (Physics of the Earth and Planet Inter., Earth and Planet. Sci. Lett., J .
Volcanol., Geothenn Res., BuU.Volcan., Tectonophysics).
Исследования в России были поддержаны неоднократно грантами
Р Ф Ф И , в Германии - DFG, в Швеции (Упсальский университет) - K V A , в
США (Университет Сев. Дакоты) - грантами DOE и EPSCOR.
Начальный
интерес
к
динамике
геологических
процессов
сформировался в студенческие годы во время занятия поисками рудных
месторождений методами газовой съемки в бескрайних степях Казахстана.
Первым конвективньпл процессом, с которым я столкнулся благодаря моим
руководителям к.г.-м.н. С.А.Воробьеву и проф. А.П.Соловову, было
перемешивание
газа
в
приповерхностных
пористых
породах
стохастическими колебаниями атмосферного давления. Решение в
случайных функциях для распределения давления и скорости фильтрации с
глубиной до сих пор является одной из самых полезных и интересных моих
работ. Переход аспирантом в Институт экспериментальной минералогии
АН СССР после недолгой работы практического геохимика приблизил
меня к расплавам, диаграммам плавкости, к динамике фазовых переходов.
На начальном этапе знакомства с веществом при высоких температурах и
давлениях большую пользу принесло общение с профессионалами из
лаборатории
магматизма:
с
А.Кузнецовым,
А.Чехмиром
и
М.Б.Эпельбаумом. Позднее К.И.Шмулович приобщил меня к вулканам и
дегазации магмы. Основным же соавтором в экспериментальной части
работы, которая приложила руки к тому, чтобы пузыри и гавайиты стали
так близки и понятны, является Т. Салова. А все основные теоретические
результаты диссертации были получены благодаря сотрудничеству с
членом-корреспондентом РАН В.П.Трубицьшым. Знакомство с ним стало
решающим шагом, без его высоких требований и превосходно развитого
физического чутья мне не удалось бы добиться успеха. За помощь при
обсуждении научных результатов автор благодарен академику РАН
В.А.Жарикову,
долгое
время
руководившему
Институтом
экспериментальной минералогии РАН и внесшему огромный вклад в науку.
Нельзя не отдать должное заграничным коллегам. Ни разу не побывав на
вулкане Этна, я благодаря сотрудничеству с проф. Pietro Armienti стал
одним из основных его исследователей методами экспериментальной
петрологии. Во время многочисленных визитов в достопочтенный
Пизанский университет, в котором прославился Галилей и работал
Э. Ферми, довелось наблюдать гавайиты на экране электронного
микроскопа и в окуляры FTIR микроанализатора с самим Pietro. Очень
полезны были беседы с его коллегами вулканологами M.Rossi, M.Pollacci,
P.Papale и их мудрым шефом Prof. Fabricio Innocenti. Изначально туманные
мысли автора по поводу сидементационной конвекции были основательно
упорядочены в ходе совместной работы с Prof. Нагго Schmeling - сначала в
Байроте, а затем во Франкфурте-на-Майне. Он же привил мне вкус к
использованию маркеров для численного моделирования адвективного
переноса. Мир гидротерм открылся мне в Грэнд Форксе, Северная Дакота в
общении с Ahmad Ghassemi. Оказалось, что аппарат, который мы с
В.П.Трубицьшым использовали для расчета отчасти мистического
пограничного слоя застывающей магматической камеры, с небольшими
модификациями годится для решения задачи об одномерной гетерофазной
конвекции рассолов в корневых частях геотермальных полей.
ОСНОВНОЕ СОДЕРЖАНИЕ РАБОТЫ.
Работа состоит из экспериментальной и теоретической частей. В
экспериментальной части (глава I и П) представлены результаты
оригинальных исследований процессов йузырения и кристаллизации
магматического расплава. Получены важные данные по нуклеации и
формированию распределений по размерам пузырей в кислых расплавах.
Кристаллизация основных фаз гавайитового расплава изучена при
постоянном переохлаждении и при сбросе давления и дегазации при
верхнекоровых давлениях. Результаты использованы для интерпретации
условий кристаллизации гавайитовой магмы вулкана Этна (Сицилия).
Теоретическая часть посвящена решению системы дифференциальных
уравнений, описывающих охлаждение, рост и оседание кристаллов в
пограничном слое, в равновесном и кинетическом приближениях (главы IIIV). Особо рассмотрена конвекция в гетерофазной магме, инициированная
потоком кристаллов или пузырей (главы VI-VII). Методом граничных
интегралов решена задача о стационарной конвекции, найдены критические
условия для ее начала. Модели гетерофазной конвекции применены к
описанию корневых частей геотермальных систем и конвекции, вызванной
пузырением магм в камерах под вулканами (глава VIII).
I.
Физика пузырения магмы.
Фазовое превращение начинается с зарождения новой фазы в среде
неравновесной, старой фазы. Зарождение и последующий рост с
конечными скоростями в условиях закрытой системы, т.е. без больших
относительньк перемещений фаз, может отвечать условиям быстро
протекающих дегазации и кристаллизации магмы в субвулканических и
вулканических условиях. Динамику этих процессов можно изучать
экспериментально в лабораторном масштабе времени.
Гомогенное зарождение. Основой теории зарождения является
представление о важном вкладе поверхностных сил, межфазовой энергии в
энергию вещества в субмикронных объемах. Повышенная энергия
вещества с большой удельной поверхностью способна компенсировать ту
разницу химических потенциалов, которая является движущей силой
фазового
превращения.
При
некотором
критическом
размере
устанавливается неустойчивое равновесие. Для критического пузыря
водного флюида справедливо соотношение между приращениями
химических потенциалов растворенной и выделившейся воды при сброседавления от равновесного Ро до текущего Pi:
Д//."^''=Д//^, Р„=*Р,
(1.1)При сбросе давления на несколько сотен бар парциальный мольный объем
растворенной в расплаве воды меняется незначительно, а флюида
8
значительно. Проинтегрировав, соответствующие члены по давлению и
учтя вклад поверхностной энергии пропорциональный кривизне пузыря,
получим
^.„-УЛР.-Р.),
.,=±,.,,.^Ж1Ш^^^^
(2.1)
Откуда нетрудно выразить размер критического зародыша через
характеристический капиллярный размер Ro и относительное падение
давления Pi/Po
"
'
(l +
При
— "
rf(l-e)(;c-l))""-'"-x'
размере
кластера,
^-YJSE^R_^=1L^X
меньшем
=^
критического,
Н20 .^ Н20 тг
(3.1)
его рост
термодинамически невыгоден, поскольку ^l п>(1 „,- Когда критический
размер оказывается достигнут в результате случайных флуктуации,
начинается стабильный рост пузыря. В рамках классической теории
нуклеации скорость гомогенной нуклеации выражается через пересыщение
флюидного компонента в расплаве, его диффузионную подвижность и
размер критического зародыша.
Анализируя выражение для гомогенной скорости нуклеации нетрудно
получить оценку времени задержки при этом процессе. Условно ее можно
принять равной времени, за которое образуется, хотя бы один пузырь в
характеристическом объеме расплава в 1 см^ за обозримое время от 1 сек
до 1 недели. Ввиду экспоненциального характера зависимости скорости
роста пузырей от пересыщения, пороговая величина сброса давления (или
эквивалентного пересыщения), необходимая для
начала гомогенной
нуклеации, слабо зависит от принятой величины максимального времени
ожидания. Наши оценки минимального понижения давления близки к
экспериментальным (Mangan and Sisson, 2000). Критический порог
пересьпцения свойственен гомогенной нуклеации не только при
пузырении, но и при других фазовых переходах. Так при росте кристалла
существует первое критическое переохлаждение (от долей градуса до 10" в
зависимости от системы), выше которого становится возможен рост по
механизму двумерного зарождения центров роста на грани. При более
низких переохлаждениях рост происходит по спиральному механизму
(движению ступени выхода спиральной дислокации по поверхности
грани), являющемуся по существу механизмом самораспространяющегося
гетерогенного зарождения (Snnagawa, 1985). Аналогичное критическое
переохлаждение существует и для гомогенного объемного зарождения
кристаллов в расплаве (Donaldson, 1981). Существование порога
гомогенного зарождения может проявляться в нелинейных эффектах при
кристаллизавди расплавов (Hort, 1995).
Дополнительной особенностью гомогенного зарождения пузырей
флюидной фазы из расплава является минимальное начальное содержание
воды, при котором возможно гомогенное зарождение. Анализ уравнения
9
скорости гомогенного зарождения показывает, что при начальном
содержании воды в кислом расплаве 2-3 мас.% (оценка зависит от
экспериментальной точности определения многих параметров) или
давлении насыщения 250-500 бар сброс давления до атмосферного не
обеспечивает условия начала гомогенной нуклеации. Для основных
маловязких
расплавов
критическое
начальное
содержание
для
осуществления гомогенного зарождения пузырей ниже и составляет около
0.5мас.%. Если отсутствовали условия для гетерогенной нуклеации,
магматическая вода может сохраняться в закаленных вулканических
стеклах.
Роль вязких напряжений. В
последнее время в вулканологии
большое
значение
стало
придаваться
явлению
фронтального
разрушения
закристаллизованной на 80-90%
магмы экструзий с частичным
расплавом,
насыщенным
летучими компонентами под
давлением на несколько сотен
Рис.1. Фронт пузырения, распространявшийся бар
выше
атмосферного
от внешней границы, высота микрофото 80 (Aldibirov
Dingwell
1996).
*"™Экспериментально
было
установлено, что при внезапном контакте пористой среды с атмосферой (в
результате обрушения конуса) дегазация и хрупкие разрушения композита
начинаются с поверхности и сопровождаются послойным отстрелом
обломков. В больших масштабах это может привести к катастрофическим
эксплозивным
извержениям
типа
Сан-Хелен.
Аналогичное
вышеописанному процессу фронтальное пузырение от поверхности
расплава наблюдалось нами в эксперименте (см. рис.1). Фронт пузырения
распространялся на некоторую глубину от поверхности и останавливался.
В экспериментальных образцах наблюдаются повторяющиеся полосы
мелких субмикронных пузьфей, параллельные внешней границе. Скорость
распространения фронта бьша достаточно велика, поэтому нам удалось
зафиксировать только продукты этого процесса, а не его стадии. Известно,
что скорость роста пузьфей на ранней стадии лимитируются вязкостью
расплава, а на более поздней - диффузионным подтоком воды из
окружающей среды (Bagdassarov et al., 1996). На «вязкой» стадии роста
можно принять, что давление флюида внутри пузыря отвечает его
содержанию в расплаве, т.е. приблизительно равно начальному давлению
Ро, при котором расплав был насыщен, а давление в расплаве на удалении
равно внешнему Pi< Ро. Рост пузырей описывается в приближении
разбиения объема смеси на сферические ячейки с пузырем внутри. Нами
10
показано, что при этих допущениях давление в расплаве Рт, равное
симметричной части тензора вязких напряжений, не зависит от радиальной
координаты (расстояния от поверхности пузыря) и составляет
р„-РЛ\-^~^),г,=2,г
'о
"о
''l
м'=р,-р,
(4.1)
*"
где г - радиус пузыря, rj - внешний радиус оболочки, а (ri"-r") - величина
пропорциональная количеству расплава, приходящемуся на один пузырь.
Значение п=2 отвечает плоскому решению, а п=3 - сферически
симметричному. Из формулы (4.1) вытекает, что при некотором
относительном размере пузыря давление становится нулевым и даже
формально отрицательным при дальнейшем его увеличении. В этом случае
могут быть достигнуты условия гомогенной нуклеации в условиях
растяжения при микро-разрывах расплава даже вопреки классической
теории нуклеации. Вязкость кислого расплава при содержаниях воды в 41.5 мас.% и при ликвидусной температуре составляет 10*-10' Пасек,
поэтому он проявляет при больших скоростях деформации вязко-упругие
свойства и может трескаться. Вследствие механических эффектов
пузырение может распространяться как фронт горения, когда множество
растущих мелких пузырей создает условие для зарождения следующих.
Приближение к свободно деформируемой поверхности облегчает условия
возникновения зоны пониженного давления. На внешней границе фронта
пузырения иногда наблюдаются деформированные пузьфи с трещинами
растяжения на концах. Вблизи поверхности слои пузырей деформируются
в пологие складки малой амплитуды. Описанное явление пока никак не
исследовано другими авторами, которые лишь констатируют влияние
вязких деформаций на форму маленьких пузырей в окрестностях больших,
начавших рост раньше (Larsen and Gardner, 2000).
Наличие гетерогенных центров зарождения в расплаве значительно
облегчает дегазацию, делает процесс пузырения более регулярным и
однородньп^! по объему. В качестве центров выделения флюидной фазы
могут выступать поверхности раздела кристалл/расплав. Микро-пузыри
менее растворимых, чем вода, флюидных компонентов, таких как
благородные газы, азот, углекислота, также могут выступать центрами
роста водных пузырей. Причем в силу зависимости химического
потенциала воды от ее доли во флюиде, микропузыри также должны
преодолеть барьер по размеру, прежде чем дальнейший рост станет
термодинамически выгодным процессом. Но вероятность преодоления
барьера гораздо выше, чем для гомогенного зарождения, поэтому
становится возможной перманентная нуклеация или несколько актов этого
процесса при умеренных пересьпцениях.
Измерение распределений пузырей по размеру. При декомпрессии
расплава пузыри перманентно зарождаются и растут, что приводит к
формированию
распределений
пузырей
по
размерам. Под
11
(дифферейциальным)
распределением
n(R)
(размерность
ед./см'')
понимается функция, связанная с интегральным распределением по
размерам Ф(Я), равным числу частиц размера не больше R в единичном
объеме магмы:
я
|л(г)Л- = Ф(Л)
о
(5.1)
Распределений пузырей по размерам, рассчитанные Торамару (Toramaru,
1995) по классической модели гомогенной зарождения, имеют ярко
выраженный экстремум, отвечающий одному акту гомогенного
зарождения. Распределения по размерам
„
пузырей в продуктах быстрой дегазации
3W1
магмы в природе (Mangan and Cashman,
"
-о-SIM
1996; Simakin et al., 1999) описывается
экпоненциальной или более сложной
' ■ ' \ / \
зависимостью, практически монотонно
убывающей (за исключением локальных
колебаний) с ростом ]размера пузырей.
я
m
w
Наблюдаемая разница в распределениях
свидетельствует о том, что природные
Рис.2.
Экспериментальные магмы гетерофазны и нуклеация пузырей
распределения пузырей в гранит­ носит гетерогенный характер.
ном расплаве по размерам в после­
К
началу
наших
исследований
довательные моменты времени
(Симакин
и
Салова,
1998)
при сбросе давления.
экспериментальных работ по изучению
динамики пузырения на уровне распределений по размерам практически
не было. И в более поздних работах (Gardner et al., 1999) приводятся
гистограммы числа пузырей, установленных в двухмерных сечениях с
изображений шлифов. Между тем истинная функция распределения
находится при пересчете на объем (стрейкоррекция), а без этой операции
кажущаяся плотность крупных пузырей' завышается, и гистограммы с
экстремумом отличаются от монотонных распределений (пример
распределений на рис.2). Нами экспериментально изучено пузырение
гранитного расплава при сбросе давления с постоянной скоростью
примерно за сутки с 1000 до 500 атм при Т=750°С. Дегазация протекала в
режиме гетерогенной нуклеации на микропузырях, наполйенных
преимущественно азотом, захваченных расплавом при наплавлении.
Состав газа в микропузырях оценен путем наблюдений на криостолике,
которые показали, что температура кристаллизации заполняющего их газа
составляет -131°С. Измерены распределения по размерам пузырей
флюидной фазы в стеклах, полученных при последовательной закалке
серии образцов в ходе нескольких экспериментов. Установлено, что при
дегазации объем выделившейся фазы отстает от равновесной величины. В
зависимости от содержания и размера воздушных пузырей в исходном
12
■
i;
стекле объем выделившейся фазы при
500 атм. либо
приближается к
равновесному (27 об.%), либо остается
вдвое меньше (13 об.%). Сигмаобразная
ол) форма
зависимости
объема
выделенного флюида от времени (см.
Рис. 3) с максимальным отклонением от
равновесия в промежуточной стадии
дегазации отмечается также и в работе
(Gardner et al., 1999). Интегральное
объемное содержание числа пузырей
(размером свыше 5 мкм) в ходе
дегазации существенно не меняется,
Рис.3.
Объемное
содержание
оставаясь на уровне 3-6 10* една мм^
пузырей в образцах в зависимости от
расплава.
Некоторое
возрастание
давления,
рассчитанное
по
распределению пузырей поразмерам. объемного содержания числа пузырей
Сплошная линия объемная доля на единицу объема стекла можно
интерпретировать как гетерогенную
нуклеацию со скоростью (1-3)10^ см"^ сек''. Конечные распределения могут
отвечать по форме экспоненциальному, степенному (фрактальному) или
унимодальному с максимумом в зависимости от объема и распределения
микропузырей воздуха в исходном стекле (расплаве).
Фрактальное распределение в генеральной совокупности. Несмотря на
наблюдаемое разнообразие и эволюцию распределений по размерам в
индивидуальных опытах, совокупность распределений из различных
экспериментов и природных данных укладывается в одну генеральную
совокупность. Для сопоставления были использованы распределения по
размерам пузырей (bubble size distribution, BSD) природных образцов,
выбранных из продуктов эксплозивной активности Mt. Etna и Vulcano
(Aeolian Islands, Italy).
Образцы были отобраны с точки зрения представительности процесса
дегазации при быстрой декомпрессии магмы в магматических фонтанах
или других формах вулканической активности. В образцах гавайитов Этны
и латитов с Вулкане процессы течения не оказывали серьезного влияния на
распределение и форму пузырей, которая была близка к сферической. Все
измеренные BSD имеют явно выгнутую кверху форму с выраженным
ростом плотности распределения пузырей размером менее 0.4-0.06 мм в
зависимости от образца. Для сравнения на общую диаграмму в
билагорифмических координатах (lg(n), 1п(г)) нанесена совокупность
экспериментальных BSD (из работы Simakin et a l , 1999).
030
13
о
<
.'•Ч;
I >
VUdOno рыП№«
С
& Е в 1 Efrx "969 « n e t
SeCKrfi I f ine i B o w ?
в
GtS3
Scor-y. г £л1 'икад
•
»
и з 52
в е к ^ B f w )9»в ftr,«t
£Т2в»
EtfwiIMe*»^»!
~ " ~- —
■';
oat
Lfe4(1ueB«'.«it
fclt«10S!*Pi^i
ein« M ' i W i
Wbota •j(;>e1rtW*a daw *«•
"
""
1
•''<s>jw^
•.J'fc
-its
>-4.
1 ^ о (i> л fin4)
Рис. 4. Сводка данных по экспериментальным и природным распределениям
пузырей по размерам (Simakin et al., 2000)
Генеральная совокупность описывается линейной функцией с
коэффициентом корреляции 0.92 и наклоном к= -2.8+0.1 (см. рис.4).
Некоторые природные BSD строго следуют этой зависимости. Другие,
содержащие крупные пузыри (с диаметром более 1 мм) имеют такой же
наклон, но сдвинуты в область меньших значений плотности. Несмотря на
особенности индивидуальных BSD, вся совокупность описывается общей
зависимостью в диапазоне размеров пузырей от 5 мкм до 1 мм. Найденная
степенная форма распределения по размерам возникает в различных
геологических и геофизических системах (Turcotte; 1992). Степенной закон,
вообще говоря, может отражать фрактальную природу некоторого процесса
физически инвариантного в пшроком диапазоне размеров. Значение
показателя экспоненты к=-2.8 является типичным для процедуры
заполнения пространства сферами, инвариантной относительно линейного
масштаба. Эта идея, высказанная нами (Simakin et al., 1999), получила
развитие в работе группы Н. Marder из университета г. Бристоль (Англия).
В результате численного моделирования процесса зарождения и роста
частиц при фазовом переходе Blower et al., 2002 показали, что при числе
актов нуклеации более пяти распределение частиц по размерам
приближается
к
степенному
фрактальному
с
к=2.3.
Другим
пространственно инвариантным процессом, который предположительно
может привести к фрактальному распределению по размерам пузырей,
является каскадное слияние всплывающих пузырей разного размера.
Согласно ОаопасЪ et al.(1996) в этом случае распределение асимптотически
приближается к степенному с к=4.99. Близкую размерность (к=4.11) мы
обнаружили у B S D воздушных пузырей в исходном стекле, вероятно
сформировавшихся при коалесценцияи микропор при наплавлении из
порошка.
14
II. Кристаллизация гавайитового расплава сопряженная с дегазацией
Декомпрессия магматического расплава, во многих случаях
адиабатическая,
без
значительных
теплопотерь,
сопровождается
выделением летучих компонентов после достижения уровня насыщения. В
свою очередь потеря летучих компонентов, как правило, ведет к
кристаллизации расплава. Это связано с тем, что вхождение летучих
компонентов, особенно воды, ведет к понижению температуры равновесия
кристалл - расплав (ликвидусной температуры) за счет понижения
активности соответствующих минералов в расплаве. Особенно сильное
влияние вода оказывает на температуру равновесия плагиоклаз - расплав.
Выделение воды вызывает кристаллизацию расплава, даже при неизменной
или растущей за счет теплоты кристаллизахщи температуре. Явление
кристаллизации,
сопряженной
с
дегазацией
было
изучено
экспериментально на примере расплава гавайитового состава.
Гавайиты Этны. Гавайит является умеренно щелочной вулканической
породой. Гавайитовые лавы - характерный продукт вулканизма во
внутриплитных условиях и в горячих точках. При извержении самых
прод)тстивных в мире вулканов Килауэа (Гавайи) и Этна (Сицилия)
преобладают гавайитовые лавы. Самый большой активный вулкан в Европе
Этна расположен в локальной зоне растяжения (тектоническом окне) у
границы Апеннинской аккреционной призмы (Doglioni et al., 2001). Его
практически непрерывная активность и необычайно высокое содержание
летучих компонентов объясняются механизмом генерации магмы с
участием флюидного компонента, высвобождаемого из погружающейся в
зоне субдукции Ионической плиты (Metrich et al., 1993). Эволюция состава
лав Этны от субщелочного к натровому щелочному и, со временем, к
щелочным породам калиевой спецификации (Schiano, 2001; Tonarini et al.,
2001) также связывается с участием коревого материала из погружающейся
Ионической плиты.
Проверка надежности расчета диаграммы плавкости по программе
M E L T S . До начала наших экспериментальных исследований диаграмма
плавкости гавайитов изучена очень незначительно, можно назвать лишь две
работы (Trigila el al., 1990; Metrich and Rutherford, 1998). Наиболее
плодотворным представляются не многочисленные определения условий
плавления гавайитов в широком диапазоне вариаций состава, содержаний
воды, летучести кислорода и давлений, а проверка достоверности
результатов расчетов по программе M E L T S (Ghiorzo and Sack, 1995).
Исходя из этого, были определены ликвидусные температуры главных фаз
гавайитов в условиях насыщения водой, характерных для вулкана Этна,
при 600 бар и содержании воды в расплаве 2.4+0.1 мас.% и летучести
15
кислорода на уровне NNO. Для определения была использована
прецизионная затравочная методика, при которой температура равновесия
пары минерал/р*асплав определялась по началу кристаллизации твердой
фазы с обрастанием затравки. Найдены температуры парных равновесий
оливина (1120±8°С), клинопироксена (1105±5°С) и плагиоклаза (1070±5°С)
Для плагиоклаза найденная таким путем температура является
температурой метастабильного продолжения ликвидуса, который примерно
на 20° вьпце температуры равновесного появления плагиоклаза при
равновесной кристаллизации гавайитового расплава после выделения
оливина и клинопироксена. Сравнение экспериментальньпс определений
ликвидусов главных фаз с результатами расчетов с M E L T S показало, что
ликвидусные температуры
клинопироксена и оливина в гавайитах
воспроизводятЬя с точностью ±5°, а расчетный ликвидус плагиоклаза
занижен, но не более чем на 20-25°. Дополнительно определена
температура ликвидуса клинопироксена при содержании воды в 1.5 мас.%,
которая также совпала с расчетной. Установлено, что программа M E L T S
не только хорошо описывает толеитовую систему, но успешно применима к
термодинамическому моделированию щелочных базальтов.
Скорости роста основных фаз гавайита.
Для вьфащивания
кристаллов из гавайитового расплава применялись затравки, т.е. кристаллы
минералов близких к интересующим нас по составу. Таким образом удается
обеспечить беспрепятственный рост изучаемых кристаллов до начала
массовой, объемной кристаллизации и установить их скорости роста. При
Т=1060°С (переохлаждение 15°) и Рн2о=600 бар скорость роста плагиоклаза
оценена в 1.410'^ см/сек, а время задержки гетерогенного зарождения
оказалось близко к нулю. При том же содержании воды и Т=1080°С
(АТ=25±5°) лейстовидный клинопироксен растет со скоростью 2-2.510'^
см/сек, а при содержании воды в 1.5 мас.% скорость роста оценена в 3-410'
* см/сек. Помимо этого установлено, что с ростом переохлаждения
происходит увеличение содержания примесей: глинозема в клинопироксене
и железа в плагиоклазе. Плагиоклаз затравок, а также новообразованной
при дегазации и перекристаллизации фазы содержит железа не более 0.60.7 мас.% (в пересчете на FeO). Лейсты, растущие при 1060°С (ДТ=10-15°),
содержат порядка 0.9-1.2 мас.% FeO. С ростом переохлаждения до 20-25°С
содержание железа растет до 1.5-1.7 мас.%. Максимальное содержание
железа до 3-3.5 мас.% наблюдается в закалочных лейстах. Знание
особенностей состава и порядка величины скорости роста кристаллов
позволило критически интерпретировать распределения по размерам
кристаллов клинопироксена и плагиоклаза из лав вулкана Этна.
Модели распределения кристаллов по размерам кристаллов. Для
изучения динамики кристаллизации и гетерофазной конвекции важно знать
размер флюидных пузьфей и кристаллов, поэтому, строго говоря,
моделирование
конвекции
с
одной
равновесной
оседающей
16
кристаллической фазой становится трехмерной задачей, где третье
измерение - это размер в функции плотности распределения кристаллов по
размерам. Самая простая постановка состоит в описании распределения по
размерам в однородной системе. Однородность достигается либо
отсутствием
перемещения
фаз
(локальная
однородность),
либо
перемешиванием объема жидкости, претерпевающей фазовый переход.
Тогда справедливо уравнение сохранения для функции плотности
распределения:
—^f-!-^——{y„n-Ddt
дг
^) = а
д г
,
(1.2)
г„
где член справа отвечает потере кристаллов из конвектирующего объема, в
закрытой системе а=0. Согласно (Martin and Nokes, 1989) скорость оседания
кристаллов из турбулентно перемешиваемого объема пропорциональна
Стоксовой скорости осаждения, а значит т = 2 . Граничным условием для
функции распределения является ее значение при нулевом радиусе, равное
отношению скорости нуклеации Jo к скорости роста ¥„
n(0)~Ji/V„.
Наибольшее применение получило упрощенное решение уравнения (1.2),
предложенное Marsh (1985). Оно описывает стационарное состояние в
пренебрежении диффузионным членом (D=0), а оседание, согласно Marsh,
описывается с а=1, ш=0. Принято, что скорость кристаллизации не зависит
от радиуса кристалла и является постоянной, как и скорость нуклеации.
Такая аппроксимация отвечает условиям кристаллизации солей из водного
раствора в техническом кристаллизаторе. Нетрудно установить, что
решением является экспоненциальная функция
«Г'-; = " „ е х р Г - - ^ Л
n,=J,/V^,
(2.2)
Поскольку экспоненциальной функцией часто можно аппроксимировать
природные распределения кристаллов по размерам, либо полностью либо
частично, формула (2.2) получила широкое распространение (Mangan, 1990;
Cashman and Marsh,1988). Временная константа TQ, неизбежная из-за
соображений размерности в уравнении (1.2), интерпретировалась Маршем
(Marsh, 1988) как эффективное время пребывания магмы в конвективно
перемешиваемом очаге в условиях постоянного притока из глубинного
источника. Задав по каким-то независимым соображениям величину
скорости роста, можно оценить То.
Другое тривиальное частное решение описывает кристаллизацию в
закрытой системе в пренебрежении дисперсией скорости роста и ее
зависимостью от размера кристаллов. Тогда решением уравнения (1.2),
удовлетворяюпцш граничному условию n(0,t)=J(t)/Vcr(t) является
и(г,0=Л'*)/^'„('*),
где время t* отвечает началу роста кристаллов радиуса х{г=\
(3.2)
V„{T)dT).
17
Это тривиальное решение, по существу, означает, что распределение
кристаллов по размерам рассматривается как запись отношения скорости
зарождения к скорости роста за все время кристаллизации данного объема
магмы. Как вытекает из анализа природных и экспериментальных данных,
для быстрых процессов кристаллизации вулканитов эта аппроксимация
более верна, чем модель Марша (уравнение (2.2)).
Распределение кристаллов плагиоклаза в лавах Этны по размеру.
Распределения кристаллов плагиоклаза по размерам в гавайитовых лавах
вулкана Этна хорошо изучены. Пример распределения по размерам (CSD)
кристаллов плагиоклаза из лав вулкана Этна, найденного путем анализа
оптических изображений шлифов, и прошедшего стереокоррекцию,
приведен на рис.5. Распределение представляет в логнормальных
координатах комбинацию практически линейного начального участка
(экспоненциального в нормальных координатах) при малых размерах
(менее примерно 0.6 мм) и плато при размерах более примерно 1.5 мм.
Сочленение происходит по негладкой кривой с несколькими резкими
колебаниями. Аналогичные распределения по размерам типичны для^
многих минералов, например, оливина (Annienti and Pareschi, 1992) из
вулканических пород.
С установленными скоростями роста плагиоклаза из гавайитового
расплава в условиях близких к природным можно оценить размеры,
которых могут достичь кристаллы плагиоклаза за счет кристаллизации,
сопряженной с дегазахщей магмы. Скорость подъема магмы в канале имеет
порядок 1м/сек (Trigila et al., 1990). С учетом того, что плотность
вспененной магмы невелика, гидростатическое давление в 600 бар
(условие насыщения по Metrich et al., 1993) отвечает глубине 2000-3000 м.
Поэтому порядок времени кристаллизации дегазируемой магмы может
составлять 2000-3000 сек. Приняв скорость роста плагиоклаза равной 1.5
-i-2 10'' см/сек, получим максимальную длину кристаллов порядка 2(1.5 -^2
10'* см/сек 2000^3000 сек) = 0.6-1.2 мм. Это очень приближенная оценка,
но представляется, что по порядку величины она верна, В этом случае
примерно экспоненциальный участок распределения кристаллов по
размерам отвечает декомпрессионному росту с увеличением скорости
нуклеации по мере потери воды. По неопубликованным данным П.
Армиенти при изучении кристаллов размером меньше 25 мкм (при
сильных оптических увеличениях) и еще более мелких (нанокристаллов) с
использованием S E M отмечаются квази- экспоненциальные участки
распределения с еще большими плотностями и склонами {dlg(n)/dR). При
столь малых размерах заведомо не работает механизм
оседания
кристаллов, предложенный Маршем.
18
Важно также рассмотреть состав
кристаллов
малых
размеров.
Самые малые, изученные на
рентгеновском
микрозонде,
i .1 \ д„
кристаллы
плагиоклазов
из
гавайитовых лав вулкана Этна
имеют состав Ап57-АПб4 (см. рис.5)
с
содержанием
железа
(пересчитанного на FeO) не
превышающим 1 % . Природные
Размер, мм
„ . „
'
кристаллы
содержат
много
Рис.5. Распределение кристаллов плагиоклаза
'^
из
гавайитов
Этны
по
размерам, ^««Ь'"« '«^леза, чем мелкие
неопубликованные данные P.Armienti (2002). закалочные лейсты И близки к
синтезированным
при
декомпрессии. По этому признаку можно заключить, что переохлаждение
(пересьпцение) природных магм по плагиоклазу не превышает, в
основном, 10-15°. Лишь эпизодически, при кристаллизации отдельных зон,
насыщенных расплавными включениями, оно, вероятно, было больше.
Таким образом, увеличение числа кристаллов, увеличение плотности
распределения к малым размерам, сопровождающееся смещением состава
кристаллов к альбиту, можно связать с падением содержания воды.
Потенциально, сопоставляя наклоны этих участков, а также сопоставляя
составы плагиоклазовых кристаллов, можно оценить скорости дегазации и
подъема расплава к поверхности.
Наиболее вероятно, что и экспоненциальные распределения по
размерам кристаллов из лавовых озер на Гавайских островах и из лав
вулкана Стромболи объясняются, как и на вулкане Этна, зарождением и
ростом кристаллов в квази-закрытых для седиментации условиях за
относительно короткое время.
Наибольшую ценность для оценки
скоростей обмена магмы в магматической камере под вулканами
представляет плато больших размеров.
19
Ш
Полное решение задачи о кристаллизации бинарной системы
эвтектического типа с оседанием кристаллов.
Первые две главы диссертации посвящены анализу кристаллизации и
дегазации
магматического
расплава
преимущественно
на
экспериментальном материале. Экспериментальные масштабы времени
порядка 3 (IHUO*) сек. отвечают масштабу времени самых быстрых
вулканических процессов в природе. В зависимости от состава вязкость
магмы варьирует от 10^ до 10* пуаз, причем часто при близликвидусных
температурах реология становится неньютоновской. Поэтому небольшие
по размеру кристаллы и пузьфи флюидной фазы, с радиусом менее 1 мм, не
успевают за это время переместиться на значимые расстояния, и
гетерогенную систему в объеме несколько кубических см можно, в первом
приближении, рассматривать как закрьпую. При более медленных
процессах, в частности при остывании расплава в камерах с линейными
размерами от нескольких десятков метров до нескольких километров
перемещения фаз становятся важными. За счет комбинации направленной и
объемной кристаллизации происходит магматическая дифференциация.
Наибольшее развитие этот подход по;^чил в работах М.Я. Френкеля и его
группы впоследствии (Арискин и др., 1987; Френкель, 1995). Однако
точных решений комбинированной задачи Стефана с направленной и
объемной кристаллизацией оседающей фазы до последнего времени не
бьшо. На первом этапе мы сохранили допущение, принятое у наших
предшественников, о существовании некоторой эффективной скорости
оседания кристаллов. Лишь в этом случае удалось получить аналитические
решения рассматриваемых проблем.
A-fl
»'»
> '
ч
Z
ч
V,
^.I't:
T
»
L
L
/B+L
A+L
А-^В
c. a
Рис.6. Схема пограничного слоя с
оседающей надэвтектической фазой (А).
Слева вверху распределение температуры,
излом отвечает границе Стефана - полного
затвердевания. Слева внизу упрощенная
диаграмма плавкости эвтектического типа.
20
Основные уравнения модели.
Нами рассмотрена для простоты
двухкомпонентная
система
(компоненты
А
и
В)
эвтектического
типа
с
кристаллическими
фазами
постоянного состава А и В. На'
фронте полного затвердевания
происходит
эвтектическая
кристаллизация фаз А и В, а
перед фронтом кристаллизуется
и оседает надэвтектическая фаза
А (см. рис.6). Перенос твердой
фазы можно описать, либо с
помощью
функции
распределения
по
размерам
n(R,z,t) кристаллов, растущих со скоростью V(z,t), с учетом зависимости
скорости осаждения от размера кристаллов Us(z,R)
El+^!lM:3t+y^,,tA = 0
dt
dz
(1.3)
dR
либо относительно ее третьего момента - объемной доли кристаллов £(zj)
и средней скорости осаяедения U,
|^a«t^^^^^^^^,^^|^^5^^r.
dt
dz
dR
dt
dz
(2.3)
p,
где Г - функция источника, определяемая ниже (уравн. (2.2.11)), р=4л/3
для сферических кристаллов, а усредненные параметры определяются как:
£(z,t) = 4л-/3^яЛ^Л
(3.3)
u:-u,.s.u,.^!pi^
(4.3)
Je " Л «ж
Уравнение сохранения для компонента расплава А, который участвует в
объемной кристаллизации, имеет вид
I-[(l.e)p,CJ+^[(l-£)p,Cu,J
= -r
(5.3)
dt
dz
где Г- функция источника, С - содержание компонента А в расплаве, U/ скорость расплава. Аналогично записывается уравнение для компонента В
j.[(l.e)p,(l-Q]+^[(l-e)p,il-QU,J
=0
(6.3)
где функция источника Гв=0, а концентрация компонента В равна (1-С).
Уравнение теплового баланса запишем в энтальпийной форме
^
=
-.^(еи,Н,+(1-е)и,Н,)
-^(е ил Н1-е )U,H,)++Я^у
Я 0
(7.3)
где Я- коэффициент теплопроводности, а энтальпия смеси Н находится
суммированием по фазам (кристаллам и расплаву).
Н= Н,е + Н,(1-е),
где
Н, «(с,,Г+ Я „ ^ л , Н, « {c^J+HJp,
(8.3)
Это уравнение теплового баланса для двухфазной среды может быть
переписано как уравнение для температуры (Симакин и Трубицын, 1995;
1997), если пренебречь разницей в теплоемкости и плотности жидкой и
кристаллической фаз
|^,,..; = „,0..МЬ
(„)
где Кт - коэффициент температуропроводности, L - теплота плавления.
21
Чтобы определить функцию источника генерации кристаллов Г,
необходимо еще одно уравнение. При достаточно большой скорости
фазового превращения система не отклоняется далеко от равновесия. В
равновесном состоянии температура расплава и его состав связаны
ликвидусными соотношениями.
В линейном приближении
С = 1-(С-Г)/ж,
(10.3)
где Тт это температура плавления (индекс m от melting) чистого
компонента А, а /я/ ликвидусный склон. Уравнение (10.3) определяет
количество твердой фазы, вьвделяющейся при понижении температуры, и
тем самым неявно определяет функцию Г.
Альтернативным является определение Г, использующее скорость V^r
роста кристаллической фазы А и величину ее поверхности Q. Допуская
независимость скорости роста от размера, получим для фазы с
распределением по размерам п(г).
r = p.jn(r)n(r)Vdr, Г = К(ДГ)
(11.3)
где скорость роста зависит от переохлаждения AT равного разности
ликвидусной и текущей телшературы:
ЦС)=Т/.(1-С)т,
(12.3)
Таким образом, для шести неизвестных функций T(Z,R,t), C(Z,R,t), Ui(Z,t),
Us(Z,R,t), IJZ.t), £(Z,t) мы получили систему из шести уравнений; трех
уравнений сохранения массы (2.3), (5.3), (6.3) для двух компонетов ( Сд и
Св) расплава и содержания кристаллов е, уравнения переноса тепла (9.3),
уравнения фазового равновесия (12.3) и соотношения между абсолютной и
относительной скоростями оседания Us (4.3). В связи с гиперболической
природой уравнений седиментации возможны разрывы в распределении
кристаллов ф), на что указывал М.Я.Френкель. Подобные явления
встречаются в газовой динамике и при инфильтрационном метасоматозе,
что одним из первых описано Д.С.Коржинским.
Анализ решения рассматриваемой проблемы проведен в приближении
квазистационарного решения при заданном тепловом потоке на фанице
полного затвердевания. Нестационарные законы сохранения бьши
преобразованы с введением новых координат связанных с движущейся
границей затвердевания t'=t, z'=z-\Uydt. При этом производная по
времени преобразуется согласно выражению:
-'^-uA
22
так что скорость перемещения границы попадает в стационарное решение,
которое качественно зависит от соотношения скорости оседания 5
надэвтектической фазы А и скорости перемещения фронта затвердевания
Ub- После манипуляции с первыми интегралами уравнений сохранения
получается важное соотношение для распределения доли кристаллической
фазы А в зависимости от расстояния до фронта полного затвердевания
.(z) = .
'-^^
(l-{7,(z)/C/J(l-C(z))
(13.3)
где В - это константа интегрирования, которая подбирается для
удовлетворения соответствующего граничного условия для г. При
скорости осаждения больше Ub
е(0) = 0, С = С^, В = С,
(14.3)
Для поглощающего пограничного слоя (U/Ub<l) фаничное условие
необходимо ставить на внешней хранице
£(«) = 5„, С = С„, Д = С„ + е„(1-у)(1-С„)
(15.3)
где здесь и далее У=С/Д4. Простое алгебраическое соотношение (13.3) для
Е(С) позволяет изучить свойства распределений ф), C(z), не находя их
точного вида. Для завершения решения необходимо решить уравнение
теплового баланса с учетом (13.3). Величина скорости фронта полного
затвердевания также является параметром, который отыскивается, исходя
из граничных условий. Распределение температуры найдено в неявном
виде как 7=Ф(Т,В), где Ф достаточно громоздкое выражение, приведенное
в диссертации. Важно, что распределение температуры представляет собой
гладкую функцию, близкую к простейшему экспоненциальному
распределению, свойственному классической задаче Стефана без
объемной кристаллизации.
Распределение объемной доли твердой фазы может быть разрывным,
поскольку оно описывается гиперболическим уравнением адвективного
переноса. Оказывается, что большое значение имеет скорость противотока
жидкости в слое с оседающими кристаллами. При пренебрежении этим
эффектом из (13.3) вытекает, что при равенстве Стоксовой скорости
оседания и фронта затвердевания Qo-S/Ub=l) решение имеет особенность.
Оно не имеет физического смысла, так как в этом случае содержание
кристаллов у фронта (при jo<l) и на бесконечном удалении (при jo>I)
стремится к бесконечности при J о стремящемся к единице. Тогда как при
Ja>l все надэвтектические кристаллы фазы А, образующиеся в погранслое
оседают и не захватываются фронтом. Напротив, приУо<7 все кристаллы
захватываются фронтом затвердевания.
23
Учет эффекта противотока жидкости. Если кристаллы образуются в
тепловом погранслое, их Стоксова скорость обычно значительно не
превосходит скорости фронта затвердевания, и вследствие влияния
противотечения и возрастания вязкости можно ожидать, что в части
термального погранслояу</, а в частау>/.
Минимальный уровень скорости осаждения, при котором появляется
интервал с кристаллами, избегающими захвата, достигается при jo=l- В
работе показано, что первая критическая скорость осаждения Si,
отвечающая этому условию равна
S,.=
S
(16.3)
pmiC,(C„-C.+Ste)
Аналогичным образом, второй критический уровень Стоксовой скорости,
обеспечивающий условие j>l во всем интервале (Се,Со) достигается при^'о
=4.791 (при Со^О.б, Се=0.3, общее выражение громоздко, но получается
тривиально).
^.791.Q
(1-С
1
(J7 3)
рга.с (С -С +Stei!-b^)
' '
"
*
(1-CJ
Таким образом, при значениях скоросщ осаждения 5 сг (O.SJ все
кристаллы фазы А захватываются фронтом затвердевания. В этом режиме
скорость фронта затвердевания достигает своего максимального значения
равного Ub-Si. При увеличении S происходит частичное осаждение
кристаллов, сопровождающееся падением скорости фронта до уровня
S2/4.791. При S >S2 появляется стационарное решение, описывающее
осаждение всей надэвтектической фазы А. Смена режима от полного
захвата до полного осаждения происходит при изменении скорости
осаждения от уровня 5/ до S2=2.9Si (при Ste=2, Со=0.6, €^=0.3).
Эффективный размер кристаллов при этом изменяется примерно в 2 раза.
При S>Si возможно множество стационарных состояний с различными
распределениями объемной доли кристаллов и различной степенью
фракционирования. Теоретический анализ стационарных состояний был
дополнен численным решением уравнений модели. Численные решения
показали, что при малых и больших скоростях осаждения быстро
устанавливается соответствующие стационарные состояния. Тогда как при
скоростях осаждения, близких к скорости перемещения границы Стефана,
возможны решения с поглощающим погранслоем, который не покидают
оседающие кристаллы надэвтектической фазы. Также возможны
разрывные решения разных типов. Более реалистичная картина
распределения кристаллов получается с учетом их переменных размеров в
процессе роста. Она рассмотрена в главе V диссертации.
24
IV. Гидротермальная конвекция типа тепловой трубы
в двухфазном флюиде
Гидротермальная система с двухфазным (жидкость-газ) флюидом
очень близка по физическому механизму и способу формального описания
к разобранной выше кристаллизации погранслоя магматической камеры с
оседающими
кристаллами
надэвтектической
фазы.
Изучение
гидротермальной конвекции имеет растущее практическое значение в
современном мире. Приповерхностные гидротермы в значительной степени
используются для производства энергии в различных странах мира. В связи
с этим моделирование конвекции водного флюида на глубинах до 2-3 км и
температурах до 300° стало рутинной инженерной практикой (Pruess, К.,
code T O I F G H I ; Oldcnhurg, С , and Pruess. К., 1993). Вместе с тем, по мере
освоения практически всех доступных приповерхностных гидротерм в
Италии (Ландорелло - Тровалли) и Японии началась эксплуатация
горизонтов с глубиной 3-4 км (Barelli et al., 2000; Kasai et al,, 2000). Ha этих
глубинах обнарз^жены рассолы различного состава. При кипении рассолов
могут возникать условия для существования обширных двухфазных зон,
причем флюиды могут быть весьма агрессивны и способны к растворению
окружающих пород. Нами теоретически рассмотрены режимы одномерного
двухфазного течения в модельной реакционной среде, состоящей из соли. В
реальности реакционный флюид реагирует с силикатными породами и
способен выщелачивать отдельные окислы металлов, оставляя инертные
компоненты (АЬОз и в меньшей степени SiO^). В природе система флюидпорода может быть далека от равновесия. Тем не менее, полученные
оценки скорости тепло- и возможного массопереноса в такой системе
представляют интерес для оценок масштабов природных процессов. С
помощью
модели рассмотрены закономерности
функционирования
реальных гидротермальных полей на Западе С Ш А .
Основные уравнения модели. Течение каждой из фаз двухфазного
флюида подчиняется модифицированному уравнению Дарси
Я,=-^—рд
Ml
oz
+p,g) ,
(1.4)
где i=l,2 номер фазы, kri(s) - относительный коэффициент
проницаемости в смеси в зависимости от объемной доли газовой фазы е, р,
плотность i фазы, qi - поток массы i фазы.
Закон сохранения энергии записывается в энтальпийной форме и
включает в себя тепловой эффект кипения
1.((Н,р,Е + Н^р/1-е))ф + Н,(1-ф))+-^{д,Н,+д^Н^)=Я^Т,
(2.4)
где ф - пористость породы, Н - энтальпия жидкости и газа, Л коэффициент теплопроводности. Законы сохранения массы выписываются
для каждого компонента в каждой фазе. Рассмотрен случай простой
25
системы
соль-вода,
явлениями
гидролиза
и
фракционирования
пренебрегают. Для растворенного вещества справедлива система уравнений
^(фер, С,)+|-(9,С,) = Г„ +Г,:
от
(3.4)
сг
1.(ф(Х-е)р,С,) + ~(я,С^) = -Т,,
(4.4)
где Cg и С/ являются содержаниями растворенного вещества в
газообразной и жидкой фазе соответственно. Го - функцией источника,
отвечающей процессу растворения и осаждения твердой соли, Л функцией источника, описывающей перераспределение соли между
газовой и жидкой фазами. Аналогичные уравнения справедливы для
содержания воды в жидкой и газообразной фазах. Условия равновесия
связывают химические потенциалы компонентов в фазах:
соли в газе и жидкости:
А„Й (С, ,PJ) = //„„ (С,, Л Т)
воды в газе и жидкости:
/ ' „ „ „ ( С „ Л П = //»„„(С,,Л7')
соли в твердой фазе и жидкости: ft'^{Р,Т) = Мш,(Q,Р,Т)
Относительные проницаемости фаз задаются в виде
k\=(s-e;r, к;=(1-^-^)",
(5.4)
(6.4)
где жидкость становится подвижной при объемной доле жидкой фазы
£>£>", а газ при 1-ё>е^. Для простоты часто принимается п=1, Sg=0, а'=0.
Граничные условия задаются не для всех переменных. В трехфазной
системе уравнения химического равновесия задают значения Cg, О, Р как
функции
температуры.
Для
температуры,
подчиняющейся
дифферешшальному уравнению второго порядка, необходимо два
граничных условия, таких как заданный тепловой поток и температура на
нижней границе кипения. Потоки qi, qg линейно зависимы и одна константа
определяет подток флюида через границу кипения снизу. М ы
рассматриваем закрытую систему, и некоторое количество воды попадает
снизу лишь для того, чтобы компенсировать осаждение соли.
Моновариантное равновесие в трехфазной системе. Нами изучено
квази-стационарное решение системы уравнений (1.4-5.4). Подобное
решение хорошо изз^ено для чисто водного флюида (Hardee, 1982;
McGuinness, 1996). Трехфазная система (соль-жидкость-газ) является
моновариантной и поэтому подобна чистому флюиду. Система уравнений
химического равновесия для трехфазного равновесия задает составы
газообразной и жидкой фаз, а также давление как функции температуры
(Cg(T), Ci(T), Р(Т)). Этого достаточно и для того, чтобы выразить плотности
фаз pi(T) и pg(T) от температуры (BischofF, 1991). Принципиальным
отличием солесодержащей системы является то, что помимо теплопереноса
в тепловой трубе происходит и массоперенос. Нетрудно отыскать
26
выражение для функции источника Го, задающей растворение и отложение
соли в процессе течения и связанной с пористостью колонны
g/ ^i-^g
Рис.7. Зависимость объемной доли
жидкости
от
температуры
при
стационарной гетерофазной конвекции в
трехфазной системе жидкостъ-пар-соль.
(7.4)
В
силу
однозначной
связи
давления и температуры при
трехфазном равновесии, можно
исключить давление, используя
уравнения
теплового
баланса
(2.4), Дарси (1.4) и с учетом
уравнений сохранения (3.4) и
(4.4). Из уравнений сохранения и
Дарси можно также получить в
явном виде выражения для потока
ОДНОЙ из фаз, скажем, qg, и
выражение градиента давления в
допущении
упрощенного
выражения для относительных
^.т^^''
^ " " " ^ '™"'"*^ " ° ' ' " ' ° ' ' «коэффициентов
проницаемости
жидкой и газообразной фаз
(к] =е,к[=\-е). В итоге получается выражение, связывающее долю
жидкой фазы и температуру в стационарном состоянии. Следует отметить,
что на топологию решения влияет тот факт, что зависимость Р(Т) имеет
экстремум при Т=596°С и заканчивается при T^SOO'C плавлением чистой
соли. Экзотические стационарные состояния е(Т) могут существовать и при
Т>596°С, но негативный наклон зависимости Р(Т) требует нагрева и
кипения сверху для этого типа решений, что не соответствует
геотермальным системам. Мы ниже рассматриваем решения для Т<Тти- На
рис.7 представлены линии зависимости е(Т), рассчитанные численно для
разных тепловых потоков через колонну. Этот график демонстрирует
важное свойство решения, имеющего две ветви, отвечающие газо­
доминирующей и растворо-доминирующеи системам. Аналогичное
свойство имеет и чисто водная система,, но потоки тепла, переносимые
флюидом с солевой нагрузкой, больше. При этом также происходит
растворение соли опускающимся флюидом и отложение ее на фронте
кипения.
Дивариантное равновесие в двухфазной системе. В рамках принятой
модели область двухфазного флюида не ограничивается трехфазной
кривой, а простирается между ней и критической кривой на РТ плоскости в
область больших давлений. В колонне двухфазного флюида в
дивариантном равновесии возможны процессы обмена солью между
газообразной и жидкой фазами, объемное кипение и конденсация. При этом
27
равновесие с солью достигается только на фронте кипения. В данном
приближении соль может попасть в систему лишь на фронте конденсации.
Формально в равновесном приближении это может происходить за счет
массообмена с вышележащим рассолом, находящимся в состоянии
термальной конвекщш. В реальности же возможно, что кислый
сконденсированный газообразный флюид при опускании будет
вьоцелачивать вмещающие породы.
В
приближении
двухфазного
равновесия
растворения
и
кристаллизащш соли не происходит, поэтому Го=0. На основании этого
легко вычислить стащюнарный поток соли в колонне
,,=-A^ = -ks \fy^
,
(8.4)
где В - константа, задаваемая потоком соли на границе колонны.
Используя это выражение и предположив для простоты, что вязкость фазы
пропорциональна ее плотности и =т,
легко получить выражение,
связывающее долю жидкой фазы с потоком соли в удобном виде
квадратного уравнения
5* =
= е(\-е)
(9.4)
(С,(ЛГ)-С,(Р.Г))(/7,(ЛГ)-рДР,Г))
Это квадратное относительно Б уравнение при заданных Р и Т, которое
имеет решение лишь при В*<1/4. Критическое значение В*= 1/4 отвечает
максимальному тепловому и массовому потоку через двухфазную колонну.
Одно решение с меньшим значением е отвечает газо-доминирующей
системе, а с большим корнем - растворо-доминирующей. На этом пути
оценен максимальный тепловой поток через колонну
АЯ,.^^ЛА-Р,)
Aft
При проницаемости 3.10"'*т^, Лр=0.% г/см^; ро = 1г/смЗ и прочих значениях
параметров, отвечающих системе, эта простая формула дает максимальное
значение теплового потока 4.3 ват/м^, близкое к численной оценке
максимального потока в трехфазной системе в 6 ват/м^. Этот поток можно
сравнить с потоком тепла при тепловой конвекции (Trubitsyn et al., 1993)
е = 0.218ДЯЖ^77, Ur = ^ , и, JJ^6^,
(11.4)
где L - вертикальная протяженность гидротермальной системы, Др(АТ) разность плотности, отвечающая перепаду температуры ДТ. Тепловой
поток в 0.9 ват/м^ получается при zir=200''C, д„=2.б г/см^, Ср^\ дж1т1о,
1=500 м, А:7=10'^ см^/сек, аг=Ъ-Ш\ /1=10"'^ м^ .
28
Прццожение
^
Fhncl
0»
"
L+V
0
№atpip<
0
Америки.
С
развиваемой
0
0
0
0
при
■dtssduttort
Она
рпащяимоп
точки
зрения,
скорость
вмещающих
достаточно
отвечает
пористости
пород
скорости
велика.
увеличению
примерно
0.&05 о б . % в год.
на
0.02-
Температура
предполагаемой
Liqiad
системам
найденной
конвекции
Huid
некоторым
модели,
рас'гворения
0
• ••.••>» ^ <»
Sealed zone
к
геотермальным
двухфазной
конвекции отвечает температуре
^
реологического
от
хрупкой
рис.
8)
.
перехода
к
пород
пластичной
А
это
(см.
значит,
что
Рис.8. Схематичное изображение двухфазной
линейная скорость компакции,
конвекции в корневой зоне гидротермальной
т.е.
' оседания
вьпнележащих
системы.
Ее
положение
совпадает
с
высокой
и
реологическим переходом во
вмещающих пород, ожидается
равной
н
е
с
к
о
л
ь
к
и
м
с
м
в
год.
породах при 350-400°С
Описываемый
процесс
можно
определить как высокотемпературный карст, к о т о р ы й , к а к и о б ы ч н ы й карст
карбонатных
Основные
пород
поверхностными
результаты
водами,
проявляется
этой г л а в ы докладывались
Стэнфордском
семинаре
Неожиданным
сюрпризом
по
геотермальным
оказалось
автором
системам
соседство
с
в
на
в
рельефе.
ежегодном
2002
докладом
году.
J . Moore
et
al.(2002), п о с в я щ е н н ы м к о н в е к ц и и с противотоком в двухфазном флюиде в
геотермальной системе Karaha-Telaga B o d a s (Индонезия) и подтверждающим
предложенную схему. О н и привели данные по ж и д к и м и г а з о в ы м в к л ю ч е н и я м
в кристаллах ангидрита, кальцита и флюорита, о б н а р у ж е н н ы м в минеральных
отложениях в породах, в м е щ а ю щ и х геотермальную систему. С а м ы е глубокие
образцы
ангидрита
эквивалента
NaCl.
содержат
Вся
рассолы
с
содержанием
совокупность
соли
минералогических
31
мас.%
данных
в
была
интерпретирована к а к проявление процесса- кипения и конденсации флюида,
сопровождающееся
погружением
кислого
агрессивного
конденсата
с
растворением в м е щ а ю щ и х пород. Д р у г и м явлением, вероятно с в я з а н н ы м с
конвекцией
агрессивного
флюида
на
глубине,
является
опускание
поверхности на многих геотермальных полях. Н а геотермальном поле Geysers
( С . К а л и ф о р н и я ) , по крайней мере, последние 20 лет, поверхность проседает со
скоростью 4.5±0.5 см/год ( M o s s o p et al., 1997). А н а л о г и ч н а я скорость оседания
3.5
см/год зафиксирована
(Wicks
et
отвечающим
al.,
2001).
на
геотермальном
Предполагается,
распределению
скоростей
что
поле Coso
(Ю.Калифорния)
источником
оседания, я в л я е т с я
деформации,
плоская
зона
контракции на глубине 4.6-3.9 к м в Coso и 4.1-3.8 к м в случае Geysers, ч т о в
обоих
случаях
отвечает
условиям
упруго-пластического
реологического
перехода.
29
V. Описание кристаллизации бинарного расплава с оседанием
кристаллов на уровне распределений кристаллов по размерам
Оперируя средними скоростями оседания кристаллических фаз в
задаче о застывании магматической камеры, М. Френкель с соавторами и
мы вслед за ними, возможно, чрезмерно упрощаем ситуацию. Размер
кристаллов, а вместе с ними и скорость оседания меняется при росте
кристаллов. Следовательно, от кинетики роста будет зависеть и
эффективность фракционирования, а также поток твердой фазы от верхнего
контакта остывающей интрузии, с которым может быть связаны
фавитационные течения.
Нами получено несколько решений задачи о кристаллизации
бинарного расплава, такой же, как в III главе. Для простоты мы
ограничились ростом одной надэвтектической фазы А перед фронтом
полного затвердевания. Для распределения кристаллов фазы А, растущих в
объеме, справедливо уравнение (1.3). Мы анализируем стационарное
решение при скорости продвижения границы затвердевания Ub
д>1(и.(г,Ю-и,)^удп
&
_
ал
.j^.
'
^' ^
где f/, = (1 - £■) ■ 5(Л, е), а Стоксова скорость осаждения равна
S = 2(p,-p,)gR'/9vp,
(2.5)
Объемное содержание надэвтектической фазы А в расплаве равно
£(2,1)=4л-/з|°л/г^Л
(3.5)
А функция источника может быть выражена через второй момент
распределения кристаллов по размерам
Г(2) = ^!гУ ^R^n(,R,z)dR
(4.5)
В стационарной модели приближенно можно положить, что в любой
момент времени кристаллы гетерогенно зарождаются только в узком
слое в окрестности нижней границы теплового погранслоя, так как
внутри погранслоя кристаллы уже возникли раньше, а в основном
объеме камеры температура еще несколько выше температуры
ликвидуса. Рост кристаллов на микровключениях начинается не
одновременно, поэтому в (R,z) пространстве решение для функции
плотности n(R,z,t) занимают некоторую полосу. Но для простоты в
настоящей работе мы рассмотрим предельный случай, когда при
понижении температуры ниже ликвидусной начинается одновременный
30
рост всех кристаллов. Е м у отвечает бесконечная скорость гетерогенной
нуклеации, и плотность распределения по размерам и пространству
wfK,z,?j имеет вид 5 - функции
r,{R,2) = S(R-R*{z))N{z),
(5.5)
где R'(z,t) - уравнение траектории одиночного оседающего кристалла, а
N ( z ) - плотность распределения кристаллов по пространству равная
первому моменту распределения по размерам
DC
*
*
|п(Л, z)dR = |<?(Л - Л (2)) • N(z)dR = N(z)
(6.5)
о
о
Траектория кристалла параметризуется временем / и может быть
представлена парой уравнений для Л и z . В движущейся системе
координат уравнение траектории записывается в виде
i^.VM. t-.u.-V,
(7.5)
В стационарном случае они сводятся к одному уравнению
dR' _ V{z)
(8.5)
dz ~U,-U,
В общем случае траекторий может быть несколько с учетом
возможности нескольких актов зарождения и возвратного (в
движущихся координатах) хода кристаллов, покидающих погранслой.
Подставляя (5.5) в (1.5) и учитывая свойства ^Л-/?*^-функции,
получим дифференциальное уравнение для изменения N(z) вдоль
траектории осаждения кристаллов (характеристик ур. 1.5)
^^t/,-t/J = - ^ ^ ,
(9.5)
dz
dz
где dUs/dz ~(dUs/dR)(dRVdz).
Таким образом, состояние погранслоя
задается распределениями четырех функций C(z), T(z), N(z), R(z). Первые
два являются решением равновесной задачи, вторые два находятся из
совместного решения пар дифференциальных уравнений (8.5) и (9.5) для
R*(Z)HN(Z).
В принятом приближении нуклеация принимается гетерогенной, а
распределение температуры и объемной скорости роста задаются
равновесным решением. Вычисляются размер кристаллов (R(z)) и их
интегральное объемное содержание (N(z), в ед./см ). Рассмотрены случаи
поглощающего пограничного слоя, когда скорость перемещения фронта Ub
велика и все кристаллы захватываются фронтом затвердевания. При этом
31
учитывался эффект противотока жидкости и влияние эффективной
вязкости, зависящей содержания кристаллов, на скорость их оседания.
Помимо этого использовано представление функции распределения в
виде суммы дельта-функций, отвечающих нескольким идеальным
дискретным актам нуклеации. Показано, что существуют решения с
частичным захватом кристаллов и их частичным осаждением. При
моделировании на уровне распределений по размерам эффект накопления
кристаллов, достигающих при своем росте скорости осаждения, равной
скорости перемещения фронта полного затвердевания, очень контрастен и
локален. В целом результаты наших исследований подтверждают
принципиальные выводы М.Я.Френкеля о том, что при дифференциации
силлов наиболее значимым является не «дождь» кристаллов от верхнего
контакта, а периодическая конвективная дестабилизация верхнего
погранслоя из-за накопления кристаллов.
Полученные решения не доведены до уровня непосредственного
петрологического приложения. Вместе с тем, подходы для дальнейшего
развития очевидны. В связи с ростом интереса к описанию процессов
затвердевания и дифференциации магматических интрузий на уровне
распределений кристаллов по размерам возможны обобщения на случай
нескольких кристаллизующихся фаз и учета фактических данных по
диаграммам плавкости и кинетики зарождения и роста.
32
V I . Стационарная седиментационная конвекция
Как рассмотрено в предыдущих главах, м,агматический расплав лишь в
редких случаях можно рассматривать как чистую жидкость. При снижении
давления происходит дегазация и сопряженная кристаллизация магмы,
точно так же, как и при ее охлаждении. В неоднородной смеси жидкости,
кристаллов и газовых пузырей градиенты плотности могут достигать
больших значений и вызвать гетерофазную конвекцию.
В том случае если можно пренебречь относительными перемещениями
фаз-включений в жидкости, течения в ней становятся эквивалентны
обычным течениям Рэлея-Тейлора. Конфигурация, когда легкий материал
помещен под тяжелым, абсолютно неустойчива. Наиболее хорошо
исследовано развитие неустойчивости Рэлея-Тейлора в среде с плоскими
горизонтальными слоями. Показано, что исходное синусоидальное
возмущение границы раздела двух слоев развивается по экспоненте,
геометрические и реологические константы слоев определяют длину волны
максимально растущей моды.
Линейный анализ устойчивости. Исследования устойчивости в общем
виде в слоистой среде с перемещающейся фазой, которая задает
стратификацию, не проведено. Трубицын и Харыбин, 1989 исследовали
частный случай гравитационной устойчивости горизонтального слоя
оседающих частиц. Проанализирована линейная устойчивость мгновенной
конфигурации с линейным уменьшением содержания частиц дисперсной
фазы книзу. Рассмотрим эти результаты подробнее.
Ниже представлена наиболее простая формулировка модели конвекции, а в
главе V I проведен анализ исходных уравнений гетерофазной конвекции.
Во-первых, это уравнение Навье-Стокса в приближении Буссинеска.
Гетерофазность фигурирует в нем только в гравитационном члене через
плотность смеси p=esP,+(l-s^pi
-VP+pg + V.\fi;\<?{U,)^W(U,Y]\=Q
(6.1)
Во-вторых, перемещение диспергированной фазы (твердой для
определенности) задается уравнением адвективного переноса
^+V-('(£/,+S)eJ = 0, S=S„
(6.2)
В-третьих, строго говоря, в гетерофазной среде справедлив аналог
уравнения несжимаемости с источником, который в приближении малого
содержания твердой фазы заменяется обычным уравнением несжимаемости
для скорости жидкости:
'V-u,=V-(Se.)»0
I
1
•
нос. НАЦИОНАЛЬНА." 1
БИБЛИОТЕКА
\
СПет«|4ург
ОЭ ЭМ «ВГ
I
'
(6-3)
..
П р и линейном анализе устойчивости рассматривается гетерофазный
слой толщиной Н, расположенный у верхней храницы столба чистой
жидкости. В качестве начального состояния принимается линейное
распределение содержания твердой фазы е, с глубиной z:
s\z,t=Q) = E,-eJH*2,
Z6(0,Z)
(6.4)
Слой оседает со скоростью S :
e = e\z-St,t), ye(St,H
+ St], e = e„.ye[0,St]
Возмущения по глубине концентрации (у) и вертикальной компоненты
скорости жидкости (Uy) с глубиной задаются в экспоненциальном виде
fr.UJ
= [Т(у), W(y)] ехрГД t) cos(far>
Уравнения развития возмущений бесконечно малой амплитуды имеют вид
^r+^-fr=0
(6.5)
dy
Ra,k'r + D'W = 0, D' = d'/dy'-k^,
(6.6)
где сидементационное число Рэлея определяется через характеристические
константы системы как
^^JPrP,)s,gH^^re,gH^
(6.7)
5//
Sv
Дисперсионный анализ зависимости Я(Ка,,к) п р и соответствующих
граничных
условиях
выполнен
полуаналитическими
методами
Трубицыным и Х а р ы б и н ы м (1987, 1989, 1991). И м и показано, ч т о
существует критическое число Рэлея, п р и котором
начинается
незатухающая сидементационная конвекция, т.е. Я>0. П р и непроницаемых
и скользких границах слоя
> = 0,Я:
U^=^^UJ^y^=0
ими установлено, что критическое число Рэлея составляет ^j*=104.76, а
волновое число критической моды равно к*=Ъ.О. Следует подчеркнуть, что
начинающаяся конвекция сосредоточена внутри слоя с градиентбм
содержания примеси. Поведение возмущения границы между этим слоем и
чистой жидкостью не рассматривалось.
Стационарная
седиментаиионная
конвекиия.
Л и н е й н ы й анализ
устойчивости свидетельствует о том, ч т о нет принципиальной разницы
между термической и седиментационной конвекцией. Термальная
конвекция происходит при сопряжении кондуктивного теплопереноса и
гравитационного течения. П р и условии, ч т о К а т < К а т " " кондуктивныи
теплоперенос ведет к диссипации возмущений плотности, и конвекция н е
34
развивается. При идеальной седиментационной конвекции (пренебрегая
гидродинамической дисперсией, которая сродни диффузии) возмущения
плотности рассеиваются адвективным переносом в частном случае
линейного неустойчивого распределения примеси в слое при Кажд<^Яжл"'Однако, строго никем не показано, что оседание дисперсной фазы способно
стабилизировать гравитационную неустойчивость типа Рэлея-Тейлора в
двухслойной системе. Тем не менее, и в этом случае существуют
критические явления при поддержании постоянного потока твердой фазы у
верхней или газообразной у нижней границы.
Для поиска условий существования и формы стационарного
седиментационного течения в двумерном приближении использован метод
функции тока. Тогда после масштабирования ( Я - масштаб длины, soмасщтаб содержания равный значению на границе, S - масштаб скорости)
уравнение (6.1) перепишется в виде:
A'(w) = -Ra.~-
(6.8)
дх
При постоянном потоке твердой фазы сверху и некоторых допущениях
внутри рассматриваемого прямоугольника размера (ах1) будут области,
содержащие кристаллы и чистую жидкость, разделяемые границей уь=/(х).
Тогда в правой части уравнения (6.8) появится производная от ступенчатой
функции Хевисайда
A'(v)=-Ra.S(y-f(x))ne.
(6.9)
Полуаналитическое решение линейного уравнения (6.9) можно
записать с использованием граничного интеграла по границе сред и
соответствующей граничным условиям функции Грина О(хь(0,уь0).х,у), где
время выступает как пгфаметр:
V'(x,y,t)= j
Q(xb{t),y,{t),x.y.)sm(arctg(^)dl
(6.10)
Приняв, что вертикальные границы симметричные, а горизонтальные
скользкие, получим выражение для функции Грина для области (а х 1)
С(оЛх.у,) = ^^ЪЪ^
яа
,/Г
.,т—)^^п^У>,
(6.11)
(пС/аЛ-ггУ
а
где A„''Sin(!^^^)sm(n/rd)
а
Теперь, когда поле скоростей известно, необходимо переписать уравнение
адвективного переноса твердой фазы в форме уравнения для границы
областей с твердой фазой и без нее. В рассматриваемом стационарном
случае:
^=^^^=r-|^+i;/|^.^=/rxj
скь
Uij,
Oxb
(6.12)
ay,,
35
Численно получены решения системы (6.11)-(6.12) для различных значений
„„ ,
седиментационного
числа
Рэлея.
Установлено,
что
стационарное
конвективное решение появляется при
достаточно медленной седементации,
отвечающей числу Рэлея Rased >105
при 0=0.7. Кристаллы, попадая в
верхний
погранслой
вдоль
всей
верхней
границы,
течением
захватываются в нисходящий поток,
который быстро достигает даш и
растекается вдоль него. Кристаллы
Рис.9. Зависимость седиментационного повидают систему у нижней границы в
числа Нусельта от седиментационного
числа Рэлея. Пунктир - параметры
нестабильного
включения
чистой
жидкости малого объема. Аспектное
отношение 1.
Результате
осаждения,
а
чистшг
жидкость поднимается в восходящем
потоке. Таким образом, получается,
что
течение
развивается
около
включения чистой жидкости в смеси.
Если включение изначально было недостаточно большого относительного
объема, оно погружалось на дно и весь объем заполняется оседающими
кристаллами. В работе определен седиментационный аналог числа
Нуссельта, отражающего эффективность конвективного переноса по
сравнению с простым оседанием. Показано, что, как и в случае тепловой
конвекции, число Нуссельта равно 1 при докритическом уровне
седиментационного параметра Рэлея, а после этого порога растет скачком,
как отношение полной площади ячейки к площади двухфазной зоны (см рис.
9).
Развитые формы седиментационнои конвекции расчитаны нами численно с
помощью программы, реализующей алгоритм, разработанный по методу
конечных элементов. Опицание алгоритма
представлено в главе V I I , пример умеренно
надкритического течения представлен на рис.
10.
Экспериментальные
данные
по
гравитационьш течениям в гетерофазных
средах. Прямой проверки
теоретических
оценок
условий
стационарной
седиментационнои
конвекции
или
Рис.10. Рассчитанное распределе-
КОНВеКТИВНОЙ
ние оседающих
НСуСТОЙЧИВОСТИ
В
СЛОе
кристаллов
в
градиентным
Восходящи поТо" ж и т
в
диспергированной фазы не проводилось.
конвективной ячейке при Raff" 500.
центре.
36
w ^
С
распределением
Известные нам эксперименты относятся к моделированию
гравитационного течения слоя, обогащенного подвижной дисперсной
фазой. Thomas et al. (1993). использовали двухслойную жидкую систему с
более вязким верхним слоем. Снизу равномерно по площади задавался
поток газовых пузырей. В силу неразрывности потока на границе при
замедлении пузырей возрастало их содержание, со временем зона с
пузырями в основании нижнего слоя расширялась. Авторы установили, что
при достижении этим пограничным слоем некоторой критической толщины
он терял устойчивость с образованием плюмов (диапиров) и вспльшал. В
качестве критерия наступления конвекции выступало равенство времени
развития неустойчивости Рэлея-Тейлора и характеристического времени
нарастания ширины слоя с пузырями. При малой относительной толщине
легкого слоя согласно Lister and Kerr (1989) максимальная скорость роста
неустойчивости а равна а, = 0.232^ ^ А характеристическое обратное
\dh
М
_S
Поэтому, согласно Thomas et al.,
h dth-k- h
1993 неустойчивость наступает при толщине слоя h, при которой 01=02.
Этот критерий равносилен введению локального числа Рэлея с масштабом
длины, равным толщине слоя. Критическое значение такого эвристического
критерия будет 4.31. На практике установлено лишь то, что расстояние
между вспльгаающими диапирами с пузырями примерно отвечает длине
волны самой быстрой ноды неустойчивости Рэлея-Тейлора (с точность 3040%). Похожее исследование провели Багдасаров и Дорфман (Bagdassarov et
al., 1996). Они изучили гравитационную неустойчивость оседающего слоя
платиновых частиц в частично расплавленном граните. Для ускорения
оседания в центрифуге создавалось ускорение lOOOg. Опыты проводились в
удлиненных (3:1) вдоль пути оседания ампулах. Авторы считают, что
переход от простого Стоксова оседания отдельных частиц к конвективным
движениям начинается при равенстве Стоксова времени, равного
отношению толщины слоя к скорости оседания, и характерного времени
роста развития неустойчивости Рэлея-Тейлора (RTI). Существенным
отличием оседания шариков в смеси кристаллов и расплава является
большая дисперсия скорости оседания, приводящая к дополнительному
механизму диссипации. Таким образом, стабилизация RTI понимается как
следствие недостатка времени для ее развития, а не в смысле теории
развития бесконечно малых возмущений. По порядку величины
седиментационного числа Рэлея (около 100) переход от простого оседания к
течению в опытах (Bagdassarov et al., 1996) происходит в области
критического значения, найденного нами.
время разрастания слоя а ,«
37
V n . Численное моделирование седиментационной конвекции
Полуаналитическими методами можно исследовать гетерофазную
конвекцию лишь в достаточно простых постановках. Для приближения к
геологически значимым случаям необходимо численное моделирование. В
работе развит метод для решения задач'^ гетерофазной конвекции,
основанный на известной технике маркеров. Дополнительные возможности
предусмотрены для описания дегазации и относительного перемещения
фаз, что не было сделано ранее в приложениях гидродинамики к наукам о
Земле.
Для решения уравнения Навье-Стокса в переменных давлениескорость был использован метод конечных элементов, в формулировке,
близкой к описанным в работах (Huges et al., 1979; Fletcher, 1984; Pelletier
et al., 1989; King et al., 1990; Флетчер, 1991). Для исключения давления
применяется метод штрафных функций (Huges et al., 1989) Давление
представляется по аналогии с уравнением состояния слабо сжимаемой
жидкости через дивергенцию скорости
p = -XV.U,
(1.7)
где А, - произвольный большой множитель со значением до Ю", После
подстановки давления из уравнения (5.1.5) в уравнение Навье-Стокса оно
становится функцией только скорости течения U
/7g + V(AV.C/)+V.r=0, r = /y(Vt/ + VC/'')
(2.7)
Для улучшения точности решения мы использовали алгоритм Узавы
(Fortin and Fortin, 1985; Pelletier et al., 1989). Этот алгоритм позволяет
получить поле скоростей с V • (/ = о даже для очень плохо обусловленных
проблем, когда вязкость жидкости меняется на много порядков. Давление
итерационно уточняется согласно модифицированной формуле (2.7) для
шага к
>o^+V(AV.[/J-Vi',,, + V-(//(V[/.+V[//)) = o
(3.7)
давление корректируется как
Р,=Р,.,-ЛУ-и,
(4.7)
Начальное приближение /*„ может быть принято равным О или, при
моделировании нестационарного процесса, обусловленного переносом
плотности, значению на предыдущем временном шаге Р„. Обычно для
достижения точности, близкой к машинному нулю, требуется 1-3
итерационных цикла.
38
Для аппроксимации компонент скорости на плоскости использовались
билинейные элементы, такие, что с каждым узлом j прямоугольной сетки
ассоциируется функция
N'=(a'jf-b',)(c)t}-d;)
где ^ и т) - нормированные координаты в элементе е (х и у, соответственно),
принимающие значения в интервале [0,1]. Коэффициенты в выражении
подбираются так, чтобы обеспечить Nj =1 в узле j и О в других узлах
элемента е. Давление описывалось с помощью пробных функций меньшей
степени, т.е. кусочно-постоянных на прямоугольном элементе функций.
Уравнения решались методом Галеркина в слабой формулировке. Эта
процедура хорошо известна и приводит к матричному уравнению
(К + Я M)U„ = «PCPj., ) + F
(5.7)
где К- матрица вязкости, М- матрица сжимаемости, F - гравитационный
член, а V(Pk-i) - член, отвечающий градиенту давления с предыдущего
шага итерации.
Для примера приведем формулу для элемента вязкой матрицы,
ассоциированную с элементом е
^ dN, дЫ„ ^ SN, <?N„
Kl= \n
^N, dN„
ск &f
dN, dN^
^
^ dN, dN„
d^ ^
(3c
dN, SN„,
дк ^
dS
Интегрирование гравитационного члена проводится с повышенной
точностью по сетке 3x3, поскольку распределение плотности описывается
с большим разрешением с помощью метода маркеров. Интеграл
представляется в виде суммы
\pgN,<m= Х 1»;;Ж,7,)5^,(#,'7/),
EL-m
(6.7)
М 3^.1 3
где веса равны произведению соответствующих коэффициентов при
численном взятии интеграла функции одной переменной:
Вектора точек интегрирования и весов в Гауссовом методе равны (
Zienkiewicz, 1983):
4 =л= [-0.7746, 0.0,0.7746]
W = [0.55556,0.88889,0.55556]
39
Решение уравнения седиментаиионного переноса с равным нулю
диффузионным членом вызывает технические трудности. Описание
переноса на Эйлеровой сетке связано с появлением неизбежной численной
диффузии. Для борьбы с численной диффузией применяются специальные
разностные схемы повышенного порядка. Для метода конечных элементов
разработана специальная модификация (Петрова-Галеркина), в которой
минимизируется численная диффузия в направлении, ортогональном
направлению течения. Но эти методы не исключают эффекта численной
диффузии полностью.
В качестве альтернативы выступает Лагранжев метод описания
переноса, уравнение
^ + Л У ( О Д =0
что эквивалентно
—^ = 0,
Dt
(7.7)
(8.7)
^
'
где плотность задана на характеристиках уравнения (7.7), а именно в
точках кривых y(t)=^(x(t)), таких что
± = и,, ^ = и,
(9.7)
Скорость перемещения для дисперсной фазы равна сумме скорости
жидкости и относительной скорости перемещения фаз С//+ S.
Проверка правильности работы программы. Конвективный код был
успешно протестирован для модификации, описывающей термическую
конвекцию. Дла проверки расчета массопереноса методом маркеров может
быть использовано решение задачи о течении со свободной деформируемой
поверхностью. Мы решили проблему расчета стационарного профиля
соляной экструзии в условиях медленного растворения сверху. Эта задача
была поставлена C.Talbot, пионером исследований растекания глубинных
солевьпс масс, выжатых в тектонические нарушения. Подобные соляные
«ледники» встречаются на юге Ирана, известны подводные экструзии на
дне Мексиканского залива. Классическое решение о нестационарной форме
растекающейся капли под действием силы тяжести получено X. Хаппертом
(Huppert, 1982). Когда соль регулярно растворяется с поверхности
выпадающими осадками и выдавливается снизу с постоянной (медленно
меняющейся) скоростью, возможно стационарное состояние. Профиль
двумерного потока (зависимость высоты h от расстояния до начала течения
х) можно описать простым решением, найденньпи нами методом тонкого
слоя:
Ь(х)=д/л/2(6-х), b = Q, 0<x<b,
(10.7)
где безразмерный подток соли задается через линейный масштаб L и
скорость растворения 5: Q = Q/SL.
40
Для
численного
решения этой проблемы мы
использовали
формально
гетерофазньш
подход:
солевая
жидкость
выдавливалась в невесомую
и
маловязкую,
но
Рис.11.
Профиль
соляной
экструзии
с несжимаемую
жидкость,
поверхностным
растворением
рассчитанный Д^^ расчета поверхности
численно с использованием модели двухслойной экструзии
применялась
среды с оседающими маркерами. Сплошная линия описанная выше программа,
аналитическое решение в приближении тонкого основанная
на
технике
^'^°'^■
маркеров,Полученный
численно при сетке 30x15 стационарный профиль (см. рис. 11) оказался
близок к теоретическому, описываемому формулой (10.7). Длина
стационарного течения определяется соотношением баланса, так как
растворение компенсируется оттоком через верхнюю границу и
отклоняется от точного значения менее чем на 1 % . Аспектное отношение
численного решения также соответствует аналитическому решению с
точностью до 1 % . Численное решение отличается по форме поверхности
Над питающей струей, посколько при аналитическом решении течение у
правой границы принято горизонтальным.
V I I I . Численные модели конвекции с пузырением в магматической
камере
Группа М.Я. Френкеля (позднее А. Арискина и Е. Коптев-Дворникова)
развила метод полуколичественйого моделирования динамики становления
пластовых интрузий (силлов). Оперируя скоростями оседания, как
подгоночными параметрами, и описывая конвекцию, как перемешивание в
заданном ими интервале глубин, авторы в деталях описали химическую и
минералогическую зональность изученных силлов. Своей деятельностью
они задали высокие требс(вания к следующему шагу в этих исследованиях.
Описание «онвекции должно быть как минимум двумерным, а
кристаллизации - на уровне распределений кристаллов по размерам.
Менее изученной, но более важной с практической точки зрения,
является проблема дегазации магмы в субвулканических камерах. Режим
дегазации, в конечном счете, определяет степень эксплозивности
извержения и его потенциальную опасность. Дегазация магмы в
вулканическом процессе может происходить как при подъеме ее в
магмоводе к поверхности и извержении, так и при длительной выдержке в
камере под вулканом. В первом случае она может достигнуть
катастрофического по скорости и последствиям характера, тогда как во
втором летучие (в основном вода с небольшим содержанием SO2, СОг, Нг и
41
др. (Symonds et al., 1994)) могут выделяться длительное время, достигая по
трещинам поверхности при температурах, близких к магматическим.
Стационарная высокотемпературная дегазация отмечается, например, на
вулкане Кудрявый, Курильские острова (где газы выделяются уже в
течение порядка 100 лет при температуре до 940''С) и на вулкане СатцумаИводзима, Япония (в течение более 1000 лет при температуре до 870''С
(Tkachenko et а!., 1992; Hedenquist, 1994).
Описание дегазаиии в модели. Для моделирования конвекции магмы в
условиях дегазации и сопряженной с ней кристаллизацией была
использована программа, описанная в предыдущей главе. Главной
особенностью, впервые реализованной в геофизических приложениях,
было описание в технике маркеров (характеристик) процесса дегазации в
приближение малой степени пузырения. Уравнения баланса для смеси
могут быть упрощены. Мы предполагаем, что давление в камере может
быть приближенно представлено гидростатическим давлением в слое
чистого расплава, что верно при малом содержании пузырей:
P = P, + P,gz
(1.8)
где Ро - это фоновое давление, которое определяется механическим
равновесием камеры во вмещающих породах. Латеральными вариациями
давления вследствие вязких напряжений мы также пренебрегаем, что
верно при мальпс числах Рейнольдса. На данном этапе Ро принимается
постоянным, регулируемым сопряженными механическими процессами,
такими как перетекание расплава и перемещение блоков по трещинам. Эта
формулировка является другим крайним случаем, если принимать во
внимание приближение постоянного объема. В последнем случае даже
незначительные перемещения пузырей или двухфазного расплава
способны вызвать неоправданно большое изменение Рд (см., например,
Michaelides, 1989; Bagdassarov, 1994). Однако расширение камеры и
деформация вмещающих пород могут нивелировать этот эффект.
С использованием уравнений баланса масс для компонентов
гетерофазного расплава, исключив дивергенцию скорости расплава (не
равную 0) и подставив давление в соответствии с гидростатической
аппроксимацией, получим
± = |:+(f/,+S*)V£ = r
(2.8)
где Г функция источника, а S* - абсолютная скорость перемещения
флюидной фазы
S* = S(\-s)
Г'П+П
42
(3.8)
(4.8)
(^,.5X1-.).
^5.8)
А.+г
tyW)__C^
_(1^£)_
(,g)
A„+z 2(1-C,)^ Д1 + ^А„)'
где A=pPa/pi, ho=PJptg. Безразмерный параметр A является отношением
характеристических плотностей газовой и жидкой фаз, /?„ - общее давление
в едйнрщах гидростатического эквивалента мощности камеры.
В (2.8) входит скорость перемещения газовой фазы 5* с учетом
эффекта противотока, что имеет место и в других задачах с относительным
перемещением фаз (Simakin et al., 1994).
Первое слагаемое в источниковом члене /} учитывает расширение
пузырей в гидростатическом поле давления. Уравнение (2.8) только с этим
членом имеет решение вдоль траектории E(Z)= eo/(ho+z). Второй член
описывает собственно пузырение.
При полном обороте частицы жидкости в конвективной ячейке
скорость жидкости Ui меняет знак, и интеграл от второго члена равен 0.
Перенос пузырей со скоростью 5* ведет к вертикальному транспорту
летучих компонентов и потере их через верхнюю границу.
Различные случаи конвекции с дегазацией.
Магматические камеры под вулканами обладают зональностью по
составу. Зональность проявляется в смене состава последовательно
извергаемых продуктов в одном извержении, зональности продуктов
палеоизвержений, таких как туфы. Анализ состава кристаллов и
расплавных включений в лавах также свидетельствует о широкой
распространенности явлений смешения магм различного состава в
магматических камерах. Стратификация магмы, особенно по содержанию
летучих компонентов, оказывает большое влияние на режим конвекции.
Нами численно промоделированы различные случаи конвекции в
стратифицированных по летучему компоненту объемах расплава.
Когда содержание летучих компонентов близко к насыщению, любые
перемещения расплава могут вызвать пузырение, усиливающее контраст
плотности и ускоряющее вызвавшее его движения. Положительная
обратная связь между композиционной конвекцией и пузырением может
быть причиной катастрофической неустойчивости. Численно прослежено
развитие неустойчивости слоя магмы, насыщенного во всем объеме, с
учетом изменения гидростатического давления с глубиной. Конвекция
инициируется за счет элемента расплава с пузырями, помещаемого в центре
прямоугольной области. В результате моментально формируется
центральное восходящее течение и нисходящее течение по краям.
Центральная колонна с пузырями постепенно растекается у кровли, область
43
вспузыренного расплава изменяет форму от античной вазы до рюмки. В
очень короткое время верхняя и нижняя части расплава меняются местами,
лишь малая часть придонного расплава задерживается, что обеспечивает
дальнейшую конвекцию. Моделирование проводилось в приближении
постоянства давления в магме у кровли, что подразумевает синхронное
увеличение объема камеры или удаление избыточного объема двухфазной
жидкости через дайки. С учетом прочности вмещающих пород быстрая
дегазация и увеличение
объема
камеры может
привести к
катастрофическому разрушению камеры.
Более
правдоподобной
является
двухслойная начальная
конфигураиия в системе с верхним
слоем, потерявшим воду и нижним
слоем, насыщенным на уровне
границы раздела слоев. Подобная
конфигура1Ц1я может сложиться при
внедрении свежей порции более
примитивного
(не
прошедшего
кристаллизационную
Г'^гГ7^Ш^-§-^''?ЩЩЩ,
дифференциацию)
расплава.
•^^^'' '"^''"^■'' ' '^'■'■■''^'^^^^^^
Внедряющийся расплав может быть
более горячим, более плотным, но с
Рис.12.
Развитие
конвективной большим
содержанием
воды,
неустойчивости в двуслойной системе. Конвективное перемешивание и
Интенсивность заливки отражает степень
пузыренш. нижнего слоя в восходящем падение давления в предыдущих
течении. Крушше маркеры отражают извержениях снижают содержание
положение жидкости верхнего слоя, мелкие воды в верхнем слое. На контакте
маркеры отвечают жидкости нижнего слоя верхний слой будет нагреваться, а
без пузырей.
нижний
охлаждаться
и
кристаллизоваться
с
дополнительным накоплением воды. Со временем плотности верхнего и
нижнего слоев могут сблизиться за счет кристаллизационной
дифференциации. Система становится неустойчивой, и любое возмущение
(в нашем случае «капля» тяжелого расплава) ведет к возмущению границы
раздела, пузырению. При развитии возмущения слои, как и в первом
случае, меняются местами, объем системы резко возрастает, что также
должно вести к разрушению камеры, как при извержении супервулкана.
44
Плавная дегазаиия магмы в
'какере под вулканом предполагает
специальный механизм конвекции
расплава. Казахая (Kazahaya et al.,
1994) предположил, что содержание
воды в объеме расплава меньше
уровня
насыщения.
Расплав
поднимается по каналу, по которому
происходят извержения, до уровня
насыщения и дегазации, где расплав
пузырится,
частично
кристаллизуется, теряет пузыри и
опускается в камеру. Скорость
РисЛЗ, Всплывание пофаничного слоя „отери летучих
ограничивается
между двумя слоями в стратифици<^г
' „ -^
„
^ ^
сечением
«трубы».
Не
рованнои магматической камере. Предпола-
'^•'
Tisi№, что более плотный нижний слой анализировалось
возможное
насыщен водой на уровне раздела, влияние скорости потери пузырей
Двухфазный
расплав
пересекает на скорость конвекции. В «трубе»
недонасыщенный
верхний слой
и степень дегазации может быть очень
растекается у кровли.
высока, а содержание пузырей
достигать уровня пены. Наша
модель позволяет рассчитать конвекцию и в такой постановке при условии,
что степень пузырения мала и действуют принятые приближения. В
качестве альтернативного рассмотрен другой сценарий постепенной потери
летучих
в
двухслойной
гравитационно
стабильной
системе.
Предполагается, что локально верхний слой насыщается водой на границе с
нижним слоем (например, за счет выделения пузырей из нижнего слоя при
кристаллизации или при колебании общего давления в камере). В
результате конвективно-дегазационной неустойчивости этот погранслои
периодически всплывает сквозь ненасыщенный флюидом расплав и,
достигая кровли, теряет пузыри и питает флюидом высокотемпературные
фумароль!.
45
ОСНОВНЫЕ ЗАЩИЩАЕМЫЕ ПОЛОЖЕНИЯ
1.
Впервые показано, что дегазация магматических расплавов,
протекающая при малом активационном пороге нуклеации, происходит
равномерно в пространстве с распределением растущих пузырей по
размеру,
приближающемуся
к
степенному,
фрактальному.
Экспериментально установлено, что при гомогенной нуклеации и
достаточно большом пересыщении возможен фронтальный, аналогичный
фронту горения, режим дегазации магматических расплавов.
2.
Дегазация магмы при подъеме к поверхности вызывает ее
кристаллизацию, связанную
с ростом ликвидусных
температур
кристаллических фаз. В свою очередь содержание кристаллических фаз
влияет на динамику извержения. Согласно нашим экспериментальным
кинетическим данным по кристаллизации гавайитов за счет дегазации,
установлено, что при извержении Этны образуются кристаллы плагиоклаза
и клинопироксена с максимальными размерами до 0.6-0.8 мм,
ограничивающими сверху экспоненциальные участки в распределениях по
размерам.
3.
Получено полное точное решение задачи о кристаллизации
бинарного расплава эвтектического типа с осаждением кристаллов в
пограничном слое с постоянной эффективной Стоксовой скоростью
осаждения и учетом противотока жидкости. При значении скорости
осаждения кристаллов, соизмеримой со скоростью перемещения фронта
полного затвердевания, возможно появление разрывных решений с
частично поглощающим (догоняющим) кристаллы погранслоем.
4.
Предложена модель кристаллизации погранслоя с оседающими
кристаллами в приближении распределения кристаллов по размерам и в
зависимости от расстояния до фронта полного затвердевания. Получены
решения с использованием разложения распределения по дельта-функциям
на траекториях, что эквивалентно конечному числу кратковременных актов
нуклеации. Показано, что замедление осаждения кристаллов у границы
полного затвердевания приводит к возрастанию их числа и формированию
узкой зоны, обогащенной оседающим минералом в согласии с
предсказаниями М.Я.Френкеля.
5.
Формализм для описания гравитационного оседания частиц в
гетерофазной реакционной среде переменного состава применен к
решению задачи о течении в гидротермальной среде с гетерофазным
(жидкость - газ) флюидом в приближении двухкомпонентной системы
соль-вода. Показано, что теплоперенос ускоряется в несколько раз по
сравнению с чисто водным флюидом при аналогичных перепадах
температурах.
Подобный
механизм
(одномерной)
конвекции
в
46
гетерофазном флюиде может иметь место в корневых частях активных
геотермальных систем в Калифорнии, Италии и Японии.
6.
Проанализирован режим гетерофазной конвекции при
постоянном потоке дисперсной фазы (газа снизу, твердой фазы сверху).
Показано, что такая двухфазная конвекция может быть рассмотрена как
частный вид естественной конвекции подобный тепловой конвекции. Для
начала этой конвекции должно быть достигнуто критическое значение
седиментационного
числа
Рэлея,
полуаналитическими
методами
установлена зависимость седиментационного числа Нусельта от
седиментационного числа Рэлея.
7.
Метод маркеров обобщен для решения проблемы адвективного
переноса с реакцией в гетерофазной жидкости. На его основе разработан
численный код для моделирования гетерофазной конвекции, сопряженной с
дегазацией и кристаллизации. Проанализированы различные сценарии
развития конвекции в результате пузырения и кристаллизации в
двухслойной магматической камере под вулканом с нижним слоем,
содержащим больше, чем верхний, растворенных летучих компонентов.
ОСНОВНЫЕ РАБОТЫ, ОПУБЛИКОВАННЫЕ ПО ТЕМЕ
ДИССЕРТАЦИИ
1.
Жариков В.А., Симакин А.Г. и Эпельбаум М.Б. (1991).
Моделирование возможности возникновения гранитоидных магм при
взаимодействии базальтовых расплавов с веществом коры// Вестник МГУ,
серия 4. Геология. №2. С. 1-15.
2.
Симакин А.Г. и Трубицын В.П. (1993). Структура верхнего
пограничного слоя кристаллизующегося интрузива при седиментации
кристаллов// Физика Земли. №4. С.20-29.
3.
Жариков В.А., Эпельбаум М.Б., Боголепов М.В. и Симакин
А.Г.(1994) Процессы гранитообразования// Экспериментальные проблемы
геологии/ Ред. Жариков В.А. М.: Наука. С.83-104.
4.
Симакин А.Г. и Чевычелов В.Ю.
(1994). Кинетика
кристаллизации гранитного расплава: приложение к моделированию
образования гранитов фации эндоконтакта// Экспериментальные проблемы
геологии/ Ред. Жариков В.А. М.: Наука. С.121-141.
5.
Симакин А.Г. и Чевычелов В.Ю. (1995). Эксперименталь-ное
изучение кристаллизации Fsp из гранитного расплава с различным
содержанием воды// Геохимия. №2. С. 216-237.
6.
Симакин А.Г. и Трубицын В.П. (1995). Эволюция структуры
остывающей магматической камеры// Физика Земли. №2. С. 40-52.
47
7.
Сгтакин AT. и Эпельбаум М.Б. (1995). Поведение летучих'
компонентов при плавлении водосодержащих пород на контакте с
базитовыми расплавами// Петрология. №3. С. 349-358.
8.
Симакин А.Г., Трубицын В.П. (1996). Струйный режим
осаяадения кристаллов и конвекция расплава в магматических камерах//
Физика Земли. №7. С. 1-7.
9.
Симакин AT. и Трубицын В.П. (1997). Захват фронтом
затвердевания кристаллов растущих и оседающих в магматической
камере// Физика Земли. №10. С. 3-13.
10. Симакин AT,
Трубицын В.П. и Харыбин Е.В. (1998).
Распределение по размерам глубине для кристаллов, осаждающихся в
застывающей магматической камере// Физика Земли. №8. С. 30-37.
11. Симакин А.Г., Армиенти П. И Салова Т.П. (2000). Сопряженная
дегазация и кристаллизация: экспериментальное изучение при плавном
снижении давления// Геохимия. №6. С. 579-591.
12. Симакин А.Г. и Салова Т.П. (2000). Динамические микро­
режимы пузырения гранитного расплава// Физика Земли (№4) 57-63.
13. Симакин А.Г.
и
Трубицын В.П.
(2000). Эффекты
композиционной конвекции при дегазации стратифицированной магмы//
Физика Земли. №11. С. 57-69.
14. Симакин А.Г. и Салова Т.П. (2001). Экспериментальные данные
по эволюции распределения пузырей по размеру при плавной дегазации
водонасьпценного гранитного расплава// Геохимия. №3. С. 294-304.
15. Симакин А.Г, Салова Т.П (2003). Кристаллизация плагиоклаза
из гавайитового расплава в эксперименте и в вулканическом канале//
Петрология. №12. С. 98-109.
16. Симакин А.Г, Салова Т.П., Армиенти П. (2003). Кинетика
роста клинопироксенз из водосодержащего гавайитового расплава//
Геохимия. №12. С. 1275-1285.
17. SimakinA.G., Trubitsyn V.P. andSchmeling Н. (1994). Structure of
the upper boundary layer of a solidifying intrusion with crystal sedimentation//
Earth and Planet. Sci. Lett. Vol.126. P. 333-349.
18. Simakin A.G., Schmeling H. and Trubitsyn V.P. (1997). Convection
in melts due to sedimentative crystal flux from above// Phys. Earth Planet. Inter.
Vol.102. P. 185-200.
19. Simakin A.G.. Armienti P.. EpeVbaum M.B. (1999) Coupled
degassing and crystallization: experimental study at continuous pressure drop,
with application to volcanic bombs// Bull. Volcanol. Vol.61. P.275-287.
48
20. Simakih A.G. an/frubitsyn V.P. (2000). Diferentiation of twocomponent melt solidifying in a magma chamber// Computational seismology
and Geodynamics. Vol.31. P.31-51.
21. Simakin A.G. and Botchamikov R. (2001). Degassing of stratified
magma by compositional convection// J.Volcanol. and Geothenn. Res. Vol.105.
P. 207-224.
22. Simakin A. and A. Ghassemi A. (2003). Salt loaded heat pipes:
steady-state Operation and related heat and mass transport// Earth and Planet.
Sci. Lett. Vol. 215. P.411-424.
23. Simakin A., Salova, T and P., Armienti.P. (2005) Experimental
study of Hawaiitic melt crystallization at low pressure// Europ.J.Mineral (in
press).
НЕКОТОРЫЕ АБСТРАКТЫ ДОКЛАДОВ
НА Н А У Ч Н Ы Х СОВЕЩАНИЯХ
1. Computer simulation of structure of contact granites: an attempt at
the deficit of nucleation data, V All-Union Symposium "Kinetics and
Dynamics of geochemical process", 23-25 May 1989, Chemogolovka.
2. Experimental investigation of Fsp crystallization from a granitic
melt with application to the structwe of contact rocks, with V.U.Chevychelov
and M.B.Epel'baum; Third International Symposium on Experimental
Mineralogy, Petrology and Geochemistry. Edinburgh, UK, 5-7 April 1990, p.
33.
3. Numerical modelling of the formation of the vein-greisen ore
deposit of the orthoraagmatic genesis, with V.A.Zharikov, M.B.Epelbaum,
G.P.Zaraisky and V.S.Kudyaev; 8 lAGOD Symposium august 12-18 1990,
Ottawa, Canada Program with abstracts, 1990, p. 142.
4. Experimental investigation of granitic melt crystallization, with
V.U.Chevychelov; Abstracts of the X I I (last) All-Union conference on the
Experimental Mineralogy, 24-26 September, 1991, Miass.
5. Simakin A.G. A Simple model for treatment of experimental data
on crystal growth from the melt // Intern. Mmeral. Assoc. 16th General
Meeting, 4-9 September 1994 Pisa, Italy, Abstracts, p.378.
6. Chevychelov V.Ju., Simakin A. Kinetics of Fsp crystal growth from
the granite melt with variable water content. Experimental study. // Intern.
Mineral. Assoc. 16th General Meeting, 4-9 September 1994 Pisa, Italy,
Abstracts, p.378.
49
7. Simakin A.G., Schmelmg H., Trubitsyn V. Purely sedimentative
convection with constant crystal .flux from above // Europ. Union Geosc. 8
Strasburg, France, 9-13 April 1995 Terra abstracts, p.300.
8. Simakin A.G. Magma crystallization with Crystal Settling in Upper
Boundary layer: Mixed approach // V.M.Goldschmidt Copnference, March 31
- April 4, 1996, Heidelberg, Germany; J.of Conference Abstracts 1 (1996),
p.571,
9. Simakin, A.; Botchamikov, R. Effects of compositional convection
at degassing of stratified magma. Geophysical Research Abstracts Volume 2,
2000 EGS 25the General Assembly, Nice.
10. A. Simakin A. and A. Ghassemi Salt loaded heat pipes: steadystate operation and related heat and mass transport. Proceedings twentyseventh Workshop Geothermal Reservoir Engineering January 28-30,2002.
50
'^' ^ 9 5 70
РНБ Русский фонд
2006-4
20730
Документ
Категория
Без категории
Просмотров
0
Размер файла
2 528 Кб
Теги
bd000101575
1/--страниц
Пожаловаться на содержимое документа