Содержание к диссертации
Введение
1 Обзор литературы по гидродинамическим моделям лавовых течений 8
1.1 Введение 8
1.2 Изотермические модели 13
1.3 Неизотермические модели 22
1.4 Течения с фазовыми переходами 26
1.5 Выводы 31
2 Постановка задачи о неизотермическом течении лавы на осесимметричной искривленной поверхности 35
2.1 Уравнения движения 35
2.2 Система уравнений при малых углах наклона подстилающей поверхности 9 40
2.3 Система уравнений при конечных 9 42
3 Стационарные изотермические течения 45
3.1 Течения при малых 9 45
3.2 Течения при конечных 9 48
3.3 Развитие модели для описания вязкопластических течений . 49
3.4 Выводы 52
4 Автомодельное решение задачи о росте лавового купола на конической подстилающей поверхности 53
4.1 Автомодельное решение для конечных при степенном законе массоподвода 53
4.2 Обсуждение полученных решений 57
4.3 Автомодельное решение для конечных при экспоненциальном законе массоподвода 61
4.4 Неавтомодельные режимы течения при малых 9 63
4.5 Выводы 68
5 Трехмерные изотермические течения лавы на неосесимметричной конической поверхности 70
5.1 Введение 70
5.2 Автомодельные решения при малых 9'{ip) 72
5.3 Автомодельные решения при конечных 9'(ср) 77
5.4 Стационарные решения при конечных 9{ip) 82
5.5 Структура течения вблизи переднего фронта 83
5.6 Выводы 90
6 Неизотермические течения лавы на осесимметричной кони ческой поверхности 92
6.1 Система уравнений на масштабе жерла вулкана 92
6.2 Уравнения для течения на значительном удалении от жерла 94
6.3 Численный метод 94
6.4 Результаты расчетов 97
6.5 Выводы 102
7 Асимптотические модели лавовых течений с поверхностной солидификацией 103
7.1 Постановка задачи 104
7.2 Приближение тонкого слоя 107
7.3 Асимптотические решения при S^oo 110
7.4 Численные решения при S ~ 1 117
7.5 Обсуждение 119
7.6 Выводы 123
Заключение 125
Литература
- Неизотермические модели
- Система уравнений при малых углах наклона подстилающей поверхности 9
- Автомодельное решение для конечных при экспоненциальном законе массоподвода
- Автомодельные решения при конечных 9'(ср)
Введение к работе
При экструзивном вулканическом извержении расплавленные магматические породы под действием избыточного давления в очаге поступают из жерла вулкана, что приводит к формированию течения тонкого слоя остывающей лавы на подстилающей поверхности. Лава, как правило, представляет собой смесь силикатного расплава, кристаллов и газовых пузырьков, причем концентрация кристаллической и пузырьковой фракций может существенно изменяться в процессе развития течения. Аномально большая вязкость лавы и, как следствие, малые скорости течения, сильно нелинейная зависимость вязкости от температуры и усложненная реология лавы определяют специфические свойства лавового течения в условиях экструзивного вулканического извержения [1]. Поведение и структура лавовых течений, закон распространения переднего фронта, развитие неустойчи-востей существенно зависят от свойств и расхода лавы, поступающей из жерла вулкана, топографии подстилающей поверхности, а также свойств окружающей среды.
Величина потока тепла на свободной поверхности лавового течения, как правило, достаточно велика, что приводит к формированию тонкого слоя аморфного твердого вещества вблизи свободной границы. Более медленный процесс кристаллизации во всем слое движущейся лавы постепенно приводит к полной солидификации и остановке течения. Передний фронт течения в результате остывания и солидификации останавливается. В некоторых случаях это происходит даже раньше, чем прекращается извержение из жерла вулкана.
Интерес к исследованию лавовых течений обусловлен рядом причин. Среди них - необходимость оценить прямую опасность, которую представляет распространение лавы при экструзивных вулканических извержениях, а также угрозу периодических извержений на таких вулканах, как Этна на острове Сицилия и Килауэа на Гавайских островах. Более разрушительными являются пирокластические потоки, наблюдавшиеся, например, при извержениях на вулканах Мирапи (Индонезия, 1994 г.) и Унзен (Япония, 1991-1995 гг). Пирокластическое течение представляет собой поток частиц разрушенного лавового купола на крутом склоне. Лавовые течения могут также приводить к внезапным наводнениям, оползням и грязевым потокам.
Еще одна причина, побуждающая ученых к исследованию лавовых потоков, это необходимость интерпретации данных о древних извержениях на Земле и других планетах. Данные об особенностях свободной поверхности и законе распространения переднего фронта течения позволяют оценить интенсивность извержения и реологию лавы. Среди других целей моделирования лавовых течений отметим описание формирования рудных месторождений в результате древних высокотемпературных лавовых течений, так называемых "коматитов" (komatiites), а также возникновение таких необычных явлений, как очень большие риолитовые течения на поверхности Земли, огромные тонкие лавовые купола на Венере и базальтовые течения длиной до 100 км на поверхности Земли, Луны и Марса.
При изучении лавовых течений большую роль играют как упрощенные полуаналитические модели, так и модели, основанные на прямом численном моделировании. Упрощенные изотермические модели течений сильновязкой ньютоновской или вязкопластической жидкости были построены для описания медленных течений без учета остывания и солидификации. Эти модели послужили основой для обобщения и учета зависимости вязкости лавы от температуры и интенсивности солидификации. Результаты, полученные в рамках моделей с учетом проплавлення подстилающей поверхности, позволили объяснить такие явления, как вовлечение подсти-
лающих скальных пород в лавовое течение и отложение металлических руд. Результаты экспериментов с лабораторными аналогами расплавленной магмы удовлетворительно согласуются с решениями, полученными в рамках упрощенных теоретических моделей.
Известные в литературе упрощенные асимптотические модели, основанные на приближении тонкого слоя, учитывают такие специфические свойства лавы, как зависимость вязкости от температуры и скорости солидификации, величины и объемной доли кристаллов, наличие предела текучести, неустойчивость передней кромки течения и формирование так называемых "пальцев" на переднем фронте. Однако подавляющее большинство работ посвящено течениям на горизонтальной либо наклонной плоской подстилающей поверхности. Оставлены без внимания эффекты криволинейно-сти подстилающей поверхности, приводящие к существенному искривлению переднего фронта течения и возникновению локальных максимумов скорости движения переднего фронта. Учет влияния неплоской геометрии подстилающей поверхности и эффектов теплопереноса на распространение лавового течения является одной из основных целей данного исследования.
Цель работы
Целью настоящей работы является построение асимптотических моделей неизотермического течения тонкого слоя сильновязкой жидкости с массо-подводом на криволинейной подстилающей поверхности, исследование динамики нестационарного течения с учетом сильно нелинейной зависимости вязкости жидкости от температуры, перетекания жидкости в поперечном к потоку направлении и образования твердого приповерхностного слоя в результате остывания и солидификации. Полученные в рамках построенных моделей решения могут быть использованы для описания таких явлений, как лавовый поток на склонах вулканического конуса и рост лавового купола на искривленном дне кратера вулкана при экструзивных извержениях.
Неизотермические модели
Сравнительно меньшее число работ посвящено неизотермическим пленочным течениям сильновязкой тяжелой жидкости в приложении к описанию лавовых течений. В [46] представлена модель неизотермического тонкого слоя вязкой ньютоновской жидкости на горизонтальной плоскости в приложении к течениям лавы. Эта модель представляет собой обобщение уравнений [24] на случай неизотермического течения. На свободной и подстилающей поверхностях задано нулевое граничное условие для температуры. Таким образом, предполагается, что интенсивность потока тепла на верхней и нижней границах имеет одинаковый асимптотический порядок. Эта модель является упрощенной, так как известно, что интенсивность потока тепла на свободной поверхности за счет излучения существенно превышает поток тепла в подстилающую поверхность за счет теплопроводности [2]. В модели [46] предполагается, что вязкость обратно пропорциональна температуре. После осреднения по поперечной координате получена система, состоящая из параболического уравнения для формы свободной поверхности и гиперболического уравнения для осредненной температуры. Полученные уравнения решались численно для случаев растекания фиксированного объема жидкости и течения с постоянным массоподводом.
В [47] представлено развитие модели [46] на случай течения жидкости, у которой от температуры зависит не только вязкость, но и плотность. Получены численные решения для толщины слоя и закона распространения переднего края течения. Проведено сравнение с экспериментом [48] и обоснована возможность приложения полученных результатов к описанию лавовых течений.
В [49] решена задача о стационарном неизотермическом течении тонкой пленки вязкой жидкости с учетом термокапиллярных эффектов на подогретой наклонной плоскости. Получено численное решение для формы свободной поверхности. В окрестности контактной линии уравнения исследованы аналитически. Проведено исследование линейной устойчивости относительно основного стационарного решения. Получены диапазон неустойчивых волновых чисел и максимальный коэффициент нарастания в зависимости от влияния силы тяжести, термокапиллярных эффектов и скорости движения контактной линии. В [50] приведен подробный обзор работ по моделям жидких пленок на горизонтальной либо наклонной плоскости, учитывающим поверхностное натяжение, термокапиллярные эффекты, зависимость физических свойств жидкости от температуры, влияние тепло-и массоподвода, фазовые превращения, влияние Ван-дер-Ваальсовых сил, взаимодействие с поверхностно-активными веществами, а также влияние вращения подстилающей поверхности.
В [51] выведены асимптотические уравнения трехмерного неизотермиче ского тонкого слоя ньютоновской жидкости на горизонтальной плоскости. В этой работе, являющейся существенным развитием неизотермической модели [46], проведено детальное обобщение уравнений изотермического тонкого слоя на горизонтальной плоскости [24] на случай неизотермического течения. В граничных условиях на свободной и подстилающей поверхностях задан поток тепла как известная функция температуры границы. Проведен сравнительный асимптотический анализ уравнений в зависимости от модифицированного числа Пекле. Кратко рассмотрен случай течения с солидификацией вблизи подстилающей поверхности. Полученные системы уравнений исследованы численно. Указано, что найденные решения могут быть использованы в геофизических приложениях, а также для описания течений смазочно-охлаждающих жидкостей в конструкциях ядерных реакторов.
В [52] модель осесимметричного роста вязкопластического лавового купола на горизонтальной плоскости [37] обобщена на случай неизотермического течения. Аналогичные модели неизотермического тонкого слоя ньютоновской жидкости применялись в [53] при исследовании влияния термокапиллярных эффектов на распространение капли вязкой жидкости по твердой поверхности ив [51] при исследовании течений охлаждающих жидкостей в конструкциях ядерных реакторов. При построении модели использовано предположение, что перенос тепла поперек слоя происходит достаточно быстро, так что течение можно считать вертикально-изотермическим. Приближение вертикально-изотермического тонкого слоя использовано также, например, в [54] при изучении движений тонкой пленки жидкости с растворенными в ней поверхностно-активными веществами.
Система уравнений при малых углах наклона подстилающей поверхности 9
Рассматривается течение во внутренней области при малых углах наклона образующей поверхности к горизонтали на масштабе длины Li = D (в формулах (2.1.1)-(2.1.3) положим L = Li). Здесь и далее в случаях, когда это необходимо, используем индекс і для параметров внутренней области и индекс е для параметров внешней области. Главные члены асимптотических разложений искомых функций при малых 9 имеют порядки:
Остальные функции имеют порядок единицы. Индекс 0 ниже опущен. Полагая в качестве характерного масштаба компоненты скорости v величину Q/Lf, где Q - характерный масштаб расхода массоподвода, из уравнения неразрывности имеем Q/ Lf = є/ U. Последнее равенство позволяет выразить величины масштаба скорости U и малого параметра Е{ через определяющие размерные параметры: Считаем Re = о(єі ). Подставляя (2.2.1) в (2.1.1)-(2.1.3) и оставляя главные члены, имеем:
Отметим, что динамическое уравнение в проекции на продольное направление содержит проекцию силы тяжести, производную касательного напряжения и самоиндуцированный продольный градиент давления, который связан через граничное условие с неизвестной формой свободной поверхности. В рамках сделанных предположений полученная система уравнений описывает течения не только при малых, но и при конечных и умеренно больших числах Рейнольдса.
Используя (2.2.2), получаем выражения для параметров к\ и / через определяющие размерные параметры: В табл. 2 приведены порядки величин параметров fci, / и є І для двух типов лавы: кислой и основной [84]. При вычислении значений параметров использованы данные [52].
Рассмотрим течение во внешней области на поверхности конуса при конечных углах наклона образующей конуса к горизонтали на масштабе длины Le = Re (в формулах (2.1.1)-(2.1.3) положим L = Le). Главные члены асимптотических разложений искомых функций при конечных 9 имеют порядки:
Остальные функции имеют порядок единицы. Индекс 0 ниже опущен. Для масштаба скорости U и малого параметра ее имеем выражения через определяющие размерные параметры: X?\2/3_L Р9) (2.3.2)
Динамическое уравнение полученной системы содержит лишь проекцию силы тяжести и производную касательного напряжения и не содержит самоиндуцированного продольного градиента давления. Это является принципиальным отличием от системы уравнений для малых 9 (2.2.3). Аналогично (2.2.3), система (2.3.3) описывает течения не только при малых, но и при умеренно больших числах Рейнольдса.
Граничное условие при хе = 0 должно быть заменено условием асимптотического сращивания с решением во внутренней области [86]. Последнее условие состоит в том, что асимптотика решения во внешней области при хе — 0 должна совпадать с асимптотикой решения во внутренней области при ХІ — оо, переписанной во внешних переменных.
Используя (2.3.2), получаем выражения для параметров к\ и к Порядки величин параметров fci, / и ее приведены в табл. 2 (при вычислении значений параметров использованы данные [52]). Исходя из оценки порядка величины параметра fci, слагаемым с диссипативной функцией в уравнении притока тепла (2.3.3) пренебрегаем.
В рамках рассматриваемой модели все безразмерные параметры считаются известными, кроме безразмерной величины потока тепла Nu на свободной поверхности. В целях иллюстрации рассмотрим все асимптотические возможности.
Автомодельное решение для конечных при экспоненциальном законе массоподвода
Рассматривается стационарное изотермическое течение во внутренней области, описываемое уравнениями (2.2.3). Положим в уравнениях (2.2.3) Т = 1, /І = 1 и д/dt = 0. Интегрируя динамические уравнения системы (2.2.3) с учетом граничных условий, получаем выражения:
Подставляя полученное выражение для и в уравнение неразрывности в интегральной форме, получаем уравнение для формы свободной поверхности:
Рассмотрим решения уравнения (3.1.1) в случае, когда твердая поверхность есть часть сферы размерного радиуса Re = Ьгєг , тогда О(х) = х. Уравнение (3.1.1) имеет однопараметрическое семейство решений (свободный параметр /г(0)), что является проявлением эллиптического характера исходной системы уравнений Навье-Стокса. В работах, посвященных вдуву газа с поверхности тела, обтекаемого сверхзвуковым потоком (см., например, [87]), также рассматриваются стационарные уравнения тонкого слоя с самоиндуцированным продольным градиентом давления и возникают однопараметрические семейства решений. По аналогии с [87], для отбора единственного решения необходимо поставить дополнительное граничное условие при х — оо. При х — 0 решение (3.1.1) имеет следующую асимптотику: При х — оо существует единственное решение (3.1.1), удовлетворяющее условию h(x) Ои h (x) — 0: и однопараметрическое семейство решений, удовлетворяющих условию h(x) оо и h (x) — оо:
Возрастающие решения, имеющие асимптотику (3.1.3), при х — оо выходят на линии постоянного уровня Y = const в декартовой системе координат (X, Y) на масштабах порядка радиуса сферы. Убывающие решения, для которых h(xr) = 0 при конечном хг, приобретают физический смысл в рамках стационарной постановки при дополнительном условии в точке хг. Таковым может быть наличие точечного стока в точке хг. Вместе с тем, эти кривые могут быть интерпретированы как решения нестационарной задачи в квазистационарной постановке, в которой скорость массоподвода зависит от времени как от параметра. В этом случае отбор единственной кривой из семейства должен производиться из интегрального закона сохранения массы, то есть условия равенства объема движущейся жидкости в текущий момент времени суммарному объему жидкости, поступившей из канала за все время течения. Семейства возрастающих и убывающих решений разделены сепаратрисой с асимптотикой (3.1.2).
В качестве примера на Рис. 3.1 представлены семейства решений уравнения (2.3) в системах координат (ж, г]) (а) и (X, Y) (6) (є = 1/27), найденные численно методом Рунге-Кутта четвертого порядка при условии, что на выходе из канала задан профиль скорости массоподвода в форме Пуазейля vm(x) = 1-х2.
В случае, когда твердая поверхность есть конус, имеющий малый угол наклона образующей к горизонтали и затупленный по сфере радиуса Rei уравнение для формы свободной поверхности при х Хо принимает вид XQ 1 - координата точки, в которой сфера переходит в конус):
Как и в случае сферы, семейство интегральных кривых уравнения (3.1.4) имеет сепаратрису, разделяющую возрастающие и убывающие решения. Возрастающие решения и сепаратриса имеют следующие асимптотики:
Рассматривается стационарное изотермическое течение во внешней области на поверхности конуса при точечном массоподводе в вершине с безразмерным объемным расходом Q, описываемое уравнениями (2.3.3). Положим в уравнениях (2.3.3) Т = 1, /І = 1 и д/dt = 0. Интегрирование динамических уравнений системы (2.3.3) с учетом граничных условий приводит к формулам:
После подстановки в уравнение неразрывности в интегральной форме выражения для и, имеем единственное решение для формы свободной поверхности: Формула (3.2.1) справедлива в рамках гипотезы тонкого слоя при 9 1 и х 1. Асимптотика решения (3.2.1) при х — 0 в случае, когда поверхность конуса переходит в сферу, имеет вид:
Это решение можно срастить лишь с единственным решением во внутренней области (на масштабе порядка радиуса канала), убывающим на бесконечности и имеющим асимптотику (3.1.2) (2 на Рис. 3.1, а). Таким образом осуществляется отбор единственного решения из семейства решений во внутренней области для построения равномерно пригодного решения в обеих областях [86]. Решение для течения на конусе в размерных переменных [рдх sin29_ не зависит от характерного масштаба длины, выбранного при обезразме-ривании.
Обобщим полученные стационарные решения на случай течения вязкопла-стической жидкости с реологической моделью [23]. Это реологическое соотношение учитывает наличие предела текучести и эффекты увеличения либо уменьшения вязкости при сдвиге, характерные для лавы [37]. В рамках модели [23] справедливо следующее размерное соотношение между компонентами тензора вязких напряжений и тензора скоростей деформаций:
Здесь В - предел текучести, г и е - вторые инварианты тензоров вязких напряжений и скоростей деформации, тп - положительное действительное число. Введем обозначение у = Нв{х) для поверхности, на которой достигается предел текучести. Ниже этой поверхности имеет место течение, а над ней жидкость движется в "затвердевшем" состоянии с постоянной в вертикальном сечении продольной скоростью.
После интегрирования асимптотических уравнений тонкой пленки на поверхности сферы большого радиуса во внутренней области (на масштабе порядка радиуса канала) аналогично тому, как это сделано в разд. 3.1, получаем систему уравнений для неизвестных функций h и Кв (обозначения соответствуют введенным ранее): Система (3.3.1) имеет единственное решение, удовлетворяющее условию h(x) — 0, h (x) — 0, 1ів(х) — 0 при х — оо. Это решение имеет асимптотики:
Отметим, что асимптотика для формы свободной поверхности не зависит от предела текучести Вив случае течения бингамовской жидкости (т = 1) совпадает с соответствующей асимптотикой (3.1.2), полученной выше для случая ньютоновской жидкости.
Как и в случае течения ньютоновской жидкости, система (3.3.1) имеет семейство убывающих решений, которые в рамках модели стационарного течения приобретают физический смысл при условии точечного стока в точке, где толщина пленки обращается в ноль. Семейство этих решений можно также интерпретировать как решение нестационарной задачи в квазистационарной постановке. При этом отбор единственного решения из семейства, как и в случае ньютоновской жидкости, должен производится из интегрального закона сохранения массы. Существенным отличием от случая ньютоновской жидкости является отсутствие возрастающих решений (в системе координат, связанной со сферой).
Автомодельные решения при конечных 9'(ср)
Рассмотрим нестационарное изотермическое течение на произвольной поверхности с малым углом наклона образующей к горизонтали 9, описы ваемое системой (2.2.3). В уравнениях (2.2.3) положим Т = 1 и /І = 1. Предполагается, что форма свободной поверхности является гладкой на оси симметрии:
Интегрируя динамическое уравнение и уравнение неразрывности в (2.2.3) по поперечной координате, получаем эволюционное уравнение для формы свободной поверхности:
В случае течения на горизонтальной плоскости ( д(х) = 0) уравнение (4.4.2) переходит в уравнение (1.2.1), для которого в [24] были найдены автомодельные решения уравнения при степенной зависимости объема движущейся жидкости от времени. В случае д(х) ф 0, даже при отсутствии масштаба длины в граничных условиях, уравнение (4.4.2) не имеет автомодельных решений, однако, как показано ниже, его решения имеют автомодельные асимптотики. Аналогичная ситуация возникала в [33] при изучении течений на плоскости с малым углом наклона к горизонтали от точечного (линейного) источника.
В качестве примера рассмотрим течение на поверхности конуса ( д(х) = const). В этом случае уравнение (4.4.2) не имеет автомодельного решения при точечном массоподводе (L Li), поэтому оно решалось численно для течения с массоподводом из конечного отверстия (L = Li) на поверхности в форме конуса, затупленного по сфере радиуса є-1 3. Для этого случая имеем "д(х хо) = х и "д(х хо) = хо, где хо - координата точки перехода сферы в конус.
Уравнение (4.4.2) решалось численно конечно-разностным методом. Была введена прямоугольная сетка с шагом по времени Д и шагом по пространственной координате Ах и составлена двухслойная неявная схема на шеститочечном шаблоне с использованием центральных разностей для аппроксимации dh/dx. Свободный член в уравнении (4.4.2) терпит разрыв в точке Жо, поэтому была выбрана специальная сетка, в которой точка разрыва совпадает с узлом сетки. При использовании неявной схемы дополнительного условия в этой точке не требуется [92]. Граничным условием при х = 0 является (4.4.1), а условие h = 0 при х — оо сносится на некоторое достаточно большое X, до которого передний фронт течения заведомо не дойдет за расчетное время. Схема имеет погрешность аппроксимации О (At, Ах), так как левое граничное условие аппроксимируется с первым порядком точности. Для того, чтобы тип уравнения не менялся в пределах расчетной области, при t = 0 задана ненулевая толщина слоя на всем отрезке [0,Х] по формуле h(0,x) = 10 3ехр(—ж2/25). Система линейных уравнений с трехдиагональнои матрицей, возникающая на каждом шаге по времени, решалась методом прогонки [93]. Для матрицы схемы выполняется условие преобладания диагональных элементов, поэтому метод прогонки устойчив [92]. Приближение тонкого слоя, в рамках которого было получено уравнение (4.4.2), не выполняется в области бесконечных градиентов h в малой окрестности передней границы течения. По аналогии с [37], для того, чтобы избежать разрыва решения на передней границе, в схему была введена искусственная вязкость в виде малой аддитивной постоянной в коэффициенте при старшей производной. В процессе расчетов был найден диапазон значений малой постоянной, в котором имеет место слабая зависимость решения от величины этой константы.
Результаты расчетов в случае, когда массоподвод непрерывно выходит на стационарный режим по закону проиллюстрированы на Рис. 4.3 (О (ж) - функция Хевисайда). Численные расчеты проведены для XQ = 1.2. При t — оо решение совпадает с решением соответствующей задачи в квазистационарной постановке, полученным в разд. 3.1. Найденные решения описывают распространение лавы по подстилающей поверхности при повторном экструзивном извержении вулкана.
Решение уравнения (4.4.2) в случае, когда объем движущейся жидкости растет со временем пропорционально Р , имеет автомодельную асимптотику при оо. Движение в автомодельном режиме определяется балансом сил тяжести и сил вязкого трения (продольный градиент давления со временем становится внепорядковым членом в динамическом уравнении (2.2.3)), и асимптотика совпадает с автомодельным решением задачи при конечных углах наклона образующей к горизонтали, полученным в разд. 4.1. Для этого решения закон движения переднего фронта течения (4.2.3) в безразмерной форме имеет вид (С\ - константа, определяемая граничными условиями
