ФИЗИКА ЗЕМЛИ, 2007, № 9, с. 79-86
УДК 550.311
МНОГОТЕМПЕРАТУРНАЯ МОДЕЛЬ КОМПАКЦИИ МАГМАТИЧЕСКОГО РАСПЛАВА В АСТЕНОСФЕРЕ (ЧИСЛЕННЫЙ ПОДХОД)
© 2007 г. В. В. Пак
Тихоокеанский океанологический институт ДВО РАН им. В.И. Ильичева, г. Владивосток
Поступила в редакцию 21.08.2006 г.
Разработана численная модель компакции флюида в вязком скелете с учетом фазового перехода. Температуры фаз различны. Решение находится методом асимптотического разложения относительно несжимаемого варианта, позволяющим устранить ряд вычислительных проблем, связанных со слабой сжимаемостью скелета. Для каждого приближения задача решалась методом конечных элементов. Исследовался процесс двумерной компакции магматического расплава в астеносфере под разломной зоной для одно- и двухтемпературного вариантов. Из-за пониженного порового давления в этой области концентрируется магматический поток. Магма с более высокой температурой, поступающая с нижних горизонтов, вызывает локальное нагревание скелета и интенсивное плавление его легкоплавкой компоненты. В двухтемпературной модели под разломной зоной образуется концентрационная аномалия магмы.
Показаны принципиальные ограничения, существенно усложняющие получение аналогичных результатов в рамках однотемпературной модели и обоснована необходимость применения многотемпературного варианта.
РАС8: 91.40.Jk
ВВЕДЕНИЕ
В настоящее время для описания процесса миграции магмы в частично подплавленной астеносфере широко используется модель пористой флюидонасыщенной вязкодеформируемой среды, позволяющая описать не только фильтрацию жидкой фазы, но и вязкие деформации скелета, фазовый переход и другие возможные процессы взаимодействия фаз [Каракин, 1988; McKenzie, 1984; Ribe, 1985].
Большой интерес представляет исследование этого процесса под разломными зонами. Т.к. наличие разлома в литосфере вызывает падение порового давления, магма должна устремляться в эту область, создавая концентрационную аномалию [Каракин, 1988]. По мнению многих исследователей, таким путем происходит концентрация магмы и тепла в узких рифтовых зонах океана с образованием мощных источников базальтового вулканизма [Каракин, 1988; Теркот и др., 1985]. Аналогичный механизм концентрации расплава может приводить к внутриплитовому вулканизму и вулканизму зон перехода океан - континент [Безверхний, Пак, 2003; Ефимов и др., 1986].
Выбор факторов изучаемого процесса, которые необходимо включить в модель, имеет принципиальное значение. Т.к. в рассматриваемом процессе существенную роль играет горизонтальное перераспределение магмы, необходимо использовать двумерные модели [Spiegelman, 1996]. Одно-
мерные модели, для которых возможно аналитическое исследование [Каракин, 1988; McKenzie, 1984; Ribe, 1985], здесь неприменимы.
Большинство исследований было проведено для изотермических моделей, применимых только при быстрых движениях расплава [Barcilon et al., 1986]. Но, если исследуется частичное плавление астеносферы с последующей миграцией магмы, то необходимы модели с учетом температуры фаз и фазового перехода, реализация которых имеет значительные трудности [Каракин, 1988; Трубицын и др., 2003].
Во всех разработанных неизотермических моделях компакции магмы предполагалось, что температуры скелета и флюида одинаковы. Равенство температур скелета и флюида означает, что относительная скорость движения флюида пренебрежимо мала по сравнению со скоростью межфазного теплообмена скелет-флюид. Это далеко не очевидно для всей мантии, т.к. собственная скорость флюида значительно больше, чем у скелета, а теплообмен происходит только через межфазную границу контакта с мантийной породой, теплопроводность которой чрезвычайно низкая [Теркот и др., 1985]. В этом случае температуры флюида и мантии могут быть различными.
Хотя возможность построения многотемпературных моделей теоретически обоснована [Нигма-туллин, 1987], они еще практически не использовались. Основная причина заключается в значитель-
ном усложнении модели: уравнение теплового баланса должно определяться отдельно для каждой фазы, кроме этого, в них должны быть введены дополнительные члены, характеризующие межфазный теплообмен. В результате чего, в модели появляется ряд дополнительных параметров (например, площадь межфазной поверхности), определение которых весьма затруднительно. Все это создает значительные трудности для решения задачи.
В настоящей работе построена численная модель многофазной двухтемпературной среды. Проведено численное моделирование компакции магмы под разломной зоной для одно- и двухтемпературного вариантов. Проведен сравнительный анализ полученных результатов. Показаны принципиальные ограничения однотемпературной модели и необходимость применения многотемпературного варианта.
СИСТЕМА УРАВНЕНИЙ.
НАЧАЛЬНЫЕ И КРАЕВЫЕ УСЛОВИЯ
Рассмотрим двумерную область, заполненную пористой флюидонасыщенной вязкодеформиру-емой средой, находящейся в поле силы тяжести. Предположим, что твердая фаза разделена на две компоненты ("ультраосновную" и "базальтовую") с плотностями ри, рь, скорости которых равны.
Так как на геологических временах инерционные члены в уравнениях движения скелета пренебрежимо малы, движение скелета описывается уравнениями ползущего течения (уравнениями Стокса):
П1(и,, + и,,г - иккЪ1])} + П2иккЪ1] -
(1)
- Рг - £5г-2рк«к[ 1- «к( Тк - Т* )] = 0,
где Ъу - символ Кронекера, £ - ускорение силы тяжести, П1, П2 - коэффициенты сдвиговой и объемной вязкости скелета, иг - скорость скелета, р - давление, р1 = ри(1 - ф) + рьф, р2 - плотности, соответственно, твердой и жидкой фаз, ф - объемная концентрация базальта в мантии, т1, т2 - объемные концентрации скелета и флюида (т1 + т2 = 1), а1, а2 - коэффициенты объемного терморасширения скелета и флюида, Т1, Т2 - температуры скелета и
флюида (по Кельвину), Т* , Т* - температуры скелета и флюида в начальный момент (по Кельвину).
Движение жидкой компоненты (расплава) описывается следующими уравнениями:
Рг = - ОУг - Р2т2[ 1- а2 (Т2- Т*)] 8Ъа;
В = П/,
(2)
Аналогично работам [Теркот и др., 1985; Трубицын и др., 2003], использовалась обратная степенная зависимость от концентрации т2 для объемной вязкости и гидравлического сопротивления:
П 2 =
П1
В =
Во
(т2- то) (т2- то)
где т0 = 2%, I =1, к = 2.
Межфазное взаимодействие магматического расплава с мантийной породой реализовано в модели в упрощенном варианте:
• плавлению подвергается только базальтовая составляющая;
• сжимаемость частиц твердой и жидкой фаз пренебрежимо мала, поэтому плотности фаз во всей расчетной области задавались постоянными.
Кроме того, аналогично [Трубицын и др., 2003; ШЬе, 1985], в уравнениях сохранения масс можно пренебречь разностью плотностей между базальтовой компонентой скелета и расплавом, т.е. считать, что фазовый переход (частичное плавление скелета и кристаллизация расплава) происходит без изменения объемов фаз.
Тогда уравнения массопереноса имеют следующий вид:
Э [ ( 1 - ф ) т 1 ] д?
Э(ф т1)
= -[(1 - ф)т1 ик]к; J,
(3)
т2
J,
■ЭГ = " (т2ик + ^кЬ +
где J12 - скорость фазового перехода между базальтовой компонентой и расплавом.
Используя тождество т1 + т2 = 1, сложим уравнения (3) и, после некоторых преобразований, получим следующую систему уравнений:
1 к, к
Ук, к = 0;
фкик ( 1 - т) р2 ;
Эф д ?
(4)
дт , 12
э? = "(тик + ук ),к + Р7;
где Уг - скорость фильтрации для расплава, к2 -проницаемость скелета, В - коэффициент гидравлического сопротивления, г|/ - вязкость флюида.
где т = т2.
Уравнения теплового баланса для каждой фазы должны содержать дополнительные члены, характеризующие межфазный энергообмен [Нигматул-лин, 1987]. Тепловые эффекты химических реакций, генерация тепла от радиогенных источников, а также межфазную работу сил вязкого трения, давления и т.д., в модели не учитывались.
2
Предполагалось, что компоненты скелета имеют одинаковую температуру.
Тогда тепловой баланс описывается следующими уравнениями:
д Т
Р1 Сх( 1- т)= р1 С1 (1- т)икТ 1;к + + (ХгТ 1 к )к + (1- т )т }кв ук + 621 ^ 12( ^21 С1 Т1 ) ;
дТ 2
Р2 С2 т "ду = Р2 С2 (тик + Ук) Т2, к +
+ (^Т2, к),к + 6 12+ ^ 12 ( и 12 С2 Т2 ) 5
(5)
где Сг- - удельные теплоемкости фаз, Хг - коэффициенты теплопроводности фаз, 6- - приток тепла через межфазную границу из у-той в г-тую фазу, и,г - удельная внутренняя энергия, переходящая из у-той В 7 -тую.
Кроме этого, должны выполняться следующие соотношения [6]:
^ 12¡12 612 + 621,
и,,- = и 5 = и (р, Т5(р)) = СгТ5(р)5
(6)
где ¡12 - скрытая теплота плавления базальта, Т^р) - температура солидуса "базальтовой" компоненты.
Для замыкания системы уравнений, необходимо определить функциональную зависимость для скорости фазового перехода J12. В работах [Каракин, 1988; Spiegelman, 1996] предлагается один из вариантов определения J12 через степень плавления Р, которая, при отсутствии массообмена, совпадает с пористостью, а в общем случае, зависит от динамики системы в целом и не может быть представлена в виде однозначной функции. В [Spiegelman, 1996] для описания декомпрессионного плавления мантии
дР
используется зависимость: J12 = р1и2 --— и предла-
д х2
гается упрощенная формула вычисления Р.
В предлагаемой модели использовались формулы из работы [Нигматуллин, 1987], применимые к системам с малыми неравновесностями. Эти формулы получены на основе принципа линейности, принципа симметрии (Кюри) и принципа взаимности (Онзагера):
612 = Кб ( Т1- Т2 ) + в 112 ( Т1- Т5) 5 (7)
Jl2 = в( Т1- Т2 ) + К^12 (Т1- Т8) . Подставляя в (5) соотношения (6) и (7), получим дТ
р1 С (1- т)= -Р1С1Т 1( 1- т)икТ 1, к + + (А,1 Т1,к),к + Кб(Т2 - Т1) + в112(Т2 + Т5 -2Т1) + + ^¡\2 (Т8 - Т1) + С1( Т8 - Т1) х
х[в(Т2-Т1) + КJ¡ 12(Т5- Т1)];
д Т 2
Р2 С2т "дТ = Р2С2( тик + Угк) Т2, к + (8)
+ (^Т2, к),к + Кб( Т1 - Т2 ) + в¡12 ( Т1 - Т5) + + С2( Т8 - Т2 )[в( Т1 - Т2) + КJ¡ 12 (Т1 - Т5)].
В итоге, общая система уравнений образуется из уравнений (1), (2), (4) и (8).
КРАЕВЫЕ УСЛОВИЯ
1. На верхней границе заданы следующие краевые условия:
• жесткое сцепление скелета с вышележащей толщей пород: игтг = 0, где тг - единичный вектор, касате
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.