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

Моделирование фазовых переходов и массопереноса в магматических камерах

ВАК РФ 25.00.10, Геофизика, геофизические методы поисков полезных ископаемых

Автореферат диссертации по теме "Моделирование фазовых переходов и массопереноса в магматических камерах"

Российская академия наук Институт физики Земли им. О.Ю.Шмидта

На правах рукописи УДК 550.311, 552.11,551.21

Симякйн Александр Геннадьевич

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

Специальность 25.00.10 -геофизика, геофизические методы поисков полезных ископаемых

Автореферат диссертации на соискание ученой степени доктора физико-математических наук

Москва 2005

Работа выполнена в Институте экспериментальной минералогии РАН

Официальные оппоненты:

Доктор геолого-минералогических наук Шмулович К.И. (ТОМ РАН) Доктор физико-математических наук Каракин А.В. (ВНИИ Геосистем) Доктор физико-математических наук Геншафт Ю.С. (ИФЗ РАН)

Ведущая организация: Институт геохимии и аналитической химии им. В.И.Вернадского РАН

Защита состоится "23" ноября в 11 часов на заседании диссертационного совета Д.002.001.01 в Институте физики Земли им. О.Ю.Шмидга РАН (123810 Москва, БГрузинская, 10).

С диссертацией можно знакомиться в библиотеке Института физики Земли

РАН.

Автореферат разослан " " октября 2005 г.

Ученый секретарь Диссерта кандидат физ.-мат. наук

АЛ.Трубицын

Ш6-Ч

20730

ВВЕДЕНИЕ

Необходимость и актуальность исследований. Наша живая планета находится в постоянном движении. Мантийной конвекцией сдвигаются континенты, продолжаются процессы плавления и перемещения расплавленной тверди, периодически магма выплескивается наружу, выступая наглядным доказательством подземного жара. Логика познания глобальных геодинамических процессов, закономерностей эволюции магматических камер, будь то субвулканические магматические очаги, лавовые озера, палеоинтрузии, несущие руду, или линзы частично расплавленных пород, подстилающих активные современные геотермальные системы - ведет к рассмотрению конвекции в гетерофазных средах с фазовым переходом. Особенное практическое значение в последние время признается за необходимостью предсказания причины, даты и места извержения супервулканов. В исторически обозримом прошлом (75 тыс. лет назад) извержение такого супервулкана в Индонезии, которое превосходит на порядки зарегистрированные современные извержения, послужило причиной многолетней зимы и вымирания более 90% жившего тогда человека разумного. Есть кандидаты на подобные извержения и в настоящее время. Перед учеными, изучающими геотермальные системы, также есть возможность внести свой вклад в решение энергетической проблемы будущего и преодолеть углеводородный голод, толкающий целые нации на путь разбоя и войн.

Тепловая и композиционная конвекция в магматических и связанных с ними гидротермальных системах, как правило, сопряжена с фазовыми переходами, с гетерогенизацией растворов и расплавов и процессами гравитационной сепарации фаз. Эти явления в последнее время стали очень модным объектом для исследований. Если предыдущие декады можно было бы назвать эпохой «двойной диффузионной» конвекции, то в настоящее время большое внимание придается пузыренгао, кристаллизации природных магм, процессам гетерогенизации солевого флюида. В нашей стране застывание интрузий первоначально рассматривалось B.C. Голубевым и В.Н. Шараповым по аналогии с направленной кристаллизацией металлических слитков с диффузионным перераспределением примесей без учета конвекции и осаждения кристаллов. На полуколичественном уровне моделирование становления магматических интрузий с учетом гравитационной сегрегации кристаллических фаз было предпринято группой МЯ.Френкеля в 70-90ые годы 20 века. Но конвективные процессы при этом рассматривались не на строго физическом уровне. В настоящее время сопряженные процессы пузырения и кристаллизации в ходе течения магмы от очага к поверхности рассматриваются О. Мельником и А. Барминым. При теоретическом описании классическая термическая кондрцгпмя и готярофячныр течения

РОС НАЦИОНАЛЬНА. БИБЛИОТЕКА СПетсрвург

о» ЩУ

противопоставляются. Последние рассматриваются как одна из форм течений Рэлея-Тейлора в многослойной системе, состоящей из чистых жидкостей, расположенных в пространстве в последовательности, ведущей к гравитационной неустойчивости (легкая жидкость расположена под тяжелой). Тот факт, что в гетерофазной системе гравитационная неустойчивость достигается распределением переменной доли диспергированной фазы (газообразной, твердой) в объеме жидкости, игнорируется. А диспергированная фаза может перемещаться в поле силы тяжести, что придает этому фактору плотности сходство с переносом тепла.

Основная и ель настоящего исследования состояла в изучении динамики процессов фазового перехода в расплавах на локальном уровне, в погранслое течения, процессов конвекции в гетерогенных средах с перемещающимися фазами и приложении теоретических результатов к геотермальному и вулканическому процессам. В соответствии с поставленной целью решались следующие конкретные задачи:

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

2. Изучение кристаллизации, сопряженной с декомпрессией и дегазацией магматического расплава на примере гавайитовых магм характерных для вулкана Этна, крупнейшего действующего вулкана Европы. Содержательный анализ механизмов формирования экспоненциального распределения по размерам кристаллов магматических минералов.

3. Проведение исчерпывающего анализа решений задачи о пограничном слое застывающего пластообразного магматического тела в равновесном приближении бинарного расплава с диаграммой плавкости эвтектического типа.

4. Приложение модели межфазного гравитационного течения в среде с фазовыми переходами к расчету тепло и массопереноса в гидротермальной двухфазной флюидной системе соль-вода.

5. Количественный анализ конвекции в гетерогенной жидкости при заданном потоке дисперсной фазы на границах. Нахождение критических условий возникновения конвекции в терминах седиментационного числа Рэлея по аналогии с термической конвекцией.

6. Разработка алгоритма и написание программы для численного моделирования процессов гетерогенной конвекции в магматической системе. Проверка работы программы на тестовых примерах.

7. Использование программы для расчета процесса конвекции, вызванной пузырением и сопряженной кристаллизацией в магматических камерах под вулканами.

Научная новизна.

1. Впервые проанализированы различные микрорежимы пузырения магматического расплава. Равномерная гетерогенная нуклеация при пониженном пороге активационного барьера фазового перехода приводит к заполнению пространства пузырями разных генераций с приближением распределения по размеру к фрактальному, с фрактальной размерностью 2.5-2.8. Режим фронтального пузырения возможен из-за существенных вязких напряжений у растущих пузырей малых, близких к микронным, размеров.

2. Впервые комплексно проанализирована диаграмма плавкости и кинетика роста главных кристаллических фаз гавайитов при низких давлениях водяного флюида. Рассмотрение скоростей роста и особенностей состава кристаллов плагиоклаза и клинопироксена по микропримесям дало возможность связать экспоненциальный участок функции распределения этих минералов в гавайитовых лавах вулкана Этна с процессом дегазации и сопряженной кристаллизации лавы в жерле вулкана Этна.

3. Впервые получено точное решение задачи Стефана, сопряженной с объемной кристаллизацией оседающей фазы в двухкомпонентной двухфазной системе. В зависимости от скорости оседания кристаллов и продвижения границы полного застывания возможны различные варианты распределения оседающих кристаллов, включая разрывные близкие по природе к ударным волнам, о существовании которых было впервые заявлено М.Я.Френкелем.

4. Впервые математически описано стационарное одномерное течение в двухфазном флюиде на примере системы КаС1-Н20. Оценены потоки тепла и массы через колонну флюида, ограниченную снизу границей кипения солевого флюида, а сверху конденсации газообразной фазы.

5. Впервые получены решения задачи о кристаллизации в пограничном слое с оседающими кристаллами на уровне распределений кристаллов по размерам.

6. Впервые показано, что помимо кратковременного течения типа Рэлея-Тейлора в гетерофазных жидкостях с относительным перемещением фаз возможна перманентная конвекция при постоянном потоке оседающей (всплывающей) фазы. Найдено минимальное значение седиментационного числа Рэлея, при которым возможна стационарная конвекция, и которое оказалось близко к Определенному ранее В,П.Трубицынам и Е.Харыбиным (1989) методом анализа устойчивости бескнечно малых возмущений.

7. Разработан метод решения задачи моделирования конвекции при существенно адвективном переносе диспергированной фазы с внутренним источником. Впервые в геофизических приложениях маркерная методика применена для решения на характеристиках задач о перемещении границы двухфазных областей, ограниченных фронтом растворения.

8. Впервые получены численные решения, описывающие процессы конвекции, сопряженной с пузырением и кристаллизацией в магматических очагах под вулканами.

Практическое значение. Экспериментальные данные и модели дегазации магматических расплавов могут быть практически применены в вулканологии. Они дают представление о том, насколько неравновесным может быть процесс дегазации, и как кинетика пузырения может привести к катастрофическим извержениям, позволяют интерпретировать динамику современных и палеоизвержений по распределению по размерам пузырей и кристаллов в лавах. Модель одномерного гравитационного течения в колонне гетерофазного реагирующего флюида была применена для оценки скорости тепло- и массопереноса в корневых частях геотермальных систем. Теория седиментационной конвекции применима ко многим конвективным системам как технического, так и геологического значения. Программу для расчета конвекции, вызванной дегазацией и сопряженной кристаллизацией магмы, можно использовать для расчета флюидного потока на некоторых вулканах андезитового состава Курильской и Японской островной дуг. Как показано автором, этот код также можно использовать для численного моделирования солевых экструзий в условиях перманентного поверхностного растворения.

Исходные материалы и методы исследований. Экспериментальная часть работы базируется на более 150 экспериментах с гранитным и гавайитовым водосодержащим расплавом. Экспериментальные образцы охарактеризованы сотнями снимков в отраженных электронах, десятками ИК спектров, многими сотнями анализов на рентгеновском микроанализаторе. Расчетные результаты получены только с использованием программ автора.

Построение диссертации

Работа состоит из введения, 8 глав и заключения, изложенных на 183 стр., включает 75 рис., 6 табл., список литературы из 226 наименований.

Аппробаиия работы

Результаты работы докладывались на международных конференциях в Москве (1991), Пизе (1994), в Страсбурге (1995), в Гейдельберге (1996), в Ницце (2000), в Стэнфорде (2002). На всесоюзных совещаниях в Черноголовке (1998), в Москве (1999-2004), в Миассе (1991). Материалы работы использовались в лекциях для студентов университетов г. Пиза (2002), Флоренция (2002), Мадрид (1999).

По теме диссертации опубликовано более 40 работ с тезисами, включая 6 полноценных публикаций в рецензируемых международных журналах (Physics of the Earth and Planet Inter., Earth and Planet. Sei. Lett., J. Volcanol., Geotherm Res., Bull.Volcan., Tectonophysics).

Исследования в России были поддержаны неоднократно грантами РФФИ, в Германии - DFG, в Швеции (Упсальский университет) - KVA, в США (Университет Сев. Дакоты) - грантами DOE и EPSCOR.

Начальный интерес к динамике геологических процессов сформировался в студенческие годы во время занятия поисками рудных месторождений методами газовой съемки в бескрайних степях Казахстана. Первым конвективным процессом, с которым я столкнулся благодаря моим руководителям к.г.-м.н. С.А.Воробьеву и проф. А.П.Соловову, было перемешивание газа в приповерхностных пористых породах стохастическими колебаниями атмосферного давления. Решение в случайных функциях для распределения давления и скорости фильтрации с глубиной до сих пор является одной из самых полезных и интересных моих работ. Переход аспирантом в Институт экспериментальной минералогии АН СССР после недолгой работы практического геохимика приблизил меня к расплавам, диаграммам плавкости, к динамике фазовых переходов. На начальном этапе знакомства с веществом при высоких температурах и давлениях большую пользу принесло общение с профессионалами из лаборатории магматизма: с А.Кузнецовым, А.Чехмиром и М.Б.Эпельбаумом. Позднее КИ.Шмулович приобщил меня к вулканам и дегазации магмы. Основным же соавтором в экспериментальной части работы, которая приложила руки к тому, чтобы пузыри и гавайиты стали так близки и понятны, является Т. Салова. А все основные теоретические результаты диссертации были получены благодаря сотрудничеству с членом-корреспондентом РАН В.П.Трубицыным. Знакомство с ним стало решающим шагом, без его высоких требований и превосходно развитого физического чутья мне не удалось бы добиться успеха. За помощь при обсуждении научных результатов автор благодарен академику РАН В.А.Жарикову, долгое время руководившему Институтом экспериментальной минералогии РАН и внесшему огромный вклад в науку. Нельзя не отдать должное заграничным коллегам. Ни разу не побывав на вулкане Этна, я благодаря сотрудничеству с проф. Pietro Armienti стал одним из основных его исследователей методами экспериментальной петрологии. Во время многочисленных визитов в достопочтенный Пизанский университет, в котором прославился Галилей и работал Э. Ферми, довелось наблюдать гавайиты на экране электронного микроскопа и в окуляры FTIR микроанализатора с самим Pietro. Очень полезны были беседы с его коллегами вулканологами M.Rossi, M.Pollacci, P.Papale и их мудрым шефом Prof. Fabricio Innocenti. Изначально туманные мысли автора по поводу сидементационной конвекции были основательно упорядочены в ходе совместной работы с Prof. Harro Schmeling - сначала в Байроте, а затем во Франкфурте-на-Майне. Он же привил мне вкус к использованию маркеров для численного моделирования адвективного переноса. Мир гидротерм открылся мне в Грэнд Форксе, Северная Дакота в общении с Ahmad Ghassemi. Оказалось, что аппарат, который мы с В.П.Трубицыным использовали для расчета отчасти мистического пограничного слоя застывающей магматической камеры, с небольшими модификациями годится для решения задачи об одномерной гетерофазной конвекции рассолов в корневых частях геотермальных полей.

ОСНОВНОЕ СОДЕРЖАНИЕ РАБОТЫ.

Работа состоит из экспериментальной и теоретической частей. В экспериментальной части (глава I и П) представлены результаты оригинальных исследований процессов йузырения и кристаллизации магматического расплава. Получены важные данные по нуклеации и формированию распределений по размерам пузырей в кислых расплавах. Кристаллизация основных фаз гавайитового расплава изучена при постоянном переохлаждении и при сбросе давления и дегазации при верхнекоровых давлениях. Результаты использованы для интерпретации условий кристаллизации гавайитовой магмы вулкана Этна (Сицилия). Теоретическая часть посвящена решению системы дифференциальных уравнений, описывающих охлаждение, рост и оседание кристаллов в пограничном слое, в равновесном и кинетическом приближениях (главы III-V). Особо рассмотрена конвекция в гетерофазной магме, инициированная потоком кристаллов или пузырей (главы У1-УП). Методом граничных интегралов решена задача о стационарной конвекции, найдены критические условия для ее начала. Модели гетерофазной конвекции применены к описанию корневых частей геотермальных систем и конвекции, вызванной пузырением магм в камерах под вулканами (глава VIII).

I. Физика пузырения магмы.

Фазовое превращение начинается с зарождения новой фазы в среде неравновесной, старой фазы. Зарождение и последующий рост с конечными скоростями в условиях закрытой системы, т.е. без больших относительных перемещений фаз, может отвечать условиям быстро протекающих дегазации и кристаллизации магмы в субвулканических и вулканических условиях. Динамику этих процессов можно изучать экспериментально в лабораторном масштабе времени.

Гомогенное зарождение. Основой теории зарождения является представление о важном вкладе поверхностных сил, межфазовой энергии в энергию вещества в субмикронных объемах. Повышенная энергия вещества с большой удельной поверхностью способна компенсировать ту разницу химических потенциалов, которая является движущей силой фазового превращения. При некотором критическом размере устанавливается неустойчивое равновесие. Для критического пузыря водного флюида справедливо соотношение между приращениями химических потенциалов растворенной и выделившейся воды при сбросе-давления от равновесного Р0 до текущего Р1:

АМ"т20=АМ"/°, Р0=*Р, (1.1)-

При сбросе давления на несколько сотен бар парциальный мольный объем растворенной в расплаве воды меняется незначительно, а флюида

значительно. Проинтегрировав, соответствующие члены по давлению и учтя вклад поверхностной энергии пропорциональный кривизне пузыря, получим

А*. -Л). V,(2.1)

Откуда нетрудно выразить размер критического зародыша через характеристический капиллярный размер Ro и относительное падение давления Pi/P0

о __К_ ,= v5£_ и -1Z. г--5. i"t i^

(1 + S(l - «)(* -1))""-°" _ х' b Рс' Р. К }

При размере кластера, меньшем критического, его рост

H20 ^ Н20 т/термодинамически невыгоден, поскольку ц fl>ji m. Когда критическии

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

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

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

флюидного компонента в расплаве, его диффузионную подвижность и

размер критического зародыша.

Анализируя выражение для гомогенной скорости нуклеации нетрудно получить оценку времени задержки при этом процессе. Условно ее можно принять равной времени, за которое образуется, хотя бы один пузырь в характеристическом объеме расплава в 1 см3 за обозримое время от 1 сек до 1 недели. Ввиду экспоненциального характера зависимости скорости роста пузырей от пересыщения, пороговая величина сброса давления (или эквивалентного пересыщения), необходимая для начала гомогенной нуклеации, слабо зависит от принятой величины максимального времени ожидания. Наши оценки минимального понижения давления близки к экспериментальным (Mangan and Sisson, 2000). Критический порог пересыщения свойственен гомогенной нуклеации не только при пузырении, но и при других фазовых переходах. Так при росте кристалла существует первое критическое переохлаждение (от долей градуса до 10° в зависимости от системы), выше которого становится возможен рост по механизму двумерного зарождения центров роста на грани. При более низких переохлаждениях рост происходит по спиральному механизму (движению ступени выхода спиральной дислокации по поверхности грани), являющемуся по существу механизмом самораспространяющегося гетерогенного зарождения (Sunagawa, 1985). Аналогичное критическое переохлаждение существует и для гомогенного объемного зарождения кристаллов в расплаве (Donaldson, 1981). Существование порога гомогенного зарождения может проявляться в нелинейных эффектах при кристаллизации расплавов (Hort, 1995).

Дополнительной особенностью гомогенного зарождения пузырей флюидной фазы из расплава является минимальное начальное содержание воды, при котором возможно гомогенное зарождение. Анализ уравнения

мае. ^^

скорости гомогенного зарождения показывает, что при начальном содержании воды в кислом расплаве 2-3 мас.% (оценка зависит от экспериментальной точности определения многих параметров) или давлении насыщения 250-500 бар сброс давления до атмосферного не обеспечивает условия начала гомогенной нуклеации. Для основных маловязких расплавов критическое начальное содержание для осуществления гомогенного зарождения пузырей ниже и составляет около 0.5мас.%. Если отсутствовали условия для гетерогенной нуклеации, магматическая вода может сохраняться в закаленных вулканических стеклах.

Роль вязких напряжений. В последнее время в вулканологии большое значение стало придаваться явлению

фронтального разрушения

закристаллизованной на 80-90% магмы экструзий с частичным расплавом, насыщенным

летучими компонентами под давлением на несколько сотен бар выше атмосферного (А1сШжоу, 1996).

Экспериментально было

установлено, что при внезапном контакте пористой среды с атмосферой (в результате обрушения конуса) дегазация и хрупкие разрушения композита начинаются с поверхности и сопровождаются послойным отстрелом обломков. В больших масштабах это может привести к катастрофическим эксплозивным извержениям типа Сан-Хелен. Аналогичное вышеописанному процессу фронтальное пузырение от поверхности расплава наблюдалось нами в эксперименте (см. рис.1). Фронт пузырения распространялся на некоторую глубину от поверхности и останавливался. В экспериментальных образцах наблюдаются повторяющиеся полосы мелких субмикронных пузырей, параллельные внешней границе. Скорость распространения фронта была достаточно велика, поэтому нам удалось зафиксировать только продукты этого процесса, а не его стадии. Известно, что скорость роста пузырей на ранней стадии лимитируются вязкостью расплава, а на более поздней - диффузионным подтоком воды из окружающей среды (Ва§с1а85агоу « а1., 1996). На «вязкой» стадии роста можно принять, что давление флюида внутри пузыря отвечает его содержанию в расплаве, т.е. приблизительно равно начальному давлению Р0, при котором расплав был насыщен, а давление в расплаве на удалении равно внешнему Р1< Ро. Рост пузырей описывается в приближении разбиения объема смеси на сферические ячейки с пузырем внутри. Нами

ш ш

Рис.1. Фронт пузырения, распространявшийся от внешней границы, высота микрофото 80 мкм.

показано, что при этих допущениях давление в расплаве Рт, равное симметричной части тензора вязких напряжений, не зависит от радиальной координаты (расстояния от поверхности пузыря) и составляет

= 2,3 ДРжР,-Р, (4.1)

'о "о г\ г

где г - радиус пузыря, ij - внешний радиус оболочки, a (rin-r") - величина пропорциональная количеству расплава, приходящемуся на один пузырь. Значение п=2 отвечает плоскому решению, а п=3 - сферически симметричному. Из формулы (4.1) вытекает, что при некотором относительном размере пузыря давление становится нулевым и даже формально отрицательным при дальнейшем его увеличении. В этом случае могут быть достигнуты условия гомогенной нуклеации в условиях растяжения при микро-разрывах расплава даже вопреки классической теории нуклеации. Вязкость кислого расплава при содержаниях воды в 41.5 мас.% и при ликвидусной температуре составляет 106-108 Па сек, поэтому он проявляет при больших скоростях деформации вязко-упругие свойства и может трескаться. Вследствие механических эффектов пузырение может распространяться как фронт горения, когда множество растущих мелких пузырей создает условие для зарождения следующих. Приближение к свободно деформируемой поверхности облегчает условия возникновения зоны пониженного давления. На внешней границе фронта пузырения иногда наблюдаются деформированные пузыри с трещинами растяжения на концах. Вблизи поверхности слои пузырей деформируются в пологие складки малой амплитуды. Описанное явление пока никак не исследовано другими авторами, которые лишь констатируют влияние вязких деформаций на форму маленьких пузырей в окрестностях больших, начавших рост раньше (Larsen and Gardner, 2000).

Наличие гетерогенных центров зарождения в расплаве значительно облегчает дегазацию, делает процесс пузырения более регулярным и однородным по объему. В качестве центров выделения флюидной фазы могут выступать поверхности раздела кристалл/расплав. Микро-пузыри менее растворимых, чем вода, флюидных компонентов, таких как благородные газы, азот, углекислота, также могут выступать центрами роста водных пузырей. Причем в силу зависимости химического потенциала воды от ее доли во флюиде, микропузыри также должны преодолеть барьер по размеру, прежде чем дальнейший рост станет термодинамически выгодным процессом. Но вероятность преодоления барьера гораздо выше, чем для гомогенного зарождения, поэтому становится возможной перманентная нуклеация или несколько актов этого процесса при умеренных пересыщениях.

Измерение распределений пузырей по размеру. При декомпрессии расплава пузыри перманентно зарождаются и растут, что приводит к формированию распределений пузырей по размерам. Под

(дифференциальным) распределением n(R) (размерность ед./см4)

понимается функция, связанная с интегральным распределением по

размерам Ф(Я), равным числу частиц размера не больше R в единичном

объеме магмы: я

jn(r)dr = Ф(Я) (5.1)

о

Распределений пузырей по размерам, рассчитанные Торамару (Toramaru, 1995) по классической модели гомогенной зарождения, имеют ярко выраженный экстремум, отвечающий одному акту гомогенного

зарождения. Распределения по размерам пузырей в продуктах быстрой дегазации магмы в природе (Mangan and Cashman, 1996; Simakin et al., 1999) описывается экпоненциальной или более сложной зависимостью, практически монотонно убывающей (за исключением локальных колебаний) с ростом размера пузырей. Наблюдаемая разница в распределениях свидетельствует о том, что природные магмы гетерофазны и нуклеация пузырей носит гетерогенный характер.

К началу наших исследований (Симакин и Салова, 1998) экспериментальных работ по изучению динамики пузырения на уровне распределений по размерам практически не было. И в более поздних работах (Gardner et al., 1999) приводятся гистограммы числа пузырей, установленных в двухмерных сечениях с изображений шлифов. Между тем истинная функция распределения находится при пересчете на объем (стреокоррекция), а без этой операции кажущаяся плотность крупных пузырей' завышаётся, и гистограммы с экстремумом отличаются от монотонных распределений (пример распределений на рис.2). Нами экспериментально изучено пузырение гранитного расплава при сбросе давления с постоянной скоростью примерно за сутки с 1000 до 500 атм при Т=750°С. Дегазация протекала в режиме гетерогенной нуклеации на микропузырях, наполненных преимущественно азотом, захваченных расплавом при наплавлении. Состав газа в микропузырях оценен путем наблюдений на криостолике, которые показали, что температура кристаллизации заполняющего их газа составляет -131°С. Измерены распределения по размерам пузырей флюидной фазы в стеклах, полученных при последовательной закалке серии образцов в ходе нескольких экспериментов. Установлено, что при дегазации объем выделившейся фазы отстает от равновесной величины. В зависимости от содержания и размера воздушных пузырей в исходном

Рис.2. Экспериментальные

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

стекле объем выделившейся фазы при 500 атм. либо приближается к равновесному (27 об.%), либо остается вдвое меньше (13 об.%). Сигмаобразная форма зависимости объема

выделенного флюида от времени (см. Рис. 3) с максимальным отклонением от равновесия в промежуточной стадии дегазации отмечается также и в работе (Gardner et al., 1999). Интегральное объемное содержание числа пузырей (размером свыше 5 мкм) в ходе дегазации существенно не меняется, оставаясь на уровне 3-6 106 ед-на мм3 расплава. Некоторое возрастание объемного содержания числа пузырей на единицу объема стекла можно интерпретировать как гетерогенную нуклеацию со скоростью (1-3)10Э см"3 сек"1. Конечные распределения могут отвечать по форме экспоненциальному, степенному (фрактальному) или унимодальному с максимумом в зависимости от объема и распределения микропузырей воздуха в исходном стекле (расплаве).

Фрактальное распределение в генеральной совокупности. Несмотря на наблюдаемое разнообразие и эволюцию распределений по размерам в индивидуальных опытах, совокупность распределений из различных экспериментов и природных данных укладывается в одну генеральную совокупность. Для сопоставления были использованы распределения по размерам пузырей (bubble size distribution, BSD) природных образцов, выбранных из продуктов эксплозивной активности Mt. Etna и Vulcano (Aeolian Islands, Italy).

Образцы были отобраны с точки зрения представительности процесса дегазации при быстрой декомпрессии магмы в магматических фонтанах или других формах вулканической активности. В образцах гавайитов Этны и латитов с Вулкано процессы течения не оказывали серьезного влияния на распределение и форму пузырей, которая была близка к сферической. Все измеренные BSD имеют явно выгнутую кверху форму с выраженным ростом плотности распределения пузырей размером менее 0.4-0.06 мм в зависимости от образца. Для сравнения на общую диаграмму в билагорифмических координатах (lg(n), ln(r)) нанесена совокупность экспериментальных BSD (из работы Simakin et al., 1999).

Рис.3. Объемное содержание пузырей в образцах в зависимости от давления, рассчитанное по распределению пузырей поразмерам. Сплошная линия объемная доля

I у -2 * x - 1 2 I OA*

„ ..................VU*»« ри/пм с SE в* Efr>*'969 *гл1

• See* ■ f это 'ЭО4,-«? e «CSG Lb»1.«ee*-jOt

N, < Swv. ? sv-sa «00592 * 8? fell* 10Sb

• x' eesa env »ew »mw # et«

^ > er2 SV Etrw 11Мв*>14*{ .

О i!> л ftfrti)

Рис. 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) мы обнаружили у BSD воздушных пузырей в исходном стекле, вероятно сформировавшихся при коалесценцияи микропор при наплавлении из порошка.

II. Кристаллизация гавайитового расплава сопряженная с дегазацией

Декомпрессия магматического расплава, во многих случаях адиабатическая, без значительных теплопотерь, сопровождается выделением летучих компонентов после достижения уровня насыщения. В свою очередь потеря летучих компонентов, как правило, ведет к кристаллизации расплава. Это связано с тем, что вхождение летучих компонентов, особенно воды, ведет к понижению температуры равновесия кристалл - расплав (ликвидусной температуры) за счет понижения активности соответствующих минералов в расплаве. Особенно сильное влияние вода оказывает на температуру равновесия плагиоклаз - расплав. Выделение воды вызывает кристаллизацию расплава, даже при неизменной или растущей за счет теплоты кристаллизации температуре. Явление кристаллизации, сопряженной с дегазацией было изучено экспериментально на примере расплава гавайитового состава.

Гавайиты Этны. Гавайит является умеренно щелочной вулканической породой. Гавайитовые лавы - характерный продукт вулканизма во внутриплитных условиях и в горячих точках. При извержении самых продуктивных в мире вулканов Килауэа (Гавайи) и Этна (Сицилия) преобладают гавайитовые лавы. Самый большой активный вулкан в Европе Этна расположен в локальной зоне растяжения (тектоническом окне) у границы Апеннинской аккреционной призмы (Doglioni et al., 2001). Его практически непрерывная активность и необычайно высокое содержание летучих компонентов объясняются механизмом генерации магмы с участием флюидного компонента, высвобождаемого из погружающейся в зоне субдукции Ионической плиты (Metrich et al., 1993). Эволюция состава лав Этны от субщелочного к натровому щелочному и, со временем, к щелочным породам калиевой спецификации (Schiano, 2001; Tonarini et al., 2001) также связывается с участием корового материала из погружающейся Ионической плиты.

Проверка надежности расчета диаграммы плавкости по программе MELTS. До начала наших экспериментальных исследований диаграмма плавкости гавайитов изучена очень незначительно, можно назвать лишь две работы (Trigila et al., 1990; Metrich and Rutherford, 1998). Наиболее плодотворным представляются не многочисленные определения условий плавления гавайитов в широком диапазоне вариаций состава, содержаний воды, летучести кислорода и давлений, а проверка достоверности результатов расчетов по программе MELTS (Ghiorzo and Sack, 1995). Исходя из этого, были определены ликвидусные температуры главных фаз гавайитов в условиях насыщения водой, характерных для вулкана Этна, при 600 бар и содержании воды в расплаве 2.4+0.1 мас.% и летучести

кислорода на уровне NNO. Для определения была использована прецизионная затравочная методика, при которой температура равновесия пары минерал/расплав определялась по началу кристаллизации твердой фазы с обрастанием затравки. Найдены температуры парных равновесий оливина (1120±8°С), клинопироксена (1105±5°С) и плагиоклаза (1070±5°С) Для плагиоклаза найденная таким путем температура является температурой метастабильного продолжения ликвидуса, который примерно на 20° выше температуры равновесного появления плагиоклаза при равновесной кристаллизации гавайитового расплава после выделения оливина и клинопироксена. Сравнение экспериментальных определений ликвидусов главных фаз с результатами расчетов с MELTS показало, что ликвидусные температуры клинопироксена и оливина в гавайитах воспроизводятся с точностью ±5", а расчетный ликвидус плагиоклаза занижен, но не более чем на 20-25°. Дополнительно определена температура ликвидуса клинопироксена при содержании воды в 1.5 мас.%, которая также совпала с расчетной. Установлено, что программа MELTS не только хорошо описывает толеитовую систему, но успешно применима к термодинамическому моделированию щелочных базальтов.

Скорости роста основных фаз гавайита. Для выращивания кристаллов из гавайитового расплава применялись затравки, т.е. кристаллы минералов близких к интересующим нас по составу. Таким образом удается обеспечить беспрепятственный рост изучаемых кристаллов до начала массовой, объемной кристаллизации и установить их скорости роста. При Т=1060°С (переохлаждение 15°) и Рнго=600 бар скорость роста плагиоклаза оценена в 1.4-10"5 см/сек, а время задержки гетерогенного зарождения оказалось близко к нулю. При том же содержании воды и Т=1080°С (ДТ=25±5°) лейстовидный клинопироксен растет со скоростью 2-2.510"5 см/сек, а при содержании воды в 1.5 мас.% скорость роста оценена в 3-410" 6 см/сек. Помимо этого установлено, что с ростом переохлаждения происходит увеличение содержания примесей: глинозема в клинопироксене и железа в плагиоклазе. Плагиоклаз затравок, а также новообразованной при дегазации и перекристаллизации фазы содержит железа не более 0.60.7 мас.% (в пересчете на FeO). Лейсты, растущие при 1060°С (ДТ=10-15°), содержат порядка 0.9-1.2 мас.% FeO. С ростом переохлаждения до 20-25°С содержание железа растет до 1.5-1.7 мас.%. Максимальное содержание железа до 3-3.5 мас.% наблюдается в закалочных лейстах. Знание особенностей состава и порядка величины скорости роста кристаллов позволило критически интерпретировать распределения по размерам кристаллов клинопироксена и плагиоклаза из лав вулкана Этна.

Модели распределения кристаллов по размерам кристаллов. Для изучения динамики кристаллизации и гетерофазной конвекции важно знать размер флюидных пузырей и кристаллов, поэтому, строго говоря, моделирование конвекции с одной равновесной оседающей

кристаллической фазой становится трехмерной задачей, где третье измерение - это размер в функции плотности распределения кристаллов по размерам. Самая простая постановка состоит в описании распределения по размерам в однородной системе. Однородность достигается либо отсутствием перемещения фаз (локальная однородность), либо перемешиванием объема жидкости, претерпевающей фазовый переход. Тогда справедливо уравнение сохранения для функции плотности распределения:

д n(r,t) д ... -пг"

,, ' +—(V„n-D--) = а-, (1.2)

dt дг д г г0

где член справа отвечает потере кристаллов из конвектирующего объема, в закрытой системе а=0. Согласно (Martin and Nokes, 1989) скорость оседания кристаллов из турбулентно перемешиваемого объема пропорциональна Стоксовой скорости осаждения, а значит т=2. Граничным условием для функции распределения является ее значение при нулевом радиусе, равное отношению скорости нуклеации Jo к скорости роста Vcr n(0)-J(/Vcr. Наибольшее применение получило упрощенное решение уравнения (1.2), предложенное Marsh (1985). Оно описывает стационарное состояние в пренебрежении диффузионным членом (D=0), а оседание, согласно Marsh, описывается с а=1, т=0. Принято, что скорость кристаллизации не зависит от радиуса кристалла и является постоянной, как и скорость нуклеации. Такая аппроксимация отвечает условиям кристаллизации солей из водного раствора в техническом кристаллизаторе. Нетрудно установить, что решением является экспоненциальная функция

п(г) = по )> "о = Jо^с, (2.2)

vcrxa

Поскольку экспоненциальной функцией часто можно аппроксимировать природные распределения кристаллов по размерам, либо полностью либо частично, формула (2.2) получила широкое распространение (Mangan, 1990; Cashman and Marsh,1988). Временная константа т0, неизбежная из-за соображений размерности в уравнении (1.2), интерпретировалась Маршем (Marsh, 1988) как эффективное время пребывания магмы в конвективно перемешиваемом очаге в условиях постоянного притока из глубинного источника. Задав по каким-то независимым соображениям величину скорости роста, можно оценить to.

Другое тривиальное частное решение описывает кристаллизацию в закрытой системе в пренебрежении дисперсией скорости роста и ее зависимостью от размера кристаллов. Тогда решением уравнения (1.2), удовлетворяющим граничному условию n(0,t)=J(t)/Vcr(t) является

n(r,t)=J(t*)/VJt*), (3.2)

где время t* отвечает началу роста кристаллов радиуса г (г = | Vcr(t)dr ).

Это тривиальное решение, по существу, означает, что распределение кристаллов по размерам рассматривается как запись отношения скорости зарождения к скорости роста за все время кристаллизации данного объема магмы. Как вытекает из анализа природных и экспериментальных данных, для быстрых процессов кристаллизации вулканитов эта аппроксимация более верна, чем модель Марша (уравнение (2.2)).

Распределение кристаллов плагиоклаза в лавах Этны по размеру. Распределения кристаллов плагиоклаза по размерам в гавайитовых лавах вулкана Этна хорошо изучены. Пример распределения по размерам (CSD) кристаллов плагиоклаза из лав вулкана Этна, найденного путем анализа оптических изображений шлифов, и прошедшего стереокоррекцию, приведен на рис.5. Распределение представляет в логнормальных координатах комбинацию практически линейного начального участка (экспоненциального в нормальных координатах) при малых размерах (менее примерно 0.6 мм) и плато при размерах более примерно 1.5 мм. Сочленение происходит по негладкой кривой с несколькими резкими колебаниями. Аналогичные распределения по размерам типичны для многих минералов, например, оливина (Armienti and Pareschi, 1992) из вулканических пород.

С установленными скоростями роста плагиоклаза из гавайитового расплава в условиях близких к природным можно оценить размеры, которых могут достичь кристаллы плагиоклаза за счет кристаллизации, сопряженной с дегазацией магмы. Скорость подъема магмы в канале имеет порядок 1м/сек (Trigila et al., 1990). С учетом того, что плотность вспененной магмы невелика, гидростатическое давление в 600 бар (условие насыщения по Metrich et al., 1993) отвечает глубине 2000-3000 м. Поэтому порядок времени кристаллизации дегазируемой магмы может составлять 2000-3000 сек. Приняв скорость роста плагиоклаза равной 1.5 -5-2 10'5 см/сек, получим максимальную длину кристаллов порядка 2 (1.5 +2 10'5 см/сек 2000^-3000 сек) = 0.6-1.2 мм. Это очень приближенная оценка, но представляется, что по порядку величины она верна, В этом случае примерно экспоненциальный участок распределения кристаллов по размерам отвечает декомпрессионному росту с увеличением скорости нуклеации по мере потери воды. По неопубликованным данным П. Армиенти при изучении кристаллов размером меньше 25 мкм (при сильных оптических увеличениях) и еще более мелких (нанокристаллов) с использованием SEM отмечаются квази- экспоненциальные участки распределения с еще большими плотностями и склонами (dlg(n)/dR). При столь малых размерах заведомо не работает механизм оседания кристаллов, предложенный Маршем.

Размер, мм

Важно также рассмотреть состав

кристаллов малых размеров.

Самые малые, изученные на

рентгеновском микрозонде,

кристаллы плагиоклазов из

гавайитовых лав вулкана Этна

имеют состав Ап57-Апб4 (см. рис.5)

с содержанием железа

(пересчитанного на БеО) не

превышающим 1%. Природные

„ . _ ' кристаллы содержат много

Рис.5. Распределение кристаллов плагиоклаза г

из гавайитов Этны по размерам, МеНЬше железа> чем мелкие неопубликованные данные Р.АгпнепЙ (2002). закалочные лейсты и близки к

синтезированным при

декомпрессии. По этому признаку можно заключить, что переохлаждение (пересыщение) природных магм по плагиоклазу не превышает, в основном, 10-15°. Лишь эпизодически, при кристаллизации отдельных зон, насыщенных расплавными включениями, оно, вероятно, было больше. Таким образом, увеличение числа кристаллов, увеличение плотности распределения к малым размерам, сопровождающееся смещением состава кристаллов к альбиту, можно связать с падением содержания воды. Потенциально, сопоставляя наклоны этих участков, а также сопоставляя составы плагиоклазовых кристаллов, можно оценить скорости дегазации и подъема расплава к поверхности.

Наиболее вероятно, что и экспоненциальные распределения по размерам кристаллов из лавовых озер на Гавайских островах и из лав вулкана Стромболи объясняются, как и на вулкане Этна, зарождением и ростом кристаллов в квази-закрытых для седиментации условиях за относительно короткое время. Наибольшую ценность для оценки скоростей обмена магмы в магматической камере под вулканами представляет плато больших размеров.

III Полное решение задачи о кристаллизации бинарной системы эвтектического типа с оседанием кристаллов.

Первые две главы диссертации посвящены анализу кристаллизации и дегазации магматического расплава преимущественно на экспериментальном материале. Экспериментальные масштабы времени порядка 3 (102-104) сек. отвечают масштабу времени самых быстрых вулканических процессов в природе. В зависимости от состава вязкость магмы варьирует от 102 до 106 пуаз, причем часто при близликвидусных температурах реология становится неньютоновской. Поэтому небольшие по размеру кристаллы и пузыри флюидной фазы, с радиусом менее 1 мм, не успевают за это время переместиться на значимые расстояния, и гетерогенную систему в объеме несколько кубических см можно, в первом приближении, рассматривать как закрытую. При более медленных процессах, в частности при остывании расплава в камерах с линейными размерами от нескольких десятков метров до нескольких километров перемещения фаз становятся важными. За счет комбинации направленной и объемной кристаллизации происходит магматическая дифференциация. Наибольшее развитие этот подход получил в работах М.Я. Френкеля и его группы впоследствии (Арискин и др., 1987; Френкель, 1995). Однако точных решений комбинированной задачи Стефана с направленной и объемной кристаллизацией оседающей фазы до последнего времени не было. На первом этапе мы сохранили допущение, принятое у наших предшественников, о существовании некоторой эффективной скорости оседания кристаллов. Лишь в этом случае удалось получить аналитические

Основные уравнения модели. Нами рассмотрена для простоты двухкомпонентная система (компоненты А и В) эвтектического типа с кристаллическими фазами

постоянного состава А и В. На' фронте полного затвердевания происходит эвтектическая

кристаллизация фаз А и В, а

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

Рис.6. Схема пограничного слоя с .

оседающей надэвтектической фазой (А). и надэвтектическая фаза

Слева вверху распределение температуры, А (см. рис.6). Перенос твердой

излом отвечает границе Стефана - полного фазы можно описать, либо с

затвердевания. Слева внизу упрощенная помощью функции

диаграмма плавкости эвтектического типа, распределения по размерам

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

I.

/ А-И.

А+В /

с» с

п(Я,г,1) кристаллов, растущих со скоростью У(2Л), с учетом зависимости скорости осаждения от размера кристаллов и ¡(г,К)

(1.3)

дt дг дЯ

либо относительно ее третьего момента - объемной доли кристаллов £(г, ¡) и средней скорости осаждения V,

(2.3)

5/ & дЯ д( дг р,

где Г - функция источника, определяемая ниже (уравн. (2.2.11)), р=4тг/3 для сферических кристаллов, а усредненные параметры определяются как:

= 4л73 (3.3)

+ (43)

£ и/гаК

Уравнение сохранения для компонента расплава А, который участвует в объемной кристаллизации, имеет вид

±[(1.Е)р1С]+^:[(1-Е)р1Си,] = -г (5.3)

ОТ ог

где Г- функция источника, С - содержание компонента А в расплаве, Ц -скорость расплава. Аналогично записывается уравнение для компонента В

^.[(1-е)р,(1-С)]+-^[(1-е)р,(1-С)и1] = 0 (6.3)

где функция источника Гв=0, а концентрация компонента В равна (1-С). Уравнение теплового баланса запишем в энтальпийной форме

—и,Н,+{\-е )и,Н,) + Х— (7.3)

^ = и,Н, + (1 - * )и,Н,) + Я0

где Я- коэффициент теплопроводности, а энтальпия смеси Н находится суммированием по фазам (кристаллам и расплаву).

Н = Н,е + Н,( 1-е), где

Н, *(с,,Г+ Н„)р,, Н, * (с^ТЧ Я0>, (8.3)

Это уравнение теплового баланса для двухфазной среды может быть переписано как уравнение для температуры (Симакин и Трубицын, 1995; 1997), если пренебречь разницей в теплоемкости и плотности жидкой и кристаллической фаз

д/т /-.о-.- &Т деи, где Кт - коэффициент температуропроводности, Ь - теплота плавления

(ГсгЫ-'тс^+Ь^ (9.3)

Чтобы определить функцию источника генерации кристаллов Г, необходимо еще одно уравнение. При достаточно большой скорости фазового превращения система не отклоняется далеко от равновесия. В равновесном состоянии температура расплава и его состав связаны ликвидусными соотношениями.

В линейном приближении

Сш\-(Х*-Т)1т, (10.3)

где ТтА это температура плавления (индекс m от melting) чистого компонента А, а /я/ ликвидусный склон. Уравнение (10.3) определяет количество твердой фазы, выделяющейся при понижении температуры, и тем самым неявно определяет функцию Г.

Альтернативным является определение Г, использующее скорость Vcr роста кристаллической фазы А и величину ее поверхности Q. Допуская независимость скорости роста от размера, получим для фазы с распределением по размерам п(г).

Г = р, jn(r)Cl(r)Vdr, V = V(AT) (11.3)

где скорость роста зависит от переохлаждения AT равного разности ликвидусной и текущей температуры:

WQ-T/.Q-Qm, (12.3)

Таким образом, для шести неизвестных функций T(Z,R,t), C(Z,R,t), Ui(Z,t), Us(Z,R,t), r(Z,t), s(Z,t) мы получили систему из шести уравнений: трех уравнений сохранения массы (2.3), (5.3), (6.3) для двух компонетов ( СА и Св) расплава и содержания кристаллов е, уравнения переноса тепла (9.3), уравнения фазового равновесия (12.3) и соотношения между абсолютной и относительной скоростями оседания Us (4.3). В связи с гиперболической природой уравнений седиментации возможны разрывы в распределении кристаллов e(z), на что указывал М.Я.Френкель. Подобные явления встречаются в газовой динамике и при инфильтрационном метасоматозе, что одним из первых описано Д.С.Коржинским.

Анализ решения рассматриваемой проблемы проведен в приближении квазистационарного решения при заданном тепловом потоке на границе полного затвердевания. Нестационарные законы сохранения были преобразованы с введением новых координат связанных с движущейся границей затвердевания t'=t, z'=z- jubdt. При этом производная по времени преобразуется согласно выражению:

так что скорость перемещения границы попадает в стационарное решение, которое качественно зависит от соотношения скорости оседания 5 надэвтектической фазы А и скорости перемещения фронта затвердевания иь. После манипуляции с первыми интегралами уравнений сохранения получается важное соотношение для распределения доли кристаллической фазы А в зависимости от расстояния до фронта полного затвердевания

{\-и,{2)Шь)(\-С{2)) (133)

где В - это константа интегрирования, которая подбирается для удовлетворения соответствующего граничного условия для е. При скорости осаждения больше 17ь

е(0) = О, С = Се, В = Се (14.3)

Для поглощающего пограничного слоя (и/иь<1) граничное условие необходимо ставить на внешней 1ранице

= С = С„, Д = С„ + *0(1-у)(1-С0) (15.3)

где здесь и далее Простое алгебраическое соотношение (13.3) для

е(С) позволяет изучить свойства распределений ф), С(г), не находя их точного вида. Для завершения решения необходимо решить уравнение теплового баланса с учетом (13.3). Величина скорости фронта полного затвердевания также является параметром, который отыскивается, исходя из граничных условий. Распределение температуры найдено в неявном виде как 2?=Ф(Т,В), где Ф достаточно громоздкое выражение, приведенное в диссертации. Важно, что распределение температуры представляет собой гладкую функцию, близкую к простейшему экспоненциальному распределению, свойственному классической задаче Стефана без объемной кристаллизации.

Распределение объемной доли твердой фазы может быть разрывным, поскольку оно описывается гиперболическим уравнением адвективного переноса. Оказывается, что большое значение имеет скорость противотока жидкости в слое с оседающими кристаллами. При пренебрежении этим эффектом из (13.3) вытекает, что при равенстве Стоксовой скорости оседания и фронта затвердевания {¡0~^/иь=1) решение имеет особенность. Оно не имеет физического смысла, так как в этом случае содержание кристаллов у фронта (при }0<1) и на бесконечном удалении (при jo>l) стремится к бесконечности при ]а стремящемся к единице. Тогда как при ]0>1 все надэвтектические кристаллы фазы А, образующиеся в погранслое оседают и не захватываются фронтом. Напротив, при }0<1 все кристаллы захватываются фронтом затвердевания.

Учет эффекта противотока жидкости. Если кристаллы образуются в тепловом пограислое, их Стоксова скорость обычно значительно не превосходит скорости фронта затвердевания, и вследствие влияния противотечения и возрастания вязкости можно ожидать, что в части термального погранслояу'</, а в часта}>1.

Минимальный уровень скорости осаждения, при котором появляется интервал с кристаллами, избегающими захвата, достигается при }0=1. В работе показано, что первая критическая скорость осаждения 5/, отвечающая этому условию равна

---(16.3)

рш1С1)(С0-С,+81е)

Аналогичным образом, второй критический уровень Стоксовой скорости, обеспечивающий условие во всем интервале (Се,С0) достигается при у0 =4.791 (при Со-0.6, Се—0.3, общее выражение громоздко, но получается тривиально).

8>=-479Ь° п с , (17-3)

1 р ° е (1-С0), Таким образом, при значениях скорости осаждения 5 с (>0,5;) все кристаллы фазы А захватываются фронтом затвердевания. В этом режиме скорость фронта затвердевания достигает своего максимального значения равного [/¿,=5'/. При увеличении 5 происходит частичное осаждение кристаллов, сопровождающееся падением скорости фронта до уровня 82/4.791. При 5 появляется стационарное решение, описывающее осаждение всей надэвтектической фазы А. Смена режима от полного захвата до полного осаждения происходит при изменении скорости осаждения от уровня 5/ до 52=2.9-5/ (при Ле=2, Со=0.6, С£=0.3). Эффективный размер кристаллов при этом изменяется примерно в 2 раза. При 5>5у возможно множество стационарных состояний с различными распределениями объемной доли кристаллов и различной степенью фракционирования. Теоретический анализ стационарных состояний был дополнен численным решением уравнений модели. Численные решения показали, что при малых и больших скоростях осаждения быстро устанавливается соответствующие стационарные состояния. Тогда как при скоростях осаждения, близких к скорости перемещения границы Стефана, возможны решения с поглощающим погранслоем, который не покидают оседающие кристаллы надэвтектической фазы. Также возможны разрывные решения разных типов. Более реалистичная картина распределения кристаллов получается с учетом их переменных размеров в процессе роста. Она рассмотрена в главе V диссертации.

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

Гидротермальная система с двухфазным (жидкость-газ) флюидом очень близка по физическому механизму и способу формального описания к разобранной выше кристаллизации погранслоя магматической камеры с оседающими кристаллами надэвтектической фазы. Изучение гидротермальной конвекции имеет растущее практическое значение в современном мире. Приповерхностные гидротермы в значительной степени используются для производства энергии в различных странах мира. В связи с этим моделирование конвекции водного флюида на глубинах до 2-3 км и температурах до 300° стало рутинной инженерной практикой (Pruess, К., code TOUGH2; Oldcnhurg, С., and Pruess. К., 1993). Вместе с тем, по мере освоения практически всех доступных приповерхностных гидротерм в Италии (Ландорелло - Тровалли) и Японии началась эксплуатация горизонтов с глубиной 3-4 км (Barelli et al., 2000; Kasai et al., 2000). На этих глубинах обнаружены рассолы различного состава. При кипении рассолов могут возникать условия для существования обширных двухфазных зон, причем флюиды могут быть весьма агрессивны и способны к растворению окружающих пород. Нами теоретически рассмотрены режимы одномерного двухфазного течения в модельной реакционной среде, состоящей из соли. В реальности реакционный флюид реагирует с силикатными породами и способен выщелачивать отдельные окислы металлов, оставляя инертные компоненты (AI2O3 и в меньшей степени SiCh). В природе система флюид-порода может быть далека от равновесия. Тем не менее, полученные оценки скорости тепло- и возможного массопереноса в такой системе представляют интерес для оценок масштабов природных процессов. С помощью модели рассмотрены закономерности функционирования реальных гидротермальных полей на Западе США.

Основные уравнения модели. Течение каждой из фаз двухфазного флюида подчиняется модифицированному уравнению Дарси кк„(г) дР ,

Я,=—-—PA— + /3,g) , (1-4)

jU, dz

где i=l,2 номер фазы, kri(е) - относительный коэффициент проницаемости в смеси в зависимости от объемной доли газовой фазы е, р, плотность i фазы, - поток массы i фазы.

Закон сохранения энергии записывается в энтальпийной форме и включает в себя тепловой эффект кипения

^((HlPlE + HgPli(l- 8))ф + Я, (1 - + iL я, + ) = Л Г, (2.4)

где ф - пористость породы, Н - энтальпия жидкости и газа, Л -коэффициент теплопроводности. Законы сохранения массы выписываются для каждого компонента в каждой фазе. Рассмотрен случай простой

системы соль-вода, явлениями гидролиза и фракционирования пренебрегают. Для растворенного вещества справедлива система уравнений

^(фер, С;)+|-(9(С/) = Г0 +Г, : (3.4)

at аг

1.т-е)ргС,) + ~(ягСг) = -Т1 , (4.4)

где Cg и С/ являются содержаниями растворенного вещества в газообразной и жидкой фазе соответственно. Го - функцией источника, отвечающей процессу растворения и осаждения твердой соли, Л -функцией источника, описывающей перераспределение соли между газовой и жидкой фазами. Аналогичные уравнения справедливы для содержания воды в жидкой и газообразной фазах. Условия равновесия связывают химические потенциалы компонентов в фазах: соли в газе и жидкости: (С, ,Р,Т) = //„,, (Сг, Р, Т)

воды в газе и жидкости: (С,,Р,Т) = р^ (Cg,PJ) (5.4)

соли в твердой фазе и жидкости: ц\ак (Р,Т) = (С, ,Р,Т) Относительные проницаемости фаз задаются в виде к;-(*-*,')■, k;=(i-*-*g°)n> (6.4)

где жидкость становится подвижной при объемной доле жидкой фазы £>si°, а газ при l-£>eg°. Для простоты часто принимается п=1, eg°=0, ei°=0.

Граничные условия задаются не для всех переменных. В трехфазной системе уравнения химического равновесия задают значения Cg, С/, Р как функции температуры. Для температуры, подчиняющейся дифференциальному уравнению второго порядка, необходимо два граничных условия, таких как заданный тепловой поток и температура на нижней границе кипения. Потоки qi, qg линейно зависимы и одна константа определяет подток флюида через границу кипения снизу. Мы рассматриваем закрытую систему, и некоторое количество воды попадает снизу лишь для того, чтобы компенсировать осаждение соли.

Моновариантное равновесие в трехфазной системе. Нами изучено квази-стационарное решение системы уравнений (1.4-5.4). Подобное решение хорошо изучено для чисто водного флюида (Hardee, 1982; McGuinness, 1996). Трехфазная система (соль-жидкость-газ) является моновариантной и поэтому подобна чистому флюиду. Система уравнений химического равновесия для трехфазного равновесия задает составы газообразной и жидкой фаз, а также давление как функции температуры (Cg(T), Ci(T), Р(Т)). Этого достаточно и для того, чтобы выразить плотности фаз р/(Т) и pg(T) от температуры (Bischoff, 1991). Принципиальным отличием солесодержащей системы является то, что помимо теплопереноса в тепловой трубе происходит и массоперенос. Нетрудно отыскать

выражение для функции источника Го, задающей растворение и отложение соли в процессе течения и связанной с пористостью колонны

•п 1 \ ¿<Р &

С.-С.

(7.4)

В силу однозначной связи давления и температуры при трехфазном равновесии, можно исключить давление, используя уравнения теплового баланса (2.4), Дарси (1.4) и с учетом уравнений сохранения (3.4) и (4.4). Из уравнений сохранения и Дарси можно также получить в м « « « " »•• явном виде выражения для потока

Рис.7. Зависимость объемной доли одной из фаз, скажем, и жидкости от температуры при выражение градиента давления в стационарной гетерофазной конвекции в допущении упрощенного

трехфазной системе жидкостъ-пар-соль. выражения для относительных

тепИаШУР0В УК83аНЫ ЗНаЧВШЯ П0Т0К°В коэффициентов проницаемости

жидкой и газообразной фаз (к\=е, В итоге получается выражение, связывающее долю

жидкой фазы и температуру в стационарном состоянии. Следует отметить, что на топологию решения влияет тот факт, что зависимость Р(Т) имеет экстремум при Т=596°С и заканчивается при Т=800°С плавлением чистой соли. Экзотические стационарные состояния е(Т) могут существовать и при Т>596°С, но негативный наклон зависимости Р(Т) требует нагрева и кипения сверху для этого типа решений, что не соответствует геотермальным системам. Мы ниже рассматриваем решения для Т<Ттах. На рис.7 представлены линии зависимости е(Т), рассчитанные численно для разных тепловых потоков через колонну. Этот график демонстрирует важное свойство решения, имеющего две ветви, отвечающие газодоминирующей и растворо-доминирующей системам. Аналогичное свойство имеет и чисто водная система,, но потоки тепла, переносимые флюидом с солевой нагрузкой, больше. При этом также происходит растворение соли опускающимся флюидом и отложение ее на фронте кипения.

Дивариантное равновесие в двухфазной системе. В рамках принятой модели область двухфазного флюида не ограничивается трехфазной кривой, а простирается между ней и критической кривой на РТ плоскости в область больших давлений. В колонне двухфазного флюида в дивариантном равновесии возможны процессы обмена солью между газообразной и жидкой фазами, объемное кипение и конденсация. При этом

равновесие с солью достигается только на фронте кипения. В данном приближении соль может попасть в систему лишь на фронте конденсации. Формально в равновесном приближении это может происходить за счет массообмена с вышележащим рассолом, находящимся в состоянии термальной конвекции. В реальности же возможно, что кислый сконденсированный газообразный флюид при опускании будет выщелачивать вмещающие породы.

В приближении двухфазного равновесия растворения и кристаллизации соли не происходит, поэтому Г0=0. На основании этого легко вычислить стационарный поток соли в колонне

С,-С1 р,рге/(1-е) + р1р, где В - константа, задаваемая потоком соли на границе колонны. Используя это выражение и предположив для простоты, что вязкость фазы пропорциональна ее плотности и =т, легко получить выражение, связывающее долю жидкой фазы с потоком соли в удобном виде квадратного уравнения

В* =---= е()-е) (9.4)

(С,(Р,Т)~ С, (Л Т))(р, {Р ,Т)~ рг (Л Т))

Это квадратное относительно б уравнение при заданных Р и Т, которое имеет решение лишь при В*<1/4. Критическое значение В*= 1/4 отвечает максимальному тепловому и массовому потоку через двухфазную колонну. Одно решение с меньшим значением е отвечает газо-доминирующей системе, а с большим корнем - растворо-доминирующей. На этом пути оценен максимальный тепловой поток через колонну

ЛЯ ^(А-Р.) (10

Ар

При проницаемости 3.10"16т2, Лр=0.8 г/см3; ра = 1г/смЗ и прочих значениях параметров, отвечающих системе, эта простая формула дает максимальное значение теплового потока 4.3 ват/м2, близкое к численной оценке максимального потока в трехфазной системе в 6 ват/м2. Этот поток можно сравнить с потоком тепла при тепловой конвекции (ТгиЬкзуп е1 а1., 1993)

(? = 0.218ДЯ^, Щ и, = (Ц.4)

где Ь - вертикальная протяженность гидротермальной системы, Др(ДТ) -разность плотности, отвечающая перепаду температуры ДТ. Тепловой поток в 0.9 ват/м2 получается при ¿1Г=200°С, р05=2.6 г/см3, ср5= 1 дж/г/о, ¿=500 м, ¿7=10"2 см2/сек, «=З Ю 3, Л=10'15 м2 .

Fluid i,___

on о о

0 L+V Heat pipe

14 > i '»«'»' jjy t'» I

о 0 о

О 0

4 »> »» t «> «

ÜUÜ

Sealed zone

Fluid

Liqad

Приложение к некоторым

геотермальным

системам

-dissolution

.pnecpuaon

Рис.8. Схематичное изображение двухфазной конвекции в корневой зоне гидротермальной системы. Ее положение совпадает с реологическим переходом во вмещающих породах при 350-400°С

Америки. С точки зрения, развиваемой модели, скорость ' растворения вмещающих пород при найденной скорости конвекции достаточно велика. Она отвечает увеличению пористости примерно на 0.020.005 об.% в год. Температура предполагаемой двухфазной конвекции отвечает температуре реологического перехода пород от хрупкой к пластичной (см. рис. 8) . А это значит, что линейная скорость компакции, т.е. - оседания вышележащих пород, ожидается высокой и равной нескольким см в год. Описываемый процесс можно

определить как высокотемпературный карст, который, как и обычный карст карбонатных пород поверхностными водами, проявляется в рельефе. Основные результаты этой главы докладывались автором на ежегодном Стэнфордском семинаре по геотермальным системам в 2002 году. Неожиданным сюрпризом оказалось соседство с докладом J. Moore et а1.(2002), посвященным конвекции с противотоком в двухфазном флюиде в геотермальной системе Karaha-Telaga Bodas (Индонезия) и подтверждающим предложенную схему. Они привели данные по жидким и газовым включениям в кристаллах ангидрита, кальцита и флюорита, обнаруженным в минеральных отложениях в породах, вмещающих геотермальную систему. Самые глубокие образцы ангидрита содержат рассолы с содержанием соли в 31 мас.% эквивалента NaCl. Вся совокупность минералогических данных была интерпретирована как проявление процесса- кипения и конденсации флюида, сопровождающееся погружением кислого агрессивного конденсата с растворением вмещающих пород. Другим явлением, вероятно связанным с конвекцией агрессивного флюида на глубине, является опускание поверхности на многих геотермальных полях. На геотермальном поле Geysers (С.Калифорния), по крайней мере, последние 20 лет, поверхность проседает со скоростью 4.5±0.5 см/год (Mossop et al., 1997). Аналогичная скорость оседания 3.5 см/год зафиксирована на геотермальном поле Coso (Ю.Калифорния) (Wieks et al., 2001). Предполагается, что источником деформации, отвечающим распределению скоростей оседания, является плоская зона контракции на глубине 4.6-3.9 км в Coso и 4.1-3.8 км в случае Geysers, что в обоих случаях отвечает условиям упруго-пластического реологического перехода.

V. Описание кристаллизации бинарного расплава с оседанием кристаллов на уровне распределений кристаллов по размерам

Оперируя средними скоростями оседания кристаллических фаз в задаче о застывании магматической камеры, М. Френкель с соавторами и мы вслед за ними, возможно, чрезмерно упрощаем ситуацию. Размер кристаллов, а вместе с ними и скорость оседания меняется при росте кристаллов. Следовательно, от кинетики роста будет зависеть и эффективность фракционирования, а также поток твердой фазы от верхнего контакта остывающей интрузии, с которым может быть связаны гравитационные течения.

Нами получено несколько решений задачи о кристаллизации бинарного расплава, такой же, как в III главе. Для простоты мы ограничились ростом одной надэвтектической фазы А перед фронтом полного затвердевания. Для распределения кристаллов фазы А, растущих в объеме, справедливо уравнение (1.3). Мы анализируем стационарное решение при скорости продвижения границы затвердевания Ub

Sn(U,(z,R)-Ub) + y8n ü5.

dz SR '

где U, = (1 - е ) • S(R, е), а Стоксова скорость осаждения равна

Ь-Чр.-рЛ*!*'* (2-5)

Объемное содержание надэвтектической фазы А в расплаве равно

R3dR (3.5)

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

Г(г) = 4?гК ^R2n(R,z)dR (4.5)

В стационарной модели приближенно можно положить, что в любой момент времени кристаллы гетерогенно зарождаются только в узком слое в окрестности нижней границы теплового погранслоя, так как внутри погранслоя кристаллы уже возникли раньше, а в основном объеме камеры температура еще несколько выше температуры ликвидуса. Рост кристаллов на микровключениях начинается не одновременно, поэтому в (R,z) пространстве решение для функции плотности n(R,z,t) занимают некоторую полосу. Но для простоты в настоящей работе мы рассмотрим предельный случай, когда при понижении температуры ниже ликвидусной начинается одновременный

рост всех кристаллов. Ему отвечает бесконечная скорость гетерогенной нуклеации, и плотность распределения по размерам и пространству п(Я,г,1) имеет вид 5 - функции

П(Я,2) = 3(Я-Я(г))Щг), (5.5)

где К*(гД) - уравнение траектории одиночного оседающего кристалла, а Щг) - плотность распределения кристаллов по пространству равная первому моменту распределения по размерам

г)ая = - Я (г)) • Щг)сШ = (6.5)

о о

Траектория кристалла параметризуется временем / и может быть представлена парой уравнений для Я и г . В движущейся системе координат уравнение траектории записывается в виде

¿.«,.-0. (7.5,

В стационарном случае они сводятся к одному уравнению

¿¡С _ У(г)

¿2 ~и,-и„

(8.5)

В общем случае траекторий может бьггь несколько с учетом возможности нескольких актов зарождения и возвратного (в движущихся координатах) хода кристаллов, покидающих погранслой.

Подставляя (5.5) в (1.5) и учитывая свойства д(11-К^-функции, получим дифференциальное уравнение для изменения Щг) вдоль траектории осаждения кристаллов (характеристик ур. 1.5)

, (9.5)

сЬ ¿Ь

где сШз/с12 =(<Ш8Ш)(<т */с1г). Таким образом, состояние погранслоя задается распределениями четырех функций С(г), Т(г), N(2), Щг). Первые два являются решением равновесной задачи, вторые два находятся из совместного решения пар дифференциальных уравнений (8.5) и (9.5) для К*(г) и Ы(г).

В принятом приближении нуклеация принимается гетерогенной, а распределение температуры и объемной скорости роста задаются равновесным решением. Вычисляются размер кристаллов (Я(г)) и их интегральное объемное содержание (N(2), в ед./см ). Рассмотрены случаи поглощающего пограничного слоя, когда скорость перемещения фронта иь велика и все кристаллы захватываются фронтом затвердевания. При этом

учитывался эффект противотока жидкости и влияние эффективной вязкости, зависящей содержания кристаллов, на скорость их оседания.

Помимо этого использовано представление функции распределения в виде суммы дельта-функций, отвечающих нескольким идеальным дискретным актам нуклеации. Показано, что существуют решения с частичным захватом кристаллов и их частичным осаждением. При моделировании на уровне распределений по размерам эффект накопления кристаллов, достигающих при своем росте скорости осаждения, равной скорости перемещения фронта полного затвердевания, очень контрастен и локален. В целом результаты наших исследований подтверждают принципиальные выводы МЯ.Френкеля о том, что при дифференциации силлов наиболее значимым является не «дождь» кристаллов от верхнего контакта, а периодическая конвективная дестабилизация верхнего погранслоя из-за накопления кристаллов.

Полученные решения не доведены до уровня непосредственного петрологического приложения. Вместе с тем, подходы для дальнейшего развития очевидны. В связи с ростом интереса к описанию процессов затвердевания и дифференциации магматических интрузий на уровне распределений кристаллов по размерам возможны обобщения на случай нескольких кристаллизующихся фаз и учета фактических данных по диаграммам плавкости и кинетики зарождения и роста.

VI. Стационарная седиментационная конвекция

Как рассмотрено в предыдущих главах, магматический расплав лишь в редких случаях можно рассматривать как чистую жидкость. При снижении давления происходит дегазация и сопряженная кристаллизация магмы, точно так же, как и при ее охлаждении. В неоднородной смеси жидкости, кристаллов и газовых пузырей градиенты плотности могут достигать больших значений и вызвать гетерофазную конвекцию.

В том случае если можно пренебречь относительными перемещениями фаз-включений в жидкости, течения в ней становятся эквивалентны обычным течениям Рэлея-Тейлора. Конфигурация, когда легкий материал помещен под тяжелым, абсолютно неустойчива. Наиболее хорошо исследовано развитие неустойчивости Рэлея-Тейлора в среде с плоскими горизонтальными слоями. Показано, что исходное синусоидальное возмущение границы раздела двух слоев развивается по экспоненте, геометрические и реологические константы слоев определяют длину волны максимально растущей моды.

Линейный анализ устойчивости. Исследования устойчивости в общем виде в слоистой среде с перемещающейся фазой, которая задает стратификацию, не проведено. Трубицын и Харыбин, 1989 исследовали частный случай гравитационной устойчивости горизонтального слоя оседающих частиц. Проанализирована линейная устойчивость мгновенной конфигурации с линейным уменьшением содержания частиц дисперсной фазы книзу. Рассмотрим эти результаты подробнее.

Ниже представлена наиболее простая формулировка модели конвекции, а в главе VI проведен анализ исходных уравнений гетерофазной конвекции. Во-первых, это уравнение Навье-Стокса в приближении Буссинеска. Гетерофазность фигурирует в нем только в гравитационном члене через плотность смеси р=е!р3+(1-£5)р1

-чр+я+ч-^Ьи^и^Ь0 (6Л)

Во-вторых, перемещение диспергированной фазы (твердой для определенности) задается уравнением адвективного переноса

(6.2)

В-третьих, строго говоря, в гетерофазной среде справедлив аналог уравнения несжимаемости с источником, который в приближении малого содержания твердой фазы заменяется обычным уравнением несжимаемости для скорости жидкости: _

рос. НАЦИОНАЛЬНАЯ ' (6-3)

I БИБЛИОТЕКА | ,

| С-Петервург I

> 03 М ш '

МО*

При линейном анализе устойчивости рассматривается гетерофазный слой толщиной Я, расположенный у верхней границы столба чистой жидкости. В качестве начального состояния принимается линейное распределение содержания твердой фазы е, с глубиной г:

е\г^=Щ = е0-е01Н*1, 2 6 ((>,£) (6.4)

Слой оседает со скоростью Б:

£ = £•"(2-®,'). .уе(Л,Я + &], е = е„.уе.[О,Л] Возмущения по глубине концентрации (у) и вертикальной компоненты скорости жидкости (иу) с глубиной задаются в экспоненциальном виде

[у.и,]- [Г(у),Щу)]ехр(Л I)со$(кх) Уравнения развития возмущений бесконечно малой амплитуды имеют вид

(6.5)

¿у

Ла,*2Г+ = О, 1У = с1!/¿у2-кг, (6.6)

где сидементационное число Рэлея определяется через характеристические константы системы как

Ка=(РгР,)^Н! = Уе^Н1 (67)

5//

Дисперсионный анализ зависимости Я(Яа5,к) при соответствующих граничных условиях выполнен полуаналитическими методами Трубицыным и Харыбиным (1987, 1989, 1991). Ими показано, что существует критическое число Рэлея, при котором начинается незатухающая сидементационная конвекция, т.е. Я>0. При непроницаемых и скользких границах слоя

> = 0,Я: и!=д1и,1ду1 =0 ими установлено, что критическое число Рэлея составляет И,*=104.76, а волновое число критической моды равно 3.0. Следует подчеркнуть, что начинающаяся конвекция сосредоточена внутри слоя с градиентбм содержания примеси. Поведение возмущения границы между этим слоем и чистой жидкостью не рассматривалось.

Стационарная седиментаиионная конвекиия. Линейный анализ устойчивости свидетельствует о том, что нет принципиальной разницы между термической и седиментационной конвекцией. Термальная конвекция происходит при сопряжении кондуктивного теплопереноса и гравитационного течения. При условии, что Кат<Яатсп' кондуктивный теплоперенос ведет к диссипации возмущений плотности, и конвекция не

развивается. При идеальной седиментационной конвекции (пренебрегая гидродинамической дисперсией, которая сродни диффузии) возмущения плотности рассеиваются адвективным переносом в частном случае линейного неустойчивого распределения примеси в слое при Па^^ая/"'. Однако, строго никем не показано, что оседание дисперсной фазы способно стабилизировать гравитационную неустойчивость типа Рэлея-Тейлора в двухслойной системе. Тем не менее, и в этом случае существуют критические явления при поддержании постоянного потока твердой фазы у верхней или газообразной у нижней границы.

Для поиска условий существования и формы стационарного седиментационного течения в двумерном приближении использован метод функции тока. Тогда после масштабирования (Я - масштаб длины, £о-масщтаб содержания равный значению на границе, 5 - масштаб скорости) уравнение (6.1) перепишется в виде:

Д (6.8)

дх

При постоянном потоке твердой фазы сверху и некоторых допущениях внутри рассматриваемого прямоугольника размера (ах1) будут области, содержащие кристаллы и чистую жидкость, разделяемые границей уь=/(х). Тогда в правой части уравнения (6.8) появится производная от ступенчатой функции Хевисайда

Ло, 8(у - /(х))п • е, (6.9)

Полуаналитическое решение линейного уравнения (6.9) можно записать с использованием граничного интеграла по границе сред и соответствующей граничным условиям функции Грина 0(хь(0>уь(*)>х>у)> где время выступает как параметр:

Г<Ь.у,0= \ С(Хь (0,уь(1),х,у,) ът(агс1^)<11 (6.10)

Приняв, что вертикальные границы симметричные, а горизонтальные скользкие, получим выражение для функции Грина для области (а х 1)

./Г (ппу), (6.11)

яа (т/а+п / а

а

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

^ = ^ = = (6.12) с!хь и и °хь оуь

Численно получены решения системы (6.11)-(6.12) для различных значений

седиментационного числа Рэлея. Установлено, что стационарное конвективное решение появляется при достаточно медленной седементации, отвечающей числу Рэлея Rase<j >105 при а=0.7. Кристаллы, попадая в верхний погранслой вдоль всей верхней границы, течением захватываются в нисходящий поток, который быстро достигает дна и растекается вдоль него. Кристаллы покидают систему у нижней границы в

Рис.9. Зависимость седиментационного числа Нусельта от седиментационного числа Рэлея. Пунктир - параметры Результате осаждения, а чиста* нестабильного включения чистой жидкость поднимается в восходящем жидкости малого объема. Аспектное потоке. Таким образом, получается, отношение 1. что течение развивается около

включения чистой жидкости в смеси. Если включение изначально было недостаточно большого относительного объема, оно погружалось на дно и весь объем заполняется оседающими кристаллами. В работе определен седиментационный аналог числа Нуссельта, отражающего эффективность конвективного переноса по сравнению с простым оседанием. Показано, что, как и в случае тепловой конвекции, число Нуссельта равно 1 при докритическом уровне седиментационного параметра Рэлея, а после этого порога растет скачком, как отношение полной площади ячейки к площади двухфазной зоны (см рис. 9).

Развитые формы седиментационной конвекции расчитаны нами численно с помощью программы, реализующей алгоритм, разработанный по методу

конечных элементов. Опицание алгоритма g представлено в главе VII, пример умеренно Ш надкритического течения представлен на рис. 1 Ю.

К

ш 1 -

Рис.10. Рассчитанное распределение оседающих кристаллов в конвективной ячейке при Ras™ 500. Восходящий поток жидкости в центре.

Экспериментальные_данные по

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

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

диспергированной фазы не проводилось.

Известные нам эксперименты относятся к моделированию гравитационного течения слоя, обогащенного подвижной дисперсной фазой. Thomas et al. (1993). использовали двухслойную жидкую систему с более вязким верхним слоем. Снизу равномерно по площади задавался поток газовых пузырей. В силу неразрывности потока на границе при замедлении пузырей возрастало их содержание, со временем зона с пузырями в основании нижнего слоя расширялась. Авторы установили, что при достижении этим пограничным слоем некоторой критической толщины он терял устойчивость с образованием плюмов (диапиров) и всплывал. В качестве критерия наступления конвекции выступало равенство времени развития неустойчивости Рэлея-Тейлора и характеристического времени нарастания ширины слоя с пузырями. При малой относительной толщине легкого слоя согласно Lister and Kerr (1989) максимальная скорость роста

неустойчивости а равна сг, = 0.232gAp>t. А характеристическое обратное

1 dh

время разрастания слоя сг —ч

h dt

g

= —. Поэтому, согласно Thomas et al.,

/Wf h

1993 неустойчивость наступает при толщине слоя И, при которой 01=02. Этот критерий равносилен введению локального числа Рэлея с масштабом длины, равным толщине слоя. Критическое значение такого эвристического критерия будет 4.31. На практике установлено лишь то, что расстояние между всплывающими диапирами с пузырями примерно отвечает длине волны самой быстрой ноды неустойчивости Рэлея-Тейлора (с точность 3040%). Похожее исследование провели Багдасаров и Дорфман (Bagdassarov е1 а1., 1996). Они изучили гравитационную неустойчивость оседающего слоя платиновых частиц в частично расплавленном граните. Для ускорения оседания в центрифуге создавалось ускорение 1000g. Опыты проводились в удлиненных (3:1) вдоль пути оседания ампулах. Авторы считают, что переход от простого Стоксова оседания отдельных частиц к конвективным движениям начинается при равенстве Стоксова времени, равного отношению толщины слоя к скорости оседания, и характерного времени роста развития неустойчивости Рэлея-Тейлора (ЯП). Существенным отличием оседания шариков в смеси кристаллов и расплава является большая дисперсия скорости оседания, приводящая к дополнительному механизму диссипации. Таким образом, стабилизация 1Ш понимается как следствие недостатка времени для ее развития, а не в смысле теории развития бесконечно малых возмущений. По порядку величины седиментационного числа Рэлея (около 100) переход от простого оседания к течению в опытах (Bagdassarov й а1., 1996) происходит в области критического значения, найденного нами.

УП. Численное моделирование седиментационной конвекции

Полуаналитическими методами можно исследовать гетерофазную конвекцию лишь в достаточно простых постановках. Для приближения к геологически значимым случаям необходимо численное моделирование. В работе развит метод для решения задач ^ гетерофазной конвекции, основанный на известной технике маркеров. Дополнительные возможности предусмотрены для описания дегазации и относительного перемещения фаз, что не было сделано ранее в приложениях гидродинамики к наукам о Земле.

Для решения уравнения Навье-Стокса в переменных давление-скорость был использован метод конечных элементов, в формулировке, близкой к описанным в работах (Huges et al., 1979; Fletcher, 1984; Pelletier et al., 1989; King et al., 1990; Флетчер, 1991). Для исключения давления применяется метод штрафных функций (Huges et al., 1989) Давление представляется по аналогии с уравнением состояния слабо сжимаемой жидкости через дивергенцию скорости

p = -XV.U, (1.7)

где к - произвольный большой множитель со значением до 1013. После подстановки давления из уравнения (5.1.5) в уравнение Навье-Стокса оно становится функцией только скорости течения U

/7g+V(AV-[/)+V-r=0, r = /y(Vi/+Vl/r) (2.7)

Для улучшения точности решения мы использовали алгоритм Узавы (Fortin and Fortin, 1985; Pelletier et al., 1989). Этот алгоритм позволяет получить поле скоростей с V • и = 0 даже для очень плохо обусловленных проблем, когда вязкость жидкости меняется на много порядков. Давление итерационно уточняется согласно модифицированной формуле (2.7) для шага к

^^+V(AV.[/J-V^, + V-(//(V[/t+Vi//)) = о (3.7)

давление корректируется как

= (4.7)

Начальное приближение Р0 может быть принято равным 0 или, при моделировании нестационарного процесса, обусловленного переносом плотности, значению на предыдущем временном шаге Р„. Обычно для достижения точности, близкой к машинному нулю, требуется 1-3 итерационных цикла.

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

где С, и т) - нормированные координаты в элементе е (х и у, соответственно), принимающие значения в интервале [0,1]. Коэффициенты в выражении подбираются так, чтобы обеспечить =1 в узле ] и 0 в других узлах элемента е. Давление описывалось с помощью пробных функций меньшей степени, т.е. кусочно-постоянных на прямоугольном элементе функций. Уравнения решались методом Галеркина в слабой формулировке. Эта процедура хорошо известна и приводит к матричному уравнению

(К + Я М)У„ = ) + ^ (5.7)

где К- матрица вязкости, М- матрица сжимаемости, Р - гравитационный член, а ЩРк-0 - член, отвечающий градиенту давления с предыдущего шага итерации.

Для примера приведем формулу для элемента вязкой матрицы, ассоциированную с элементом е

дх. & ду ду ду дх

дН, дЫт „ дИ. дЫ. дЫ„

I т л / т |__/__п

дх ду ду ду дх дх

Интегрирование гравитационного члена проводится с повышенной точностью по сетке 3x3, поскольку распределение плотности описывается с большим разрешением с помощью метода маркеров. Интеграл представляется в виде суммы

\pgN.d О = X ^Щ^Мт.Ъ), (6.7)

Е1-т 1=1 З.Н 3

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

Вектора точек интегрирования и весов в Гауссовом методе равны ( 1983):

4=71= [-0.7746, 0.0,0.7746] = [0.55556,0.88889,0.55556]

Решение уравнения седгшентаиионного переноса с равным нулю диффузионным членом вызывает технические трудности. Описание переноса на Эйлеровой сетке связано с появлением неизбежной численной диффузии. Для борьбы с численной диффузией применяются специальные разностные схемы повышенного порядка. Для метода конечных элементов разработана специальная модификация (Петрова-Галеркина), в которой минимизируется численная диффузия в направлении, ортогональном направлению течения. Но эти методы не исключают эффекта численной диффузии полностью.

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

fi+div(Up) = Q (7.7)

что эквивалентно -^ = 0, (8.7)

Dt v '

где плотность задана на характеристиках уравнения (7.7), а именно в

точках кривых y(t)=£(x(t)), таких что

± = U„ ^ = (9.7)

dt у dt

Скорость перемещения для дисперсной фазы равна сумме скорости жидкости и относительной скорости перемещения фаз [//+ S.

Проверка правильности работы программы. Конвективный код был успешно протестирован для модификации, описывающей термическую конвекцию. Дла проверки расчета массопереноса методом маркеров может быть использовано решение задачи о течении со свободной деформируемой поверхностью. Мы решили проблему расчета стационарного профиля соляной экструзии в условиях медленного растворения сверху. Эта задача была поставлена C.Talbot, пионером исследований растекания глубинных солевых масс, выжатых в тектонические нарушения. Подобные соляные «ледники» встречаются на юге Ирана, известны подводные экструзии на дне Мексиканского залива. Классическое решение о нестационарной форме растекающейся капли под действием силы тяжести получено X. Хаппертом (Huppert, 1982). Когда соль регулярно растворяется с поверхности выпадающими осадками и выдавливается снизу с постоянной (медленно меняющейся) скоростью, возможно стационарное состояние. Профиль двумерного потока (зависимость высоты h от расстояния до начала течения х) можно описать простым решением, найденным нами методом тонкого слоя:

h(x)=>/V2(Ъ-х), b = Q, 0<х<b, (10.7)

где безразмерный подток соли задается через линейный масштаб L и скорость растворения S: Q = Q/SL.

Для численного

решения этой проблемы мы использовали формально гетерофазный подход: солевая жидкость

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

РИС.11. Профиль соляной экструзии С несжимаемую жидкость, поверхностным растворением рассчитанный Для расчета поверхности численно с использованием модели двухслойной экструзии применялась среды с оседающими маркерами. Сплошная линия описанная выше программа, аналитическое решение в приближении тонкого основанная на технике слоя• маркеров,. Полученный

численно при сетке 30x15 стационарный профиль (см. рис. 11) оказался близок к теоретическому, описываемому формулой (10.7). Длина стационарного течения определяется соотношением баланса, так как растворение компенсируется оттоком через верхнюю границу и отклоняется от точного значения менее чем на 1%. Аспектное отношение численного решения также соответствует аналитическому решению с точностью до 1%. Численное решение отличается по форме поверхности над питающей струей, посколько при аналитическом решении течение у правой границы принято горизонтальным.

VIII. Численные модели конвекции с пузырением в магматической

камере

Группа М.Я. Френкеля (позднее А. Арискина и Е. Коптев-Дворникова) развила метод полуколичественйого моделирования динамики становления пластовых интрузий (силлов). Оперируя скоростями оседания, как подгоночными параметрами, и описывая конвекцию, как перемешивание в заданном ими интервал^ глубин, авторы в деталях описали химическую и минералогическую зональность изученных силлов. Своей деятельностью они задали высокие требования к следующему шагу в этих исследованиях. Описание конвекции должно быть как минимум двумерным, а кристаллизации - на уровне распределений кристаллов по размерам.

Менее изученной, но более важной с практической точки зрения, является проблема дегазации магмы в субвулканических камерах. Режим дегазации, в конечном счете, определяет степень эксплозивности изверженйя и его потенциальную опасность. Дегазация магмы в вулканическом процессе может происходить как при подъеме ее в магмоводе к поверхности и извержении, так и при длительной выдержке в камере под вулканом. В первом случае она может достигнуть катастрофического по скорости и последствиям характера, тогда как во втором летучие (в основном вода с небольшим содержанием 802, С02, Н2 и

др. (Бутопск й а1., 1994)) могут выделяться длительное время, достигая по трещинам поверхности при температурах, близких к магматическим. Стационарная высокотемпературная дегазация отмечается, например, на вулкане Кудрявый, Курильские острова (где газы выделяются уже в течение порядка 100 лет при температуре до 940°С) и на вулкане Сатцума-Иводзима, Япония (в течение более 1000 лет при температуре до 870°С (ТкасЬепко е1 а1., 1992; Неёепяшв^ 1994).

Описание дегазаиии в модели. Для моделирования конвекции магмы в условиях дегазации и сопряженной с ней кристаллизацией была использована программа, описанная в предыдущей главе. Главной особенностью, впервые реализованной в геофизических приложениях, было описание в технике маркеров (характеристик) процесса дегазации в приближение малой степени пузырения. Уравнения баланса для смеси могут быть упрощены. Мы предполагаем, что давление в камере может быть приближенно представлено гидростатическим давлением в слое чистого расплава, что верно при малом содержании пузырей:

где Р0 - это фоновое давление, которое определяется механическим равновесием камеры во вмещающих породах. Латеральными вариациями давления вследствие вязких напряжений мы также пренебрегаем, что верно при малых числах Рейнольдса. На данном этапе Р0 принимается постоянным, регулируемым сопряженными механическими процессами, такими как перетекание расплава и перемещение блоков по трещинам. Эта формулировка является другим крайним случаем, если принимать во внимание приближение постоянного объема. В последнем случае даже незначительные перемещения пузырей или двухфазного расплава способны вызвать неоправданно большое изменение Р0 (см., например, М!сЬаеНс1е8, 1989; Bagdassarov, 1994). Однако расширение камеры и деформация вмещающих пород могут нивелировать этот эффект.

С использованием уравнений баланса масс для компонентов гетерофазного расплава, исключив дивергенцию скорости расплава (не равную 0) и подставив давление в соответствии с гидростатической аппроксимацией, получим

где Г функция источника, а Б* - абсолютная скорость перемещения флюидной фазы

Р = Р0 + р,&

(1.8)

(2.8)

5* = 5(1-^) Г=Г,+Г2

(3.8) (4.8)

Г,=

1-£)е

(5.8)

И„+2

и, А1-*) с, (1-^)

А0+г 2(1-0,) А(1 + 2/Н„У

(6.8)

где А=рР0/р!, Н0=Р„/р£. Безразмерный параметр А является отношением характеристических плотностей газовой и жидкой фаз, Ъа - общее давление в единицах гидростатического эквивалента мощности камеры.

В (2.8) входит скорость перемещения газовой фазы 5* с учетом эффекта противотока, что имеет место и в других задачах с относительным перемещением фаз (Бтакт е1 а1., 1994).

Первое слагаемое в источниковом члене /} учитывает расширение пузырей в гидростатическом поле давления. Уравнение (2.8) только с этим членом имеет решение вдоль траектории е(г)= £о/(И0+г). Второй член описывает собственно пузырение.

При полном обороте частицы жидкости в конвективной ячейке скорость жидкости С// меняет знак, и интеграл от второго члена равен 0. Перенос пузырей со скоростью 5* ведет к вертикальному транспорту летучих компонентов и потере их через верхнюю границу.

Магматические камеры под вулканами обладают зональностью по составу. Зональность проявляется в смене состава последовательно извергаемых продуктов в одном извержении, зональности продуктов палеоизвержений, таких как туфы. Анализ состава кристаллов и расплавных включений в лавах также свидетельствует о широкой распространенности явлений смешения магм различного состава в магматических камерах. Стратификация магмы, особенно по содержанию летучих компонентов, оказывает большое влияние на режим конвекции. Нами численно промоделированы различные случаи конвекции в стратифицированных по летучему компоненту объемах расплава.

Когда содержание летучих компонентов близко к насыщению, любые перемещения расплава могут вызвать пузырение, усиливающее контраст плотности и ускоряющее вызвавшее его движения. Положительная обратная связь между композиционной конвекцией и пузырением может быть причиной катастрофической неустойчивости. Численно прослежено развитие неустойчивости слоя магмы, насыщенного во всем объеме, с учетом изменения гидростатического давления с глубиной. Конвекция инициируется за счет элемента расплава с пузырями, помещаемого в центре прямоугольной области. В результате моментально формируется центральное восходящее течение и нисходящее течение по краям. Центральная колонна с пузырями постепенно растекается у кровли, область

Различные случаи конвекции с дегазацией.

вспузыренного расплава изменяет форму от античной вазы до рюмки. В очень короткое время верхняя и нижняя части расплава меняются местами, лишь малая часть придонного расплава задерживается, что обеспечивает дальнейшую конвекцию. Моделирование проводилось в приближении постоянства давления в магме у кровли, что подразумевает синхронное увеличение объема камеры или удаление избыточного объема двухфазной жидкости через дайки. С учетом прочности вмещающих пород быстрая дегазация и увеличение объема камеры может привести к катастрофическому разрушению камеры. ...... .... г- - .......... .

^шф^Ж1

гй...........

Рис. 12. Развитие конвективной неустойчивости в двуслойной системе. Интенсивность заливки отражает степень пузырения нижнего слоя в восходящем течении. Крупные маркеры отражают положение жидкости верхнего слоя, мелкие маркеры отвечают жидкости нижнего слоя без пузырей.

дополнительным накоплением воды, нижнего слоев могут сблизиться за счет кристаллизационной дифференциации. Система становится неустойчивой, и любое возмущение (в нашем случае «капля» тяжелого расплава) ведет к возмущению границы раздела, пузырению. При развитии возмущения слои, как и в первом случае, меняются местами, объем системы резко возрастает, что также должно вести к разрушению камеры, как при извержении супервулкана.

Более_правдоподобной

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

Внедряющийся расплав может быть более горячим, более плотным, но с большим содержанием воды. Конвективное перемешивание и падение давления в предыдущих извержениях снижают содержание воды в верхнем слое. На контакте верхний слой будет нагреваться, а нижний охлаждаться и

кристаллизоваться с

Со временем плотности верхнего и

Плавная дегазаиия магмы в камере под вулканом предполагает специальный механизм конвекции расплава. Казахая (КагаЬауа е1 а1., '• . : '•■ '.'У.'.'.". 1994) предположил, что содержание •■'.■■•/.;•.•'■воды в объеме расплава меньше мН. V'I > ,' / и- 'I, уровня насыщения. Расплав

V-' поднимается по каналу, по которому

«ь 11««'лй-'! происходят извержения, до уровня

насыщения и дегазации, где расплав

пузырится, частично

кристаллизуется, теряет пузыри и

опускается в камеру. Скорость

Рис.13. Всплывание пограничного слоя похери летучих ограничивается

между двумя слоями в стратифици-____Т1

' „ 3 „ ; * сечением «трубы». Не

рованнои магматическом камере. Предпола- г

гаегся, что более плотный нижний слой анализировалось возможное

насыщен водой на уровне раздела, влияние скорости потери пузырей Двухфазный расплав пересекает на скорость конвекции. В «трубе» недонасыщенный верхний слой и степень дегазации может быть очень растекается у кровли. высока, а содержание пузырей

достигать уровня пены. Наша модель позволяет рассчитать конвекцию и в такой постановке при условии, что степень пузырения мала и действуют принятые приближения. В качестве альтернативного рассмотрен другой сценарий постепенной потери летучих в двухслойной гравитационно стабильной системе. Предполагается, что локально верхний слой насыщается водой на границе с нижним слоем (например, за счет выделения пузырей из нижнего слоя при кристаллизации или при колебании общего давления в камере). В результате конвективно-дегазационной неустойчивости этот погранслой периодически всплывает сквозь ненасыщенный флюидом расплав и, достигая кровли, теряет пузыри и питает флюидом высокотемпературные фумаролЫ.

ОСНОВНЫЕ ЗАЩИЩАЕМЫЕ ПОЛОЖЕНИЯ

1. Впервые показано, что дегазация магматических расплавов, протекающая при малом активационном пороге нуклеации, происходит равномерно в пространстве с распределением растущих пузырей по размеру, приближающемуся к степенному, фрактальному. Экспериментально установлено, что при гомогенной нуклеации и достаточно большом пересыщении возможен фронтальный, аналогичный фронту горения, режим дегазации магматических расплавов.

2. Дегазация магмы при подъеме к поверхности вызывает ее кристаллизацию, связанную с ростом ликвидусных температур кристаллических фаз. В свою очередь содержание кристаллических фаз влияет на динамику извержения. Согласно нашим экспериментальным кинетическим данным по кристаллизации гавайитов за счет дегазации, установлено, что при извержении Этны образуются кристаллы плагиоклаза и клинопироксена с максимальными размерами до 0.6-0.8 мм, ограничивающими сверху экспоненциальные участки в распределениях по размерам.

3. Получено полное точное решение задачи о кристаллизации бинарного расплава эвтектического типа с осаждением кристаллов в пограничном слое с постоянной эффективной Стоксовой скоростью осаждения и учетом противотока жидкости. При значении скорости осаждения кристаллов, соизмеримой со скоростью перемещения фронта полного затвердевания, возможно появление разрывных решений с частично поглощающим (догоняющим) кристаллы погранслоем.

4. Предложена модель кристаллизации погранслоя с оседающими кристаллами в приближении распределения кристаллов по размерам и в зависимости от расстояния до фронта полного затвердевания. Получены решения с использованием разложения распределения по дельта-функциям на траекториях, что эквивалентно конечному числу кратковременных актов нуклеации. Показано, что замедление осаждения кристаллов у границы полного затвердевания приводит к возрастанию их числа и формированию узкой зоны, обогащенной оседающим минералом в согласии с предсказаниями М.Я.Френкеля.

5. Формализм для описания гравитационного оседания частиц в гетерофазной реакционной среде переменного состава применен к решению задачи о течении в гидротермальной среде с гетерофазным (жидкость - газ) флюидом в приближении двухкомпонентной системы соль-вода. Показано, что теплоперенос ускоряется в несколько раз по сравнению с чисто водным флюидом при аналогичных перепадах температурах. Подобный механизм (одномерной) конвекции в

гетерофазном флюиде может иметь место в корневых частях активных геотермальных систем в Калифорнии, Италии и Японии.

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.

7. Симакин А.Г. и Эпельбаум М.Б. (1995). Поведение летучих' компонентов при плавлении водосодержащих пород на контакте с базитовыми расплавами// Петрология. №3. С. 349-358.

8. Симакин А.Г., Трубицын В.П. (1996). Струйный режим осаждения кристаллов и конвекция расплава в магматических камерах// Физика Земли. №7. С. 1-7.

9. Симакин А.Г. и Трубицын В.П. (1997). Захват фронтом затвердевания кристаллов растущих и оседающих в магматической камере// Физика Земли. №10. С. 3-13.

10. Симакин А.Г., Трубицын В.П. и Харыбин Е.В. (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. Simakin A.G., Trubitsyn V.P. and Schmeling H. (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., Epel'baum 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.

20. Simakih A.G. and\Trubitsyn V.P (2000). Diferentiation of two-component melt solidifying in a magma chamber// Computational seismology and Geodynamics. Vol.31. P.31-51.

21. Simakin A.G. and Botcharnikov R. (2001). Degassing of stratified magma by compositional convection// J.Volcanol. and Geotherm. 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, Chernogolovka.

2. Experimental investigation of Fsp crystallization from a granitic melt with application to the structure 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 orthomagmatic genesis, with V.A.Zharikov, M.B.Epelbaum, G.P.Zaraisky and V.S.Kudyaev; 8 IAGOD 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 XII (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. Mineral. 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.

7. Simakin A.G., Schmeling 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.; Botcharnikov, 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: steady-state operation and related heat and mass transport. Proceedings twenty-seventh Workshop Geothermal Reservoir Engineering January 28-30,2002.

»? ^ 9 5 70

РНБ Русский фонд

2006-4 20730

Содержание диссертации, доктора физико-математических наук, Симакин, Александр Геннадьевич

Введение.5

Глава 1. Физика пузырения магматического расплава

1.1 Введение 14

1.2 Теоретическая оценка условий гомогенного зарождения пузырей 17

1.3 Влияние вязких напряжений на зарождение пузырей 23

1.4 Методика экспериментов по дегазации 25

1.5 Фронтальное пузырение 28

1.6 Генеральное фрактальное распределение пузырей по размерам в эксперименте и природе 31

1.7 Наблюдение динамики распределения пузырей по размеры 36

1.8 Заключительное обсуждение и выводы к 1 главе 40-43 Приложение Решение задачи о расширяющейся вязкой сферической оболочке. 43-

Глава 2. Кристаллизация гавайитового расплава сопряженная с дегазацией

2.1 Гавайиты вулкана Этна, их экспериментальная изученность 45

2.2 Исходные материалы и методика эксперимента 47

2.3 Определение параметров диаграммы плавкости гавайитов 51

2.4 Оценки скоростей роста кристаллических фаз 55

2.5 Описание продуктов декомпрессионной кристаллизации гавайитов. 57-62 Выводы к главе 2. 62-

Глава 3. Полное решение задачи о кристаллизации бинарной системы эвтектического типа с оседанием кристаллов

3.1 Введение 64

3.2 Уравнения и граничные условия 71

3.3 Аналитическое решение для стационарного состояния 77

3.4 Соотношение между содержанием кристаллов и составом расплава в стационарном состоянии 80

3.5 Распределение температуры по глубине в стационарном состоянии 80

3.6 Анализ стационарного решения в переменных: объемное содержание кристаллов - состав расплава 82

3.7 Численное решение 85

3.8 Что говорят найденные решения о режиме кристаллизации интрузий 86-88 Заключение к главе 3 88-

Глава 4. Гидротермальная конвекция в режиме «тепловой трубы» в двухфазном двухкомпонентном флюиде

4.1 Введение 90

4.2 Основные положения модели 91

4.3 Модельные уравнения 93

4.4 Моно-вариантное стационарное решение 96

4.5 Конвекция в двухфазной (L-V) системе 99

4.6 Оценки скоростей массо-переноса и растворения 103

4.7 Сопоставление с геологическими данными 104-106 Выводы к главе 4 106-

Глава 5. Описание кристаллизации бинарного расплава с оседанием кристаллов на уровне распределения по размерам

5.1 Введение

5.2 Описание модели 108

5.3 Поглощающий кристаллы температурный погранслой 113

5.4 Траектории кристаллов, выпадающих на дно камер 110

5.5 Применимость равновесного приближения для магматических камер 112

Введение Диссертация по наукам о земле, на тему "Моделирование фазовых переходов и массопереноса в магматических камерах"

Актуальность работы. Земная кора является существенно гетерофазным объектом. Фазовые переходы и относительное перемещение жидких, кристаллических и газообразных фаз в гравитационном поле, особенно сегрегация, перемещение и застывание магматического расплава составляют основу процессов дифференциации вещества Земной коры. При малой доле расплава относительное движение кристаллической и жидкой фаз подчиняется законам компакции, т.е. вязкого гравитационного течения кристаллической пористой матрицы (Khodakovskii et al., 1995; MacKenzie, 1987) или направляется наложенными девиаторными деформациями (Simakin and Petford, 2003; Simakin and Talbot, 2002a; Simakin and Talbot 2002b; Simakin and Ghassemi, 2003).

Дифференциация достигает большого масштаба в магматических системах с высоким содержанием расплава. В процессе магматической дифференциации за счет гравитационной отсадки кристаллов и отделения и всплывания газовой фазы начальный состав магмы может существенно измениться. Классическими продуктами процесса дифференциации являются расслоенные интрузии, такие как Скаергаардская интрузия в Гренландии (Wager and Brown, 1967), Йоко-Довыренский ультрабазит-базитовый массив в Забайкалье (Кислов, 1998) и многие другие. Последовательность составов пород, связанных процессами ассимиляции корового вещества и фракционной кристаллизации безотносительно их пространственного распределения описывается на качественном уровне с помощью AFC (assimilatio-fractional crystallization model) моделей (см., например, Bohrson and Spera, 2001)

Моделирование эволюции и структуры магматических камер на строгом физическом уровне - достаточно сложная задача. На начальной стадии В.С.Голубев и В.Н.Шарапов (Голубев и Шарапов, 1974; Шарапов и Черепанов, 1986) моделировали засывание интрузивов как направленную кристаллизацию при диффузионном массообмене в расплаве. Серия работ по полуколичественному моделированию формирования пластовых интрузий в ходе отсадки кристаллов из конвектирующего расплава была выполнена М.Я.Френкелем (Френкель, 1995). Отдельные изолированные попытки аналогичных расчетов были предприняты также на Западе (Oldenburg and Spera, 1991; Worster and Huppert, 1993; Jaupart and Tait, 1995). Между тем, в настоящее время признано, что интрузивные породы приобретают свой облик, главным образом, в ходе перекристаллизации в состоянии кумулуса, т.е. пористой среды с малой долей расплава (Park and Hanson, 1999). Именно в эту стадию, скорее всего, происходит формирование магматических руд и ритмической расслоенности. Поэтому набольший интерес к расчетам процессов гетерофазной конвекции кристаллизующейся магмы проявляют вулканологи, т.к. вулканические лавы прибывают в камере и изливаются при относительно небольшой степени кристалличности, и масштаб времени эволюции вулканических камер достаточно мал и сопоставим со скоростью конвективных процессов. Активность и актуальность исследований в области вулканологии очевидна в силу реальной угрозы вулканических извержений для человека.

Развитие современных вычислительных методов, а также экспериментальных данных по кинетике и динамике фазовых превращений (кристаллизации и пузырения) в магматических расплавах позволяет подойти к количественному физическому моделированию процессов массо- и теплопереноса в природных гетерофазных магмах на более высоком уровне. Подходы и методы, развиваемые для описания гетерофазных магматических камер, применимы и для геотермальных систем с многофазным флюидом, исследования которых становятся все более актуальными в связи с увеличением глубины и солености эксплуатируемых геотерм (Barelli et al., 2000; Kasai et al., 2000).

Цель работы. Основной целью данной работы было изучение физических механизмов, контролирующих процессы фазовых превращений, относительных перемещений фаз в гетерофазных магмах и сопряженной с этими процессами конвекции. Для достижения поставленной цели потребовалось решить следующие принципиальные научные задачи:

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

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

- получить точное решение обобщенной задачи Стефана в одномерной постановке, учитывающей объемный рост надэвтектической кристаллической фазы, оседающей в поле тяжести, в приближении постоянной, эффективной Стоксовой скорости осаждения,

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

ЫаСЛ-НгО), сравнить оценки тепловых и массовых потоков с наблюдаемыми величинами,

- описать направленное затвердевание' расплава с оседающими кристаллами с применением функций распределения по размерам,

- рассмотреть теоретически различные типы гетерофазной конвекции с учетом относительного перемещения дисперсной фазы (пузырей, кристаллов),

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

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

Научная новизна. Предлагаемая к защите работа является первым междисциплинарным экспериментальным и теоретически-численным исследованием процессов дегазации и кристаллизации и сопряженного массопереноса в водонасыщенной магме в микро и макро масштабах. Все приведенные в заключении результаты исследований являются оригинальными. На основе полученных результатов можно сформулировать основные защищаемые положения:

- Впервые показано, что дегазация магматических расплавов, протекающая при малом активационном пороге нуклеации, происходит равномерно в пространстве с распределением растущих пузырей по размеру, приближающемуся к степенному (фрактальному). Экспериментально установлено, что при гомогенной нуклеации и достаточно большом пересыщении возможен фронтальный, аналогичный фронту горения, режим дегазации магматических расплавов.

- Впервые экспериментально установлено, что при кристаллизации магмы гавайитов, сопряженной с дегазацией в жерле вулкана Этна, образуется экспоненциальный участок распределения кристаллов плагиоклаза и пироксена по размерам. Максимальные размеры кристаллов плагиоклаза и клинопроксена могут достигать размеров фенокристов: 0.3-0.8 мм.

- Получено полное точное решение задачи о кристаллизации бинарного расплава эвтектического типа с осаждением кристаллов в пограничном слое с постоянной эффективной Стоксовой скоростью осаждения и учетом противотока жидкости. При значении скорости осаждения кристаллов, соизмеримой со скоростью перемещения фронта полного затвердевания, возможно появление разрывных решений с частично поглощающим (догоняющим) кристаллы погранслоем.

Предложена модель кристаллизации погранслоя с оседающими кристаллами в приближении распределения кристаллов по размерам и в зависимости от расстояния до фронта полного затвердевания. Получены решения с использованием разложения распределения по дельта-функциям на траекториях, что эквавалентно конечному числу кратковременных актов нуклеации. Показано, что замедление осаждения кристаллов у границы полного затвердевания приводит к возрастанию их числа и формированию узкой зоны, обогащенной оседающим минералом в согласии с предсказаниями М.Я. Френкеля.

Формализм для описания гравитационного оседания частиц в гетерофазной реакционной среде переменного состава применен к решению задачи о течении в гидротермальной среде с гетерофазным (жидкость - газ) флюидом в приближении двухкомпонентной системы соль-вода. Показано, что теплоперенос ускоряется в несколько раз по сравнению с чисто водным флюидом при аналогичных перепадах температурах. Подобный механизм (одномерной) конвекции в гетерофазном флюиде верятно имеет место в корневых частях активных геотермальных систем в Калифорнии, Италии и Японии.

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

Метод маркеров обобщен для решения проблемы адвективного переноса с реакцией в гетерофазной жидкости. На его основе разработан численный код для моделирования гетерофазной конвекции, сопряженной с дегазацией и кристаллизаций. Проанализированы различные сценарии развития конвекции в результате пузырения и кристаллизации в двухслойной магматической камере под вулканом с нижним слоем, содержащим больше, чем верхний, растворенных летучих компонентов.

Поставленные задачи, развитые теории и методология исследования являются по существу новым научным направлением: комплексным экспериментально-теоретическим исследованием кинетики и динамики сопряженных процессов фазовых превращений и тепло и массопереноса в природных магмах.

Теоретическая и практическая значимость. Результаты работы имеют как теоретическую, так и практическую. значимость. При изучении физики дегазации было теоретически объяснено явление фронтального пузырения за счет падения давления в зоне вязкого течения вблизи пузырей субмикронного радиуса. Теоретический интерес имеет выявленное фрактальное распределение генеральной совокупности пузырей (как в экспериментальных, так и в природных образцах), оно послужило отправной точкой к моделированию эволюции распределений в ходе множества дискретных актов нуклеации в группе (K.Mader, Bristol, UK). Отмеченный в эксперименте процесс нуклеации водных пузырей на микропузырях газов с меньшей растворимостью, обладающий малым активационным порогом, является важным механизмом гетерогенной нуклеации в природных условиях, придающий ему плавный характер. Принципиальное значение имеет впервые выявленная возможность существования порогового содержания растворенной в расплаве воды, ниже которого гомогенная нуклеация становится невозможной даже при падении давления до атмосферного.

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

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

Применение подхода, использованного для описания гетерогенных магматических расплавов, к двухфазному флюиду позволило получить практически важные оценки скорости растворения пористой среды и максимального теплового потока, переносимого «тепловой трубой» с соленым двухфазным флюидом.

Теоретически значимыми являются оценки критического седиментационного числа Рэлея, полученные методом граничного интеграла для развитого течения. Наряду с ранее изученной мгновенной конвективной неустойчивостью градиентного распределения оседающих частиц в слое жидкости, наши результаты ясно демонстрируют сходство конвекции в средах с адвективным переносом дисперсионной фазы с классическим тепловым течением. Введено определение седиментационного числа Нуссельта, найдена его связь с седиментационным числом Рэлея.

Практическое значение имеет программа для расчета конвекции, сопряженной с кристаллизацией в магматической камере под вулканом. Продемонстрировано, что в двухслойной системе при достижении близости плотностей слоев происходит перемешивание с самоускоряющейся дегазацией и разрушением вулканической камеры. Этот процесс, видимо, протекал при извержении супервулканов, поэтому предложенную модель можно использовать для прогноза таких извержений и анализа возможных путей воздействия на них. Программа допускает расчет параметров конвекции в широком магмо-подводящем канале с дегазацией и частичной кристаллизацией в верхней ее части по модели (КагаЬауа е1 а1., 1994).

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

Методика исследования. В первых двух главах главные фазовые переходы в магме на микроуровне (пузырение и кристаллизация) рассмотрены на экспериментальном материале. А в дальнейшем теоретически рассмотрены режимы их протекания на макроуровне в поле силы тяжести. В диссертации для достижения цели работы применялись как математические (метод конечных элементов, граничных интегралов, тонкого слоя, метода анализа стационарного состояния, алгоритм Узавы, метод маркеров для численного моделирования конвекции в гетерофазной среде и течений со свободной поверхностью), так и физико-химические методы экспериментальной петрологии (метод затравок для изучения кинетики роста и нахождения индивидуального ликвидуса фаз, в общем случае метастабильного, методы анализа изображений для определения распределения пузырей по размерам).

Экспериментальная часть работы проделана в содружестве с М.Б.Эпельбаумом и Т.П. Саповой, которые непосредственно в ИЭМ РАН проводили запланированные автором эксперименты на установках высокого газового давления. Препараты из образцов после опытов, включая шлифы и двусторонние полированные пластинки для ИК-измерений, изготовлялись лично автором. Измерения состава кристаллических фаз и стекол проводились на рентгеновском микроанализаторе, главным образом в ИЭМ РАН, при решающем участии А.Некрасова. Небольшое количество контрольных определений содержания воды в стеклах методом ПЖ спектроскопии проводилось в университетах г. Пиза (Италия) при участии П. Армиенти и г. Ганновер (Германия) с помощью Р. Бочарникова. Численные расчеты конвекции проводились по алгоритмам, опубликованным в литературе. Оригинальным является описание процесса пузырения и относительного .перемещения пузырей Лагранжевым методом на характеристиках уравнения переноса. Также оригинальным является аппроксимация распределения кристаллов по размерам набором дельта-функций при решении задачи о кристаллизации в пограничном слое с оседающими кристаллами.

Апробация работы. Основные положения и отдельные разделы работы неоднократно обсуждались на семинарах: семинаре Института наук о Земле Уппсальского Университета (Швеция), на совместном семинаре Национальной Вулканологической Группы и геологического факультета ун-та г. Пиза (Италия). Результаты, изложенные в диссертации, докладывались на следующих международных рабочих совещаниях:

На 16 Генеральной ассамблее Международной Минералогической Ассоциации в г.Пиза, Италия, 1994 год; на 8 сессии Европейского Союза Наук о Земле в Страсбурге, 1995; на Гольдшмитовской конференции по Геохимии в Гейдельберге, 1996; на сессии Европейского Геофизического союза в Ницце, 2000; на семинаре по Геотермальной энергии в Стэнфордском Университете, 2002.

Публикации. По результатам выполненных исследований опубликовано около 30 печатных работ, в том числе 19 статей в ведущих российских и зарубежных журналах, можно отметить 10 тезисов наиболее значимых докладов.

Структура и объем работы. Диссертация состоит из краткого введения, восьми глав, разбитых на секции, и заключения. Каждой главе предшествует подробное введение. Секции, в свою очередь, разбиты на разделы и подразделы. Нумерация секций сквозная. Работа содержит 197 стр. текста и включает 6 таблиц и 75 рисунков. Список литературы содержит 236 наименований.

Личный вклад автора. Диссертация основана на результатах самостоятельных исследований автора. Экспериментальная часть работы проделана под руководством автора в содружестве с сотрудниками и механиками лаборатории магматизма ИЭМ РАН. Многие из теоретических вопросов разрабатывались совместно с научным консультантом автора член-корр. В.П. Трубицыным, который своей повышенной требовательностью к строгости и ясности решений сыграл неоценимую роль, в ходе работы над диссертацией.

Краткое содержание работы.

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

В главе I приводятся оценки скорости гомогенной нуклеации пузырей водного флюида в кислом расплаве. Представлены данные по распределению пузырей по размерам в кислых расплавах, полученные в закалочных эксперименах при сбросе давления на аппаратах высокого газового давления. Эти распределения сравниваются с природными данными по продуктам извержения итальянских вулканов.

В главе II описаны эксперименты по кристаллизации гавайитовой магмы при малом давлении воды (600 бар), отвечающим условиям насыщения водой магм вулкана Этна. Предварительно уточнена диаграмма плавкости гавайитов. Получены оценки скоростей роста магматических минералов. Экспериментальные данные позволили интерпретировать распределения по размерам кристаллов пироксена и плагиоклаза в вулканических породах Этны.

В III-V главах теоретически рассмотрены задачи об одномерных течениях в реакционных гетерофазных жидкостях с оседанием дисперсной фазы.

В III главе рассмотрено строение двухфазной зоны на фронте затвердевания в магматическом расплаве в равновесном приближении. Проанализированы решения при различном соотношении скоростей оседания (всплывания) кристаллов и скорости перемещения фронта для стацианарного случая, получены численные нестационарные решения.

В главе IV техника описания пограничного слоя в магматической камере применена для анализа режимов одномерного течения с оседанием тяжелой компоненты в двухфазном флюиде в системе соль-вода. Приводятся оценки теплового потока в стационарном состоянии таких систем. Рассмотрены геологические примеры.

В главе V приводятся результаты численного решения задачи о кристаллизации в погранслое с оседающими кристаллами в распределениях по размерам.

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

В главе VI теоретически рассмотрены течения, которые могут возникать в гетерофазных магмах при Стоксовом переносе дисперсных фаз (пузырей, кристаллов). Найдены параметры стационарных течений при поступлении дисперсной фазы на горизонтальных границах объема магмы.

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

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

Заключение Диссертация по теме "Геофизика, геофизические методы поисков полезных ископаемых", Симакин, Александр Геннадьевич

Выводы к главе 8.

1. Впервые метод маркеров применен не только к традиционному описанию перемещения границ в многослойной среде, но к решению гиперболического уравнения переноса с источниковыми членами. Это позволило описать конвекцию в гетерогенной магме с пузырением.

2. Методом математического моделирования показано, что в многослойной магматической камере возможны интересные гидродинамические явления. При устойчивой плотностной стратификации в нижнем слое магмы ожидается повышенное содержание растворенных летучих компонентов. При продолжающейся фракционной кристаллизации плотности слоев могут сравняться, что вызовет бурную конвекцию с пузырением, всплывание нижнего слоя, облегчающегося при выделении пузырей. Такая конвективно- дегазационная неустойчивость, вероятно, является механизмом извержений супервулканов на Земле, способных за короткое время выбросить в атмосферу до нескольких сотен (тысяч) кубокилометров пирокластики.

3. При сохранении устойчивой стратификации нижний слой кислой магмы (вероятно подогреваемый периодическими поступлениями основной магмы) может служить источником выделения флюидов длительное время (от 100 до 1000 лет) в некоторых вулканах зоны субдукции (например, вулкана Кудрявый, Итуруп).

ЗАКЛЮЧЕНИЕ

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

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

2) Предложена строгая самосогласованная модель объемной кристаллизации магмы с оседающими кристаллами в термальном погранслое у Стефановского фронта полного затвердевания. Получено полное множество квази-стационарных решений проблемы при различных относительных скоростях оседания кристаллов.

3) Обобщение модели объемной кристаллизации магмы к описанию одномерного течения в гетерофазном флюиде системы соль-вода типа «тепловой трубы» позволило получить оценки тепловых и массовых потоков в корневых зонах геотермальных систем Италии, Японии и США.

4) Конвективная неустойчивость жидкости с дисперсной фазой, перемещающейся в гравитационном поле, отличается от неустойчивости Рэлея-Тэйлора и способна перейти в стационарную конвекцию, подобную тепловой.

5) Разработан смешанный Эйлерово-Лагранжев алгоритм для численного моделирования конвекции с пузырением и перемещением пузырей флюидной фазы.

6) Проведено численное моделирование конвекции с пузырением в двухслойной стратифицированной магматической камере. Выявлены режимы плавной и катастрофической дегазации.

Полученные результаты дают возможность сделать следующие выводы.

1. Полученные экспериментальные данные свидетельствуют о том, что степенной закон распределения генеральной совокупности пузырей по размеру свидетельствует о распространенности в природе гетерогенной нуклеации пузырей водного состава с малым активационным порогом.

2. Экспоненциальные распределения по размерам кристаллов плагиоклаза и клинопироксена в лавах Этны связаны с кристаллизацией, сопряженной с дегазацией при подъеме магмы к поверхности в условиях нарастания скорости нуклеации при падении содержания воды в расплаве.

3. В приближении равновесной кристаллизации на фронте полного затвердевания распределение объемного содержания вкрапленников имеет разрыв. Их содержания в твердой породе и магме связанны с величинами скоростей продвижения фронта и оседания вкрапленников. В распределении оседающих кристаллов по глубине также возможны разрывы, связанные с противотоком жидкости и с зависимостью эффективной вязкости смеси от содержания кристаллов. Разрывы являются источником решений солитонного типа.

4. Вследствие зависимости вязкости смеси (как следствие и скорости оседания) от содержания кристаллов в области эффективных скоростей оседания, близких к скорости фронта полного затвердевания, появляется множество квазистационарных решений.

5. Корректное решение задачи о кристаллизации магмы в верхнем термальном погранслое возможно лишь на уровне распределений кристаллов по размерам методами, разработанными в данной диссертации.

6. Следствием двухфазной конвекции типа тепловой трубы в реактивной среде может быть резкое возрастание скорости теплопереноса вплоть до 10-13 ват/м2 при проницаемости 10"15м2, а также явления геотермального карста с оседанием поверхности со скоростью до нескольких см/год.

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

8. Продемонстрировано, что стратифицированная магматическая камера под вулканом например, андезито-базальтового типа) может служить источником катастрофического извержения объемного типа, когда плотности нижнего и верхнего слоев сравняются естественным путем, а также обеспечивать постепенную дегазацию с длительными высокими потоками магматического флюида.

БЛАГОДАРНОСТИ

Я посвящаю эту работу своим родителям, которые предопределили мой глубокий интерес к геологии, проведя золотую пору своей юности на комбинате «Маяк» на берегу чудесного озера Иртяш. Благодаря им я появился на свет среди минералогических и петрологических чудес седого Урала. Они предоставили мне свободу выбора и смирились с геологией и скитальческой перспективой моей будущей жизни.

А далее судьба всегда предоставляла мне счастливый шанс на пути к знанию. Почти случайно я попал в ИЭМ РАН, после нескольких итераций оказавшись в лаборатории магматизма, возглавляемой мудрым М.Б. Эпельбаумом, где я прошел школу молодого бойца и провел уже более 20 лет. С вечно элегантным ветераном лаборатории Т.П.Саловой мы успешно продолжаем высокие традиции изучения расплавов с летучими компонентами.

Также почти случайно, благодаря доброй воле нашего директора академика В.А. Жарикова, я познакомился с член.-корр. В.П. Трубицыным. Он помог мне преодолеть робость выпускника группы геохимических поисков геологического факультета МГУ перед людьми, одолевшими премудрости точных наук в стенах физического и математического факультетов. Под чутким руководством В.П.Трубицына среди сотрудников его лаб. Института Физики Земли РАН (разлетевшихся кто куда в урагане последней русской революции) Е. Харыбина, О. Белавиной, О. Зацепиной, А. Боброва я рос понемногу как геофизик. Наш совместный с В.П. Трубицыным визит в Байрот в 1992 году и совместная работа с X. Шмелингом повлекли незатухающий интерес к седиментационной конвекции и методу маркеров.

Ключевой оказалась помощь Л.Л.Перчука, который заполнил мною вакансию в институте геологических наук при Упсальском университете, образовавшуюся благодаря исходу Ю.Подладчикова в Голландию. Там среди древностей земли викингов я провел счастливое время, использовав предоставленную возможность проникнуть в премудрости численных методов и механику пористой вязкоупругой среды с разрывами. Безукоризненный английский язык К.Талбота, зав. лаборатории геодинамики им. Рамберга, урожденного шотландца, сыграл решающую роль в публикации и дальнейшей жизни нескольких важных работ. В ту пору сотрудники его лаборатории С. Медведев и Б. Шотт приобщили меня к пользованию пакетом программ символьной математики МАРЬЕ и стандарту изложения своих мыслей с помощью ЬаТеХ, без которых я не представляю себе жизни отныне.

Я глубоко признателен П. Армиенти, который ввел меня в мир вулканологии и познакомил с чудесной страной Италией. Именно он и Ф. Инносенти инициировали исследования пузырения и сопряженной дегазации и кристаллизации. Его интуиция и полевые данные позволили приблизиться к пониманию природы распределений пузырей по размерам. Дискуссии с сотрудниками национальной лаборатории вулканологии в г. Пиза П.Папале, С.Нери, М. Росси были очень продуктивными и запоминающимися. Интерес к вулканам, в частности, к вулкану Кудрявый (о. Итуруп) закрепился при общении с К.И. Шмуловичем, который поставил задачу о возможных режимах дегазации, не сопровождающейся извержением. Вокруг К.И.Шмуловича сложилась группа ценителей и любителей фумарол и дегазации в лице М.Коржинского и Р.Бочарникова, общение с которыми также много дало автору.

В последние годы судьба закинула меня в прерии Сев. Дакоты, где на факультете инженерной геологии я прикоснулся к миру прикладной науки, к миру геотермальной энергии благодаря финансовой поддержке и сотрудничеству с проф. А. Гассеми. Тогда абстракция уравнений модели пластовой интрузии породила оценки параметров «тепловой трубы» в корнях геотермальных систем в гетерофазном рассоле.

В эти трудные для российской науки годы я был поддержан Российским Фондом Фундаментальных Исследований. Я также признателен за гостеприимство и финансовую поддержку, оказанные мне во время научных командировок в Италию (П.Армиенти и Ф.Инносенти), в Швецию (К.Талбот), в Германию (Х.Шмелинг, Ф.Нольц), в США (А.Гассеми), в Испанию (Х.Ринкон) и Англию (Дж.Клеменс).

Библиография Диссертация по наукам о земле, доктора физико-математических наук, Симакин, Александр Геннадьевич, Черноголовка

1. Биндеман И.Н. (1995) Ретроградная везикуляция основной магмы в приповерхностных магматических камерах модель происхождения мафических включений в кислых и средних породах. Петрология. 3(№6)576-587.

2. Гардинер К.В. (1986) Стохастические методы в естественных науках. М. Мир.

3. Глико и А.О., Петрунин А.Г. (1998) Моделирование эволюции тепломассопереноса в системе "черный курильщик- магматическая камера". Физика Земли. №7,3-10.

4. Зеленко В. и Мясников В. (1985) К теории круговых движений в барботажном слое. Изв. АН СССР, Механика Жидкости и Газа, №5,108-115.

5. Кислое Е.В. (1998) Йоко-Довыренский расслоенный массив. Издательство БНЦ, Улан-Удэ, 264 стр.

6. Лойцянский Л.Г. (1987) Механика жидкости и газа. М.Наука, 840 стр.

7. Мелихов КВ., Берлинер Л.Б. и Слинько М.Г. (1985) Влияние дисперсии скорости роста кристаллов на кинетику массовой кристаллизации. ДАН СССР, 283 (№4)917.

8. Мошинский А. И. (1991) Некоторые закономерности непрерывной кристаллизации солей из раствора. Теор.основы хим. Технологии, 25 (№2) 219-226.

9. Мясников В.П., Митронов А.П., Кочергин H.A., Дилъман В.В. (1983) Структура потоков в высоком непроточном пузырьковом слое. ДАН СССР, 269 (№4) 827-830.

10. Неймарк Б.М. и Исмаил-заде А.Т. (1994) Усовершенствованная модель погружения тяжелых тел в астеносфере. Теор. проблемы геодинамики и сейсмологии. М., Наука, 56-69 стр.

11. П. Персиков Э.С. (1984) Вязкость магматических расплавов. М.Наука. 160 стр.

12. Симакин, А.Г. (1993) Некоторые кинетические закономерности роста магматических минералов. Петрология. 1 (№5) 550-556.

13. Симакин А.Г., Салова Т.П., Эпельбаум М.Б., Бондаренко Г.В. (1998) Влияние воды на структуру промежуточного диапазона в алюмосиликатном расплаве. Геохимия (8) 768-773.

14. Симакин А.Г., Салова Т.П. (2000) Динамические микро-режимы пузырения гранитного расплава. Физика Земли (4) 57-63.

15. Симакин А.Г., Салова Т.П. (2001) Экспериментальные данные по эволюции распределения пузырей по размеру при плавной дегазации водонасыщенного гранитного расплава. Геохимия (3) 294-304.

16. Симакин А.Г., Салова Т.П., Армиенти П. (2003) Кинетика роста клинопироксена из водосодержащего гавайитового расплава. Геохимия (12) 12751285.

17. Симакин А.Г., Армиенти П., Салова Т.П. (2000) Сопряженные дегазация и кристаллизация: экспериментальное изучение при плавном снижении давления. Геохимия (6) 579-591.

18. Симакин А.Г., Салова Т.П. (2003) Кристаллизация плагиоклаза из гавайитового расплава в эксперименте и в вулканическом канале. Петрология (12) 98109.

19. Симакин А.Г., Трубицын В.П., Харыбин Е.В. (1998) Распределение по размерам и глубине для кристаллов, осаждающихся в застывающей магматической камере. Физика Земли (8) 30-37.

20. Симакин Ф.Г., Трубицын В.П. (1995) Эволюция структуры остывающей магматической камеры. Физика Земли (2) 1-13.

21. Симакин А.Г. и Трубицын В.П. (1997) Захват фронтом застывания растущих и оседающих в магматической камере кристаллов. Физика Земли (10) 3-13.

22. Симакин А.Г. и Трубицын В.П. (1993) Структура верхнего пограничного слоя кристаллизующегося интрузива при седиментации кристаллов. Физика Земли (4) 2029.

23. Симакин А.Г. и Трубицын В.П. (1996) Струйный режим осаждения кристаллов и конвекция расплава в магматических камерах. Физика Земли (8) 26-32.

24. Симакин А.Г. и Трубицын В.П. (2000) Эффекты композиционной конвекции при дегазации стратифицированной магмы. Физика Земли (11) 57-69.

25. Симакин А.Г. и Трубицын В.П. (2000) Дифференциация двухкомпонентного расплава в кристаллизующейся магматической камере. Вычислительная сейсмология (31)31-51.

26. Симакин А.Г. и Худяев В.С. (1992) Механика формирования гранитных интрузий. В сб. (Ред. Г.Р.Колонин) Генетическая модель грейзеновых рудных месторождений. М.: Наука, Новосибирск. Стр.137-154.

27. СимакинА.Г.и Чевычелов В.Ю. (1995) Экспериментальное изучение кристаллизации Fsp из гранитного расплава с различным содержанием воды. Геохимия (2)216-237.

28. Теркот Д. и Шуберт Дж. Геодинамика. Геологические приложения физики сплошных сред. М. Мир, 1985,730 стр.

29. Тимошенко С.П. Курс теории упругости. Киев, Наукова думка, 1972.

30. Ткаченко С.И., Таран Ю.А., Коржинский М.А., Покровский В.Г. (1992) Газовые струи вулкана Кудрявый, Итуруп, Курильские острова. Докл. РАН, 325 (№4) 823-828.

31. Трубицын В.П. (2000) Фазовые переходы, сжимаемость, тепловое расширение, теплоемкость и адиабатическая температура в мантии. Известия РАН, Физика Земли, 36, №2.

32. Трубицын В.П. и Харыбин Е.Б. (1987) Конвективная неустойчивость режима седиментации в мантии. Изв. АН СССР, сер. Физика Земли (7) 21-30.

33. Трубицын В.П. и Харыбин Е.Б. (1988) Гидродинамическая модель дифференциации вещества в Земле. Изв. АН СССР, сер. Физика Земли (4) 83-86.

34. Трубицын В.П. и Харыбин Е.Б. (1991) Термоседиментационная конвективная неустойчивость двухкомпонентной вязкой жидкости. Изв. АН СССР, сер. Физика Земли (2) 3-17.

35. ФерриДж. и Баумгартнер Л. (1992) Термодинамические модели молекулярных флюидов при повышенных температурах и давлениях метаморфизма. В сб. Термодинамическое моделирование в геологии. Минералы, флюиды, расплавы. М.: Мир, 354-389.

36. Флетчер К. (1991) Вычислительные методы в динамике жидкостей, т.1 Основные положения и общие методы. Мир, Москва.

37. Френкель М.Я. (1995) Термическая и химическая динамика дифференциации базитовых магм. М.Наука, 239 стр.

38. Френкель М.Я., Ярошевский A.A., Арискин A.A., Бармина Г.С., Коптев-Дворников Е.В., Киреев Б.С. (1988) Динамика внутрикамерной дифференциации базитовых магм. М.: Наука. 216 с.

39. Чевычелов В.Ю., Симакин А.Г. (2003) О механизме растворения хлора в модельном водонасыщенном гранодиоритовом расплаве: использование методов ИК спектроскопии. Геохимия (4) 1-17.

40. Чехмир А., Персиков Э.С., Эпельбаум М. и Бухтияров П. (1985) Транспорт водорода через модельный магматический расплав. Геохимия (5) 594-598.

41. Чехмир А.С., Симакин А.Г. (1991) Динамические явления во флгаидно-магматических системах. М.: Наука, 141 стр.

42. Шарапов В.Н. и Черепанов А.Н. (1986) Динамика дифференциации магм. Новосибирск, Наука, 184 стр.

43. Шарапов В.Н., Акимцев В.А., Доровский В.Н., Перепечко Ю.В. и Черепанов А.Н. (2000) Динамика развития пудно-магматических систем зон спрединга. Новосибирск, 414 стр.

44. Шарков Е.В. Петрология расслоенных интрузий. Л. Наука, 1980,174 стр.

45. Эпельбаум М.Б. (1980) Силикатные расплавы с летучими компонентами. М.: Наука, 1980.242 стр.

46. Эпельбаум М.Б., Иванов М.А. и Фокеев Е.В. (1991) Многоампульная установка высокого газового давления с револьверным устройством для быстрой закалки. Очерки физико-химической петрологии. М. Наука, Вып. 17,141-144.

47. Aldibirov М., Dingwell D.B. (1996) Magma fragmentation by rapid decompression. Nature (380)146-149.

48. Armienti P., Innocenti F., Petrini R., Pompilio M., Villari L. (1988) Sub-aphyric alkali basalt from Etna: interference on depth and composition of the source magma. Rend. Soc. It. Min. Petr. (43) 877-891.

49. Armienti P., Innocenti F., Pareschi M.T., Pompilio M., Rocchi S. (1991) Crystal population density in not stationary volcanic system: estimate of olivine growth rate in basalts of Lanzarote (Canary Islands). Contr. Mineral. Petrol. (44) 181-196.

50. Armienti P.Pareschi M.T., Innocenti F. and Pompilio M. (1994) Effects of magma storage and ascent on the kinetics of crystal growth. The case of 1991-93 Mt. Etna eruption. Contrib. Mineral. Petrol. (115) 402-414.

51. Armienti P., Pareschi M.T. and Pompilio M. (1997) Lava textures and time scales of magma storage at Mt.Etna. Acta vulcanol. (9) 1-5.

52. Ariskin A.A., Barmina G.S. and Frenkel M.Ya. (1987) Computer simulation of basalt magma crystallization at a fixed oxygen fugacity. Geoch.Int. (24) 92-100.

53. Bagdassarov N.S. and Dingwell D.B.(\ 992) A rheological investigation of vesicular rhyolite. J.Volcan.Geotherm. Res. (50) 307-322.

54. Bagdassarov N.S., Dingwell D.B. and Wilding M.C. (1996) Rhyolite magma degasing: the experimental study of melt vesiculation. Bull. Volcanol. (57) 587-601.

55. Bagdassarov N.S. and Dor/man A. (1996) Modelling of melt segregation process by high-temperature centrifuging of partially molten granites -II. Rayleigh-Taylor instability and sedimentation. Geophys. J. Int. (127) 627-634.

56. Barelli A., Bertini G., Buonasorte G., Cappetti G., and Fiordelisi A. (2000) Recent deep exploration results at the margins of the Larderello-Travale geothermal system, in: The World Geothermal Congress Japan Proc., 965-970.

57. Barmin A., Melnik O., Sparks R.S.J. (2002) Periodic behaviour in lava dome eruptions, Earth Planet. Sci. Lett. (199) 173 184.

58. Bercovici D„ Richard Y. and Schubert G. (2001) A two-phase model for compaction and damage I. General Theory. J.Geophys.Res. (106,B5) 8887-8906

59. Bergantz G.W. (1995) Changing techniques and paradigms for the evolution of magmatic process. JGR (100/B9) 17.603-17.613.

60. Bergantz G.W. and Ni J. (1999) A numerical study of sedimentation by dripping instabilities in viscous fluids. Int. J. of Multiphase Flow (25) 307-320.

61. Bergantz G.W. (2000) On the dynamics of magma mixing by reintrusion: implications for pluton assembly processes J.Structural Geol. (22) 1297-1309.

62. BischoffJ.L. (1991) Densities of liquids and vapors in boiling NaCl-H20 solutions: a PVTX summary from 300 to 500°C, Amer. J. Sci. (291) 309-338.

63. Bitner D., Schmeling H. (1995) Numerical modeling of melting processes and induced diapirism in the lower crust. Geophys.J.Int. (123)59-70.

64. Blower, J.D., Keating J.P., Mader H.M. and Phillips J. C. (2001) Inferring volcanic degassing processes from vesicle size distributions, Geophysical Research Letters 28(2).

65. Botcharnikov R.E., Knyazik V.A., Steinberg A.S. and Steinberg G.S. (1995) Measurements of outgassing velosity from fumarole outlets, Kudriavy volcano, Iturup, Kuril Islands. Proceed. Symp. Volc.-Atmosph. Interactions, Honolulu, Hawaii. 27-34.

66. Brandeis G., Jaupart C. and Allegre C. (1984) Nucleation, crystal growth and thermal regime of cooling magmas. J.Geophys.Res. (89, #B12) 10.161-10177.

67. Bohrson, W. A. and Spera, F. J. (2001). Energy-constrained open-system magmatic processes II: application of energy-constrained assimilation-fractional crystallization (EC-AFC) model to magmatic systems, Journal of Petrology (42) 1019-1041.

68. Brandeis G., Jaupart C. (1986) On the interaction between convection and crystallization in cooling magma chambers. Earth and Planet. Sci.Lett. (77) 345-361.

69. Carminati E., Augier F.T., and Barba S. (2001) Dynamic modelling of stress accumulation in Central Italy: role of structural heterogeneities and rheology. Geophys. J. Int. (144) 373-390.

70. Cashman K., Marsh B. (1988) Crystal size distribution (CSD) in rocks and the kinetics and dynamics of crystallization II: Makaopuhi lava lake. Contrib. Mineral. Petrol. (99) 292-305.

71. Cashman K.V. and ManganM.T. (1994) Physical aspects of magma degassing II: Constraints on vesiculation processes from textural studies of eruptive products. In: Volatiles in Magmas; M.R. Carrol and J.R. Holloway Eds, Reviews in Mineralogy (30) 447-478.

72. Cashman K.V. ,Mangan M.T. and S. Newman S. (1994) Surface degassing and modifications to vesicle size distributions in active basalt flows. J. Vol. Geoth. Res. (61) 4568.

73. Chiodini G., Marini L., Russo M. (2001) Geochemical evidence for the existence of high-temperature hydrothermal brines at Vesuvio volcano, Italy, Geochim. Cosmochim. Acta (65) 2129-2147.

74. Clocchiatti R, Moro A.D., Gioncada A., Jorono J.L., Mosbah M„ Pinarelli L. and Sbrana A. (1994) Assessment of a shallow magmatic system: the 1888-90 eruption, Vulcano Island, Italy. Bull. Volcanol. (56) 466-486.

75. Connolly J.A.D. and Podladchikov Yu.Yu. (2000) Temperature-dependent viscoelastic compaction and compartmentalization in sedimentary basins, Tectonophysics (324) 137-168.

76. Critescu N. and Hunsche U. (1991). A constitutive equation for salt. In Proceed. 7th International Congress Rock Mechanics, Aachen, Sept. 1991,1821-1830. Balkema, Rotterdam 3.

77. Cruden A.R., Koi H. and Schmeling H. (1995) Diapiric basal entainment of mafic into felsic magma EPSL (131) 321-340.

78. Davaille A. and Jaupart C. (1993) Transient high-Rayleigh number thermal convection with large viscosity variations. JFM (253) 141-166.

79. Doglioni, C., Innocenti, F., Mariotti, G. (2001) Why Mt Etna? Terra Nova. (13) 2531.

80. Dragoni M. and Magnanesi C .(\9%9) Displacements and stresses produced by pressurized, spherical magma chamber, surrounded by a viscoelastic shell. Phys.Earth and Planet. Inter. (56) 316-328.

81. Druzhinin 0. (1995) Dynamics of concentration and vorticity modification in a cellular flow laden with solid heavy particles. Physics of Fluid (7,#9) 2132-2142.

82. Duffleld W., Bacon C.R. and Dalrymple G.B. (1980) Late Cenozoic Volcanism, Geochronology, and structure of Coso Range, Inyo County, California, J Geophys. Res. (85, N.B5) 2381-2404.

83. Dunbar N. Jacobs G. and Naney M. (1995) Crystallization processes in artificial magma: variations in crystal shape, growth rate and composition with cooling history. Contrib. Mineral. Petrol. (120) 412-425.

84. Fenn P.M. (1977) The nucleation and growth of alkali feldspars from hydrous melts. The Canad. Mineral. (15)135-161.

85. Fletcher C.A.J. (1984) Computational Galerkin Methods. Springer Verlag, NY.

86. Fortin M and Fortin AA (1985) Generalisation of Uzawa's algorithm for the solution of the Navier-Stokes equations. Comm. Appl. Numer. Methods (1) 205-210.

87. Fournier R.O. (1999) Hydrothermal processes related to movement of fluid from plastic to brittle rocks in magmatic epithermal environment, Econ. Geol. (94) 1193-1211.

88. Francalanci L., Varecamp J.C., Vougioukalakis G., Defant M.J., Innocenti F., Manetti P. (1995) Crystal retention, fractionation and crustal assimilation in a convecting magma chamber, Nisyros Volcano, Greece. Bull. Volcanol. (56) 601-620.

89. Gaonac'h S„ Lovejoy S. and Strix J. (1996a) A scaling growth model for bubbles in basaltic lava flow. EPSL (139) 395-409.

90. Gaonac'h S„ Strix J. and Lovejoy S. (1996b) Scaling effects on the vesicle shape, size and heterogeneity of lavas from Mount Etna. Jour.Vol.Geotherm.Res. (74) 131-153.

91. Gaonac'h, H., Lovejoy, S., and Schertzer, D. (2003) Percolating magmas and explosive volcanism Geophys. Res. Lett., 30 (NO. 11), 1559, doi:10.1029/2002GL016022.

92. Gardner, J.E., Hilton, M., Carroll, M.R. (1999) Experimental constraints on degassing of magma: isothermal bubble growth during continuous decompression from high pressure, Earth Planet. Sci. Lett. (168) 201-218.

93. Ghiorso M. (1985) Chemical mass transfer in magmatic processes II: Application to equilibrium crystallization, fractionation and assimilation. Contrib. Mineral. Petrol. (90) 121-141.

94. Gibb F.G.F. Henderson C.M.B. (1992) Convection and crystal settling in sills. Mineralogy and Petrology (11)540-545.

95. Giggenbach W. F. (1997) The origin and evolution of fluids in magmatic-hydrothermal systems, in: H. L. Barnes (Ed.) Geochemistry of Hydrothermal Ore Deposits, John Wiley & Sons, New York, 3rd Edition, pp. 737-796.

96. Gonnermann H.M. and Manga M. (2003).Explosive volcanism may not be an inevitable consequence of magma fragmentation. Nature (426) 432-435.

97. Grout F.F. (1918J. Two phase convection in igneous magmas. J. Geology (26) 481499.

98. Grove T.L., Baker M.B. (1983) Effects of melt density on magma mixing in calc-alkaline series lavas. Nature (305) 416-418.

99. Guerin R.Z., Billia B. and Haldenwang P. (1991) Onset of solutal convection during directional solidification of a binary alloy. Phys. Fluids A (3,#8) 1873-1879.

100. Guglielminetti M. (1986) Mofete geothermal field, Geothermics (15) 781-790.

101. McGuinness M.J. (1996) Steady state solution selection and existence in geothermal heat pipes. I: the convective case. Int. J. Heat and Mass Transfer (39, 2) 259274.

102. Hammer J.E., and Rutherford, M.J. (2002) Kinetics of decompression-induced crystallization in Silicic melt: I. An experimental Study. J.G.R., (107,ECV8) 1-24.

103. Hardee H.C. (1982) Permeable convection above magma bodies, Tectonophysics (84) 179-195.

104. Harris A. and Stevenson D. (1997) Magma budgets and steady-state activity of Vulcano and Stromboli. Geophys. Res. Lett. (24,#9) 1043-1046.

105. Hauksson E., Hutton K. and Oppenheimer D. Earthquake Swarm at Coso Junction, Eastern California, July 2001; USGS Earthquake Hazards Program Northern California http://quake.wr.usgs.gov/eqinthenews/ci09674049.

106. Hedenquist J.W., Aoki M., Shinohava H. (1994) Flux of volatiles and ore-forming metals from magmatic hydrothermal system of Satsuma-Iwodjiva volcano. Geology (22) 585588.

107. Helz R.T. (1987) Differentiation behavior of Kilauea Iki lava lake, Kilauea Volcano, Hawaii: An overview of past and current work. Magmatic processes Physicochemical Principles. Geoch. Soc., Special Publ. No.l. Ed. B.O.Mysen, 241- 257.

108. Holz F., Behrens H„ Dingwell D.B. and Johanes W. (1995) H20 solubility inhaplogranitic melts: compositional, pressure and temperature dependence. Am Mineral (80) 94-108.

109. Hort M. (1997) Cooling and crystallization in sheet-like magma bodies revisited. EPSL (76) 297-317.

110. Huenges E., Erzinger J., KuckJ., Engeser B., Kessels W. (1997) The permeable crust: Geohydraulic properties down to 9101 m depth, J. Geophys. Res. (102) 18.25518.265.

111. Huges T JR, LiuWK, and Brooks AN. (1979) Finite element analysis of incompressible viscous flows by penalty function formulation. J.Comput.Phys. (39) 1-60.

112. Huppert HE. (1982). The propagation of two-dimensional and axisymmetric viscous gravity currents over a rigid horizontal surface. J.F.M. (121) 43-58.

113. Hurle D., Jakeman E. and Wheeler A. (1983) Hydrodynamic stability of the melt during solidification of a binary alloy. Phys.Fluids (26,#3) 624-626.

114. Hurwitz S. and Navon O. (1994) Bubble nucleation in rhyolitic melts: experiments at high pressure, temperature and water content. Earth and Planet. Sci. Lett. (122) 267-280.

115. Jambon A., Lussiez P., Clocchiatti R., Weisz J. and Hernandez J. (1992) Olivine growth rates in a tholeiitic basalt: An experimental study of melt inclusions in plagioclase. Chem. Geol. (96) 277-287.

116. Jarvis R„ Woods A. (1994) The nucleation, growth and settling of crystals from turbulently convecting fluid. J.Fluid Mech. (213) 83-107.

117. Jaupart C., Allegre C. (1991). Gas content,eruption rate and instabilities of eruption in silicic volcanoes. Earth Planet. Sci.Lett. (102) 413 429.

118. Jaupart C. and Tail S. (1995) Dynamics of differentiation in magma reservoirs. J.G.R. (100/B9) 17.615-17.636.

119. Johannes W, Holtz F. (1996) Petrogenesis and experimental petrology of granitic rocks. Springer Verlag, 335 pp.

120. Johansen S.T., Boysan F. andAyers W.H. (1987) Mathematical modelling of bubble driven flows in metallurgical processes. Appl.Sci.Res. (44) 197-207.

121. Karahaya K, Shinohara H., Sato G. (1994) Excessive degassing of Izu-Oshima volcano: magma convection in a conduit. Bull. Volcanol. (56) 207-216.

122. Karsten J.L., Holloway J.R. and Delaney JR. (1982). Ion microprobe studies of water in silicate melts: temperature dependent water diffusion in obsidian. Earth and Planet. Sci. Lett. (59) 240 -248.

123. Kerr R.C., Woods A. W., Worster M.G. and Huppert H.E. (1990) Solidification of alloy cooled from above Part 1. Equilibrium growth. J.F.M. (216) 323-342.

124. Kerr R., Woods A., Worster G. and Huppert H. (1990) Solidification of an alloy cooled from above Part 2. Non-equilibrium interfacial kinetics. J.F.M. (217) 331348.

125. Kerric D.M., Jacobs G.K (1981) Modified Redlich-Kwong equation for H20, C02 and H20-C02 mixtures at elevated pressures and temperatures. Amer. J. Sci. (281) 735-767.

126. Khodakovskii G., Robinovitcz M., Ceuleneer G., Trubitsyn V.P. (1995) Melt percolation in a partially molten mantle mush: Effect of a variable viscosity.EPSL (134)267-281.

127. King S.D., Raefsky A., Hager B.H. (1990) ConMan: vectorizing a finite element code for incompressible two-dimentional convection in the Earth mantle. Physics Earth Planetary Interiors (59)195-207.

128. Kirkpatrick R.J. (1976) Towards a kinetic model for the crystallization of magma bodies. J Geophys. Res. (81) 2565-2571.

129. Koyaguchi T., Hallwarth H. and Huppert H. (1993) An experimental study on the effects of phenocrysts on convection in magmas. J.Vol. Geotherm. Res. (55)15-32.

130. Kress V.C. and Carmichael I.S.E. (1988). Stoichiometry of the iron oxidation reaction in silicate melts. American. Mineral. (73) 1267-1274.

131. Kunda W., and Foots G. (1987) A simplified model of bubble driven flow in an axissymmetric container. Appl. Sci. Res. (44) 209-224.

132. Kuo, L.C. and Kirkpatrick R.J. (1982) Small angle X-Ray scattering study of pre-nucleation behavior oftitatnium free and titanium bearing diopside glasses. Amer. Mineral. (67) 676-685.

133. Lippmann M.J., Truesdell A.H. and Gutierrez P. (1997) What will a 6 km deep well at Cerro Prieto find, in: 21th Workshop on Geothermal Reservoir Engineering. Stanford University, Stanford, California, Proc., pp. 19-28.

134. Lister J.R. and Kerr R.C. (1989) The propagation of two-dimensional and axissymmetric viscous gravity currents at fluid interface. J. Fluid Mech. (203) 215 -249.

135. Lyakhovsky V., Hurwitz S. and Navon O. (1996) Bubble growth in rhyolitic melts: experimental and numerical investigation. Bul.Volcanol (58) 19-32.

136. Magnaudet J. and Eames I. (2000) The Motion of high-Reynolds-number bubbles in inhomogeneous liquid. Annu.Rev.Fluid Mech. (23) 659-708.

137. Manga M. (1996) Waves of bubbles in basaltic magmas and lavas. J.Geophys. Res. (101) 17.457-17.465.

138. Mangan M. (1990) Crystal size distribution systematics and the determination of magma storage times: The 1959 eruption of Kilauea volcano, Hawaii. J. Vole. Geotherm. Res. (44) 295-302.

139. Mangan M.T. and Cashman K. V. (1996) The structure of basaltic scoria and reticulite and inferences for vesiculation, foam formation and fragmentation in lava fontains. Journ.Vol.Geotherm.Res. (73) 1-18.

140. Mangan M. and Marsh B. (1992) Solidification front fractionation in phenocryst-free sheet-like magma bodies// J.Geology. (100) 605-620.

141. Mangan MT and Sisson T. (2000) Delayed, disequilibrium degassing in rhyolite magma: decompression experiments and implications for explosive volcanism. EPSL (183) 441-455.

142. Marsh D.B. (1988) Crystal size distribution (CSD) in rocks and the kinetics and dynamics of crystallisation. I. Theory. Contrib. Mineral. Petrol. (99) 277-291.

143. Marsh B.D. (1989) Magma Chambers. Annu. Rev. Earth Planet Sci. (17) 439474.

144. Marsh B. and Maxey M.R. (1985) On the distribution of crystals in convecting magma. J.Volc.Geotherm.Res. (24) 95-150.

145. Martel C., Bureau, H. (2001) In situ high-pressure and high-temperature bubble growth in silicic melts. Earth Planet. Sci. Lett. (191) 115-127.

146. Martin D. and Nokes R.( 1989) A fluid-dynamical study of crystal settling in convecting magmas. J.Petrol. (30)1471-1500.

147. Maxey M.R. (1987) The motion of small spherical particles in a cellular flow field. Phys.Fluids (30,#7)1915-1928.

148. McKenzie D.P. (1987) The compaction of igneous and sedimentary rocks Journal of the Geological Society (144, №2) 299-307.

149. Metrich N. Clocchiatti R., Mosbah M., Chaussidon M. (1993) The 1989-1990 activity of Etna: magma mingling and ascent of a H2O-CI-S rich basaltic magma. Evidence from melt inclusions. J. Volcanol. Geotherm. Res. (59) 131-144.

150. Metrich N. and Rutherford M.J.(1998) Low pressure crystallization paths of H20-saturated basaltic-hawaiitic melts from Mt Etna: Implications for open-system degassing of basaltic volcanoes. Geochim. Cosmochim. Acta. (62) 1195-1205.

151. Michaelides E.E. (1989) The role of vapor in volcanic activity. J. Volcanol. Geotherm. Res. (37) 251-260.

152. Moore J.N., Adams M.C. and Anderson A.J. (2000b) The fluid-Inclusion and Mineralogical Record of the Transition from Liquid- to Vapor- Dominated conditions in The Geysers Geothermal system, California. Economic Geology (95) 1719-1737.

153. Morse S. (1986) Thermal structure of crystallizing magma with two- phase convection. Geol. Mag. (123,#3) 205-214.

154. Mossop A.M., Murray D., Owen S. andSegall P. (1997) Subsidence at the geysers geothermal field: Results and simple models, in: 21th Workshop on Geothermal Reservoir Engineering. Stanford University, Stanford, California, Proc., 1997, pp. 377-382.

155. Mourtada-Bonnefoi C.C. and Laporte D. (2004) Kinetics of bubble nucleation in a rhyolitic melt: an experimental study of the effect of ascent rate. Earth Planet. Sci. Lett. (218) 521-537.

156. Muller, E., Heide K. andZanotto, E.D. (1993): Molecular structure and nucleation in silicate glasses. J.Non Cryst.Solids (155) 56-66.

157. Mungall J.E., Bogdassarov N., Romano C., Dingwell D.B. (1996) Numerical modeling of stress generation and microfracturing of vesicle walls in glassy rocks. J. Volcan. Geotherm. Res. (73) 33-46.

158. Navon 0., Chekhmir A., Lyakhovsky V. (1998) Bubble growth in highly viscous melts: theory, experiments, and autoexplosivity of dome lavas. EPSL(160) 763-776.

159. Ni J. and Beckermann C. (1991) A volume-averaged two-phase model for transport phenomena during solidification. Metallurgical Transacitons B. (22B) 349-361.

160. Nield. D.A (1967) The thermohaline Rayleigh-Jeffreys problem. J.Fluid.Mech. (29) 545-558.

161. Nigro N., Stori M. andIdelsohn S. (1994) Two-phase flow modelling in gas-stirred liquid vessels with SUPG-stabilized equal-order interpolation. Internat. J. for Numer. Meth. in Fluids (19)1-22.

162. Oldenburg, C.M. & Spera, F. J. (1991). Numerical modeling of solidification and convection in a viscous pure binary eutectic system. International Journal of Heat and Mass Transfer (34) 2107-2121.

163. Pareschi M.T., Pompilio M. and Innocenti F. (1990) Automated evaluation of volumetric grain-size distribution density from thin-section images. Computer & Geosciences (16) 1067-1084.

164. Park Y. and Hanson 5.(1999) Experimental investigation of Ostwald-ripening rates of forsterite in the haplobasaltic system. J.Volcanol. Geotherm. Res. (90) 103-113.

165. Pearlstein A.J., Harris R.M., Terrones G. (1989) The onset of convective instability in a triply diffusive fluid layer. J.Fluid Mech (202) 443-465.

166. Pelletier D., Fortin A. and Camarero R. (1989) Are FEM solutions of incompressible flows really incompressible ? (Or how simple flows can cause headaches!). I.J.for Numer. Methods in Fluids. (9) 99-112.

167. Proussevitch A.A. ,Sahagiart-D.L. and Anderson A.T. (1993a) Dynamics of diffusive bubble growth in magmas: isothermal case. J.Geophys. Res. (98) 22.283-22.307.

168. Proussevitch A.A., Sahagian D.L. and Kutolin V. (1993b) Stability of foams in silicate melts. J.Volcan.Geotherm.Res. (59) 161-178.

169. Putirka K. Clinoperoxene + liquid equilibria to 100 kbar and 2450 K (1999) Contrib. Mineral. Petrol. (135) 151-163.

170. Rudman M. (1992) Two-phase natural convection: implications for crystal settling in magma chambers. Phys. Earth Planet. Inter. (72) 153-172.

171. Saxena S.K., Fei Y. (1987) Fluids at crustal pressures and temperatures. I. Pure species. Contrib. Mineral. Petrol. (95) 370-375.

172. Schmeling //.(1987) On the relation between initial conditions and late stages of Rayleigh Taylor instabilities. Tectonophysics. (133) 65-80.

173. Schiano, P., Clocchiatti, R., Ottolini, L., Busa, L. (2001) Transition of Mount Etna lavas from from a mantle-plume to an island-arc magmatic source. Nature (412) 900-904.

174. Shinohara H., Karahaya K., Lowenstern J.B. (1995) Volatile transport in a convective magma column: impliation for porphyry Mo mineralization. Geology. (23, #12) 1091-1094.

175. Simakin A., Trubitsyn V., Schmeling H. (1994) Structure of the upper boundary layer of a solidifying intrusion with crystal sedimentation. Earth Planet Sci. Lett.(126) 333349.

176. Simakin A.G.,Schmeling H. and Trubitsyn V.P. (1997) Convection in melts due to sedimentative crystal flux from above. Phys. Eeath. Planet. Inter. (102) 185-200.

177. Simakin A.G., Armienti P., Epel'baum M.B. (1999). Coupled degassing and crystallization: experimental study in continuous pressure drop, with application to volcanic bombs. Bull. Volcanol.(61)275-287.

178. Simakin A. and Botcharnikov R. (2001) Degassing of stratified magma by compositional convection. J.Volcanol.Geotherm.Res.(l05)207-224.

179. Simakin A., Talbot C. (2001) Tectonic pumping of pervasive granitic melts. Tectonophysics (332) 387-402.

180. Simakin A. and Petford N. (2003) Melt redistribution during the bending of a porous, partially melted layer. Melt redistribution during folding of a partially melted layer. Geophys. Res Letts. (30,11) 1564

181. Simakin A., Talbot C. (2001) Transfer of melt between microscopic pores and macroscopic veins in migmatites. Phys. Chem. Earth (A) (26) 363-367.

182. SimakinA. and Ghassemi A.(2005) Modelling deformation of partially melted rock using a poroviscoelastic rheology with dynamic power law viscosity. Tectonophysics (397) 195-209.

183. Simakin A. and Ghassemi A. (2003) Salt loaded heat pipes: steady-state operation and related heat and mass transport. EPSL (215)411 -424.

184. Simakin A. and Ghassemi A. (2002) Heat Pipe in Porous Brine Saturated Rock: Quasi-Steady-State Operation and Related Heat and Mass Transport. 27th Stanford Geothermal Workshop 28-30 January 2002.

185. Sisson T.W. and Bacon C.R. (1999) Gas-driven filter pressing in magmas. Geology, 27,613-616.

186. Smith M. (1988) Thermal convection during the directional solidification of a pure liquid with variable viscosity. J. Fluid Mech. (188) 547-570.

187. Solomatov V. and Stevenson D. (1993) Kinetics of crystal growth in a terrestrial magma ocean. J.G.R. (98,# E3) 5407-5418.

188. Sparks R.S.J. (2003) Forecasting volcanic eruptions. Earth Planet. Sci. Lett. (210)1.15.

189. Sparks S., Huppert H., Koyaguchi T. and Hallworth M. (1993) Origin of modal and rhythmic igneous layering by sedimentation in a convecting magma chamber. Nature (361) 246-249.

190. Spera F, Yuen D., Clark S., Hong H.J. (1986) Double-diffusive convection in magma chambers: single or multiple layering. Geophys.Res.Letters (13, #1)153-156.

191. Spera F, Oldenburg C.M., Yuen D.A.(\9%9) Magma zonation: effects of chemical buoyancy and diffusion. Geophys.Res.Letters (16, #12)1387-1390.

192. Sunagawa, 1. (1984) Growth of crystals in Nature. Materials of the Erth's Interior. Tokyo. TERRAPUB, 63-105.

193. Symonds R.B., Rose W.I., Bluth G.J., Gerlah T. (1994) Volcanic-gas studies: methods, results, and application. In Volatiles in magmas, Eds. M.R. Carrol and Holloway John R.// Review in Mineralogy. (30) 1-66.

194. Swanson S„ E. (1977)' Relation of nucleation and crystal-growth rate to the development of granitic textures. American Mineral. (62) 966-978.

195. Tail S„ Jaupart C. and Vergniolle S. (1989). Pressure, gas content and eruption periodicity of a shallow, crystallising magma chamber. EPSL (92) 107-123.

196. Talbot C.J. (1998) Extrusions of Hormuz salt in Iran. In Blundell D. and Scott A. editors, Lyell: the Past is the Key to the Present. Geological Society of London, Special Publications, (143) 315-334.

197. Talbot C.J. and Jarvis R.J. (1984) Age, budget and dynamics of an active salt extrusion in Iran. Journal of Structural Geology (6)521-533.

198. Tanger J.C., and Pitzer K.S.(\989) Thermodynamics of NaCl-H20: A new equation of state for near-critical region and comparisons with other equations for adjoining regions, Geochim. Cosmochim. Acta (53) 973-987.

199. Taylor JR., Wall V.J., Powneby M.l (1992) The calibration and application of accurate redox sensors. American Mineral. (77) 284-295.

200. Thomas N. Tait S. and Koyaguchi T. (1993) Mixing of stratified liquids by the motion of gas bubbles: application to magma mixing. EPSL (115) 161-175.

201. Thomas N. Jaupart C. and Vergniolle S. (1994). On the vesicularity of pumice. J.Geophys. Res. (99, #B8) 15633-15644.

202. Tkachenko S.l. and Shmulovich K.I. (1992) Liquid-Vapor equilibrium at 400 to 600°C in aqueous system containing NaCl, KC1, CaC^ and MgCl^ Translated from: Dokladi Rossiyskoy Akademii Nauk, 326: No.6, 1055-1059.

203. Tomiya A. and Takahashi E. (1995) Reconstruction of an Evolving Magma Chamber beneath Usu Volcano since the 1663 eruption. J.Petrology. (36) 617-636.

204. Tonarini S., Armienti P., D'Orazio M„ Innocenti F., Pompilio M., Petrini R. (1995) Geochemical and Isotopic Monitoring of the Mt Etna 1989-1993 Eruptive Activity: Bearing on the Shallow Feeding System. J.Volcanol. Geotherm. Res. (64) 95-115.

205. Tonarini S., Armienti P., D'Orazio M. and Innocenti F. (2001) Subduction-like .fluids in the genesis of Mt.Etna magmas:evidence from boron isotopes and fluid mobile elements. Earth Planet. Sci.Lett. (192) 471 483.

206. Toramaru A. (1989) Vesiculation process and bubble size distribution in ascending magma with constant velocity. J. Geophys. Res. (94) 17523-17542.

207. Toramaru A. (1995) Numerical study of nucleation and growth of bubbles in viscous magmas. J. Geophys. Res. (100,B2) 1913-1931.

208. Toramaru A., Ishiwatarai A., Matsuzawa M., Nakamura N. Arai S. (1996) Vesicle layering in solidified intrusive magma bodies: a newly recognized type of igneous structure. Bull Volcanol. (58) 393-400.

209. Trigila R., Spera F. and Aurisiccho C. (1990) The Mount Etna eruption: thermochemical and dynamical inferences. Contrib. Mineral. Petrol. (104) 594-608.

210. Tuncay K., and Ortoleva P. (2001) Salt tectonics as a self-organizing process: a three dimensional reaction, transport and mechanical model (preprint).

211. Turcotte D.L. (1992) Fractal and chaos in geology and geophysics. Cambridge Univ. Press,Cambridge, pp 221.

212. Turcotte L. and Schubert G. (1982) Geodynamics. Application of continuum physics to geological problems. John Wiley and Sons, NY.

213. Unruh J. and Monastero F. (2001) in: 26th Workshop on Geothermal New seismic imaging of the Coso geothermal field, Eastern California, Reservoir Engineering, .Stanford University, Stanford, California, Proc., 2001.

214. Vergniolle S. and Jaupart C. (1986) Separated two-phase flow and basaltic eruption. J.Geophys.Res. (91, B12) 12.842-12.860.

215. Wager L.R. and Brown G.M. (1968) Layered igneous rocks. Oliver and Boyd, Edinborough,.588p.

216. Walker D., Jurewicz S. and Watson E.B. (1988) Adcumulate dunite growth in a laboratory thermal gradient. Contrib. Mineral. Petrol. (99) 306-319.

217. Wallance P., Anderson A.T., Davis A. (1995) Quantification of pre-eruptive exsolved gas contents in silicic magmas. Nature. (377) 612-616.

218. Walzer U. and Hendel R. (1999) A new conveciton-fractionaiton model for evolution of the principal geochemical reservoirs of the Earth's mantle Phys. Earth, Planet. Inter. (112) 211-256.

219. WholetzK. and Heiken G. (1992) Volcanology and geothermal energy. University of California Press, Berkeley, pp 432.

220. Wicks C., Thatcher W., Monastero F. and Hasting M. (2001) Steady-State Deformation of the Coso Range, East-Central California, Inferred from Satellite Radar Interferometry. J. Geophys. Res. (106) 13769 13780.

221. Woods A.W., Koyaguchi T. (1994). Transitions between explosive and effusive volcanic eruptions. Nature (370) 641 644.

222. Worster G. and Huppert H. (1993) The crystallization of lava lakes. J.Geophys. Res. (98,B9) 15891-15900.

223. Zanotto E.D. and Gahardi A. (1988) Experimental test of the general theory of transformation kinetics: Homogeneous nucleation in a Na20 2CaO- 3SiC>2 glass. J. Non-Crystal. Solids. (104) 73-80.

224. Zavelsky, V.O. and Salova T.P. (2002) Dissolution mechanism and speciation of water in glasses of the system Na20-Al203-Si02: An NMR study. Geochem. Intern. (40, suppl.l) 34-45.

225. Zieg M.J., Marsh B.D. (2002) Crystal size distribution and scaling laws in the quantification of the igneous textures. J.Petrology. (43, N1) 85-101.

226. Zienkiewlcz OC and Morgan K. (1983) Finite elements and approximation. John Wiley & Soms, NY.