CC BY
73
Y Жукова Е.В., Занкович А.М., Коваленко И.Г., Фирсов К.М., 2012
УДК 524.575 + 535.231.11 ББК 22.66
ГИДРОСТАТИЧЕСКАЯ МОДЕЛЬ САМОГРАВИТИРУЮЩЕГО ОПТИЧЕСКИ ПЛОТНОГО МЕЖЗВЕЗДНОГО ОБЛАКА 1
Е.В. Жукова, А.М. Занкович, И.Г. Коваленко, К.М. Фирсов
Рассмотрена гидростатическая модель оптически плотного газопылевого межзвездного облака, находящегося в равновесии с собственной тяжестью и внешним излучением. Проведены численные расчеты стационарных распределений концентрации и температуры газа в облаке. Выявлены условия для возникновения конвективной неустойчивости.
Ключевые слова: межзвездная пыль, перенос излучения, газопылевые облака, взаимодействие излучения с веществом, конвективная неустойчивость, приближение Эддингтона.
Введение
Наблюдения свидетельствуют о существенной турбулизации газа в межзвездных облаках как атомарного, так и молекулярного водорода со скоростью турбулентных движений, сопоставимой со звуковой [17; 20; 24].
Причиной турбулентности могут служить конвективные потоки, генерируемые внутри облака. Вследствие неоднородного прогрева и неоднородного выхолаживания вещества внутри оптически плотной среды тепловое равновесие среды может стать неустойчивым по отношению к малым возмущениям, либо вовсе перестать существовать как равновесие гидростатическое. Все это должно приводить к возникновению течений, стремящихся перестроить среду таким образом, чтобы в ней вновь установилось тепловое равновесие в среднем, но не статическое, а динамическое.
Турбулентность в межзвездных облаках служит фактором, препятствующим звездообразованию. Однако остановить гравитационное сжатие облака может только такая турбулентность, энергия которой превышает тепловую энергию вещества облака, а это значит, что турбулентность должна быть около- или сверхзвуковой. При этом, если нет механизмов поддержания турбулентности в облаке, оно не сможет долго сопротивляться сжатию. Действительно, в сверхзвуковой турбулентности будут рождаться ударные фронты, на которых будут происходить перекачка кинетической энергии турбулентных движений в тепловую энергию облака, временный разогрев вещества облака и последующее охлаждение через высвечивание. Характерные времена охлаждения обычно много меньше времени свободного сжатия облака, поэтому облако недолго будет находиться в стадии турбулентного динамического равновесия [3].
В качестве механизма поддержания принудительной турбулентности в облаке, могут выступать эффекты взаимодействия вещества облака с излучением, а именно: (1) приток энергии за счет прогрева газа излучением от внешних источников, например, от первых молодых звезд, родившихся в облаке, и (2) потери энергии на переизлучение разогретого газа в том диапазоне длин волн, в котором межзвездный газ в облаке прозрачен.
Взаимодействие газа облака с излучением при определенных условиях может создавать условия для развития постоянно действующей конвекции. Условия для развития конвекции возникнут, если газ в центре облака будет горячее газа на периферии облака.
В случае оптически тонкого облака условия для развития конвекции не возникают. Взаимодействие с проникающим в облако разогревающим излучением есть унарный процесс, в котором отдельный фотон взаимодействует с отдельной частицей вещества. Интенсивность такого процесса прямо пропорциональна концентрации частиц. В то же время охлаждение, как правило, есть бинарный процесс, в котором в элементарном акте взаимодействия участвуют две частицы, при этом рождается уносящий энергию фотон. Примерами бинарных процессов служат рекомбинация или тормозное излучение. Комбинированным, но в конечном итоге бинарным, является случай охлаждения газа на примесных пылевых частицах. Тепловое излучение пылевых частиц есть унарный процесс, но столк-новительное взаимодействие остывшей вследствие потерь на тепловое излучение пылевой частицы с частицей газа с последующим охлаждением газа является процессом бинарным [1].
Поскольку интенсивность охлаждения для бинарных процессов пропорциональная квадрату концентрации частиц вещества, охлаждение доминирует над нагревом в областях более плотной среды. Это означает, что в центре облака, сжатого под действием собственной тяжести, доминирует остывание, и вещество оказывается холоднее, чем на периферии облака. В «спокойных» системах (нетурбулизованные облака или газовые диски нормальных галактик) именно такое конвективно-устойчивое распределение и наблюдается.
Однако, если ядро облака содержит источники излучения, а вещество облака является оптически плотным относительно прогревающего излучения, то реализуется инверсное, конвективно-неустойчивое распределение температуры и плотности.
Разогревающим для газа и пыли выступает ультрафиолетовое (У Ф)-излучение молодых звезд О- и В-класса, относительно которого вещество облака может быть оптически плотным. В то же время выхолаживание облака наиболее эффективно происходит в инфракрасном (ИК)-диапазоне, для которого вещество облака почти прозрачно. К моменту достижения периферийных областей прогревающее УФ-излучение теряет значительную часть своей интенсивности (потери происходят вследствие перекачки энергии УФ-излучения в ИК-диапазон с последующим свободным уходом фотонов из системы), поэтому периферийные области оказываются холоднее внутренних.
Идея сверхзвуковой турбулентности высказывалась еще Вайцзеккером в 40-е гг. в связи с созданием космогонии Вселенной [26; 27], для межзвездных облаков сверхзвуковая турбулентность обсуждалась начиная с 80-х гг. [4; 9; 10; 18]. Наблюдения показывают, что сверхзвуковая турбулентность есть реальный феномен [19].
На стыке 90-х и 2000-х гг. появились первые динамические модели сверхзвуковой турбулентности в межзвездной среде [17; 20; 24]. Однако во всех этих работах по численному моделированию турбулентность задается «вручную» - через действие псевдослучайных сил, генерируемых внутри расчетной области по заданному самим вычислителем закону. Динамических моделей сверхзвуковой турбулентности в облаках с учетом переноса излучения, в которых турбулентность возникала бы «самостоятельно» вследствие развития тех или иных физических неустойчивостей, на данный момент не существует. Это связано с тем, что, во-первых, турбулентность есть существенно трехмерный феномен и численно его моделировать стало по-настоящему возможным только в самое последнее время с появлением суперкомпьютеров с большой оперативной памятью и высокой производительностью. Во-вторых, расчет процесса неодномерного (в конфигурационном пространстве) переноса излучения есть сверхгромоздкая задача, требующая на самом деле решения существенно более многомерной задачи, поскольку в уравнении переноса присутствуют дополнительные размерности - направления распространения излучения и частоты.
В астрофизике к построению многомерных динамических моделей с учетом двумерного переноса излучения перешли только в последние несколько лет. В данной работе мы подготавливаем основу для создания такой оптико-гидродинамической модели.
1. Физическая модель газопылевой среды, взаимодействующей с излучением
Распределение межзвездного вещества чрезвычайно нерегулярно, и большая его часть собрана в плотные облака [3]. Из-за облаков межзвездной пыли десятки тысяч рассеянных звездных скоплений
остаются для нас невидимыми. Основным поглощающим и рассеивающим электромагнитное излучение агентом является пыль. Присутствие в межзвездной среде пыли влияет на характеристики излучения исследуемых небесных тел, на тепловой баланс межзвездного газа, причем пылинки могут содействовать как нагреву, так и охлаждению межзвездного газа. Радиационный энергобаланс определяется равенством прихода энергии за счет поглощения УФ- и видимого излучения, приходящего от источников, и собственного излучения пыли в ИК-диапазоне - расхода тепловой энергии на испускание фотонов. В большинстве случаев среда оказывается прозрачной для собственного излучения пыли, поэтому испущенные пылью фотоны беспрепятственно покидают среду, унося тепловую энергию пылинок.
1.1. Перенос излучения
Оптическая задача, с которой мы встречаемся при описании состояния плотного межзвездного облака с центральным источником - это задача о лучистом переносе энергии. Имеем набор частиц определенной формы, размера и состава, освещенный пучком данной интенсивности и частоты. Требуется найти интенсивность повсюду в облаке, для этого возникает необходимость решать уравнение переноса излучения, определяющее изменение интенсивности излучения при прохождении его через поглощающую, рассеивающую и излучающую среду [11].
Если рассмотреть слой облака, толщина которого много меньше его радиуса, то слои можно считать плоскопараллельными, а кривизной слоя пренебречь (рис. 1).
Под средой, состоящей из плоскопараллельных слоев, следует понимать среду, оптические свойства которой зависят только от одной координаты - высоты г [7]. Считаем облако зеркально-симметричным, это позволяет ограничить построение решения только для одной из половин облака по его толщине.
Газовую составляющую облака рассматриваем как идеальный совершенный газ с показателем адиабаты е. В дальнейшем полагаем е = 5/3.
Считаем, что источники сосредоточены в центре облака, то есть в плоскости г = 0, и равномерно распределены вдоль плоскости. Считаем, что источники излучают во всех направлениях изотропно.
Вместо пространственной координаты г удобно использовать оптическую толщу облака, отсчитываемую от центра облака:
г
т = да (г^г. (1)
0
Граничные условия в терминах оптической толщи естественно поставить следующим образом: на внутреннюю границу облака (г = 0) падает поток излучения 1Р0 от источника. Внешнюю
границу (г= г0) облака полагаем абсолютно черной, то есть излучение свободно покидает границу, что соответствует отсутствию источников диффузного излучения извне. Тогда интенсивность рассеянного излучения будет зависеть только от координаты г и направления излучения, характеризуемого зенитным углом I (рис. 1).
Расстояния между частицами пыли превосходят длину волны и во много раз больше радиуса пылинки. Это позволяет описывать процесс рассеяния на независимых частицах [2]. Кроме того, длина свободного пробега излучения 1/аг значительно меньше радиуса кривизны траектории луча. Поляризацию излучения не учитываем, поскольку в средах с большой оптической толщей существенно многократное рассеяние, которое является деполяризующим процессом. Поэтому «применимость скалярного уравнения переноса не вызывает сомнений» [7], во всяком случае учет поляризации излучения не приведет к существенному изменению интенсивности и потока диффузного излучения.
Таким образом, чтобы определить изменение интенсивности излучения при прохождении его через поглощающую, рассеивающую и излучающую среду, состоящую из плоскопараллельных слоев, необходимо решить уравнение
cos” dl. = . a(z)/,. + -(Z- ql.fpicos”,cos’"\)dcos e,. (z), (2)
где d(z) = t(z) + k(z) - коэффициент ослабления;
t(z) - коэффициент рассеяния; k(z) - коэффициент поглощения;
i(z) - коэффициент излучения, учитывающий внутренние и внешние источники; p(cos cos v) - индикатриса рассеяния [1], функция, определяющая распределение рассеянного излучения по углам;
Iv = Iv(z) - искомая интенсивность излучения на частоте v.
Для того чтобы для физической системы найти интенсивность излучения I, необходимо прежде определить все коэффициенты и индикатрису рассеяния.
1.1. Оптические параметры среды
Конкретные свойства коэффициентов ослабления, рассеяния и поглощения могут быть получены из расчетных данных по формулам для одной сферической частицы по теории Ми [13; 22], соответственно
а = пa2Qe (х); а = па 2Qs (х); k = na2Qa (х). (3)
2па
Здесь а - радиус частицы; х = —— - параметр Ми, а Q Q Q - факторы эффективности.
А
Реальная пылевая среда состоит из частиц разного размера и состава, поэтому при расчете характеристики светорассеяния необходимо это учитывать. Характеристикой разницы в размерах - дисперсности пыли - является функция распределения частиц по размерам п(а) [21; 25] для межзвездной среды:
n(a)da & а'35da, аН [0.05,0.25] мк, иначе п(а) = 0, (4)
где п - концентрация частиц.
Тогда суммарный коэффициент ослабления для всех пылинок равен
а = (rt(a)na Qe (a, m)da = nd(na Qe (a, m) f (a)da, (5)
00
где fa) = n(a)/nd - нормированная функция распределения
д/ (а^а = 1, (6)
0
а комплексный показатель преломления т = п - /7, в котором мнимая часть 7 характеризует поглощение света, а действительная часть п - рассеяние.
Аналогично формуле (5) для коэффициента ослабления, для коэффициентов рассеяния и поглощения имеем:
° = na($ a2qs (a m) f (a)da, (7)
0
k = n
qna 2Qa (a, m) f (a)da. (8)
Также, чтобы верно оценить вклад рассеяния во взаимодействие излучения со средой, необходимо определить индикатрису рассеяния. Поскольку для рассеяния коротковолнового излучения на пыли преобладает рассеяние вперед, то в качестве индикатрисы рассеяния УФ-излуче-ния в качестве аппроксимации можно выбрать модельную индикатрису Хеньи-Гринстейна [5]:
)=(1 - gHG -g: cos - у, (9)
где параметр 0 ungGHnu 1 характеризует степень вытянутости вперед индикатрисы: чем он больше и ближе к единице, тем сильнее вытянутость, а в случае Щ,:нЦ = 0 имеем изотропное рассеяние.
Объемные коэффициенты рассеяния и поглощения и индикатриса рассеяния согласно уравнению переноса (18) используются как входные параметры при расчете поля излучения в атмосфере. Проблема определения численных значений этих компонент решается путем вычислений оптических характеристик пылинок по заданным комплексному показателю преломления m, функции распределения частиц по размерам f(a) и концентрации частиц n а они, в свою очередь, получаются из экспериментальных измерений. Формулы (5), (7), (8) полностью решают задачу определения оптических характеристик ансамбля сферических частиц нескольких сортов.
0
1.3. Взаимодействие излучения с веществом облака
Межзвездная среда заполнена разреженным ионизованным и нейтральным газом и пылевой средой. Электромагнитное излучение при распространении в такой среде испытывает поглощение и рассеяние.
Основное поглощение в оптическом и УФ-диапазонах (0,1-0,5 мкм) обусловлено межзвездной пылью. В более коротковолновой далекой УФ- и рентгеновской области спектра (х < 0,1 мкм) основное поглощение связано с нейтральным газом [12] (водородом и другими химическими элементами). Для жестких рентгеновских фотонов и гамма-лучей при наличии свободных электронов в среде преобладающим становится комптоновское рассеяние.
Мы будем рассматривать взаимодействие вещества облака с ближним ультрафиолетовым излучением (0,2-0,38 мкм), а значит процессы ослабления излучения веществом можно рассматривать на пыли. Также необходимо учесть, что поглощая коротковолновое ультрафиолетовое излучение, пылинка переизлучает в непрерывном инфракрасном спектре (за счет этого пылевые частицы обычно холоднее газа). Большая часть поглощенного излучения звезд нагревает пылинку, а собственное инфракрасное излучение - охлаждает, унося тепловую энергию переизлученны-ми фотонами. Из баланса прихода и расхода энергии можно определить температуру пыли.
Скорость нагрева пылинки [1]
dE+ \ r d
— = ^ J d, (10)
Jv - средняя интенсивность по всем направлениям; t = cos ”. Эффективность поглощения сравнительно мало меняется с длиной волны, и можно переписать уравнение к виду
dE+
— Р к(I.)J, (11)
dt
где v* - некоторая эффективная частота поглощаемого пылью излучения звезд. Величина v*, а следовательно, и '{* зависит от размера пылинки, ее строения и химического состава. Охлаждается пылинка излучением со всей своей поверхности. Для сферических частиц скорость охлаждения равна
dF' 0
— = ф В, (Td )dv р к*, asJd', (12)
где Bv(Td) - функция Планка;
ГSB - постоянная Стефана-Больцмана;
к^ - коэффициент поглощения собственного излучения пылинки.
Пылинка как равновесная термодинамически система излучает по закону Кирхгофа-Планка [8]. Итого, полная энергия, поглощаемая веществом в тепловом равновесии, равна полной излучаемой им энергии, то есть выполняется условие лучистого равновесия [5]
J,^ а =ф в, (та )&. (1з)
0 0
В нашем рассмотрении, где мы учитываем два эффекта - нагрев УФ-излучением и выхолаживание излучением ИК-диапазона спектра, уравнение баланса лучистой энергии можно переписать в виде
k (V*)/ = к* о8БТ%. (14)
Заметим, что кроме излучения звезд, в межзвездной среде присутствуют и другие источники нагрева, например, космические лучи. Источником охлаждения газа за счет столкновений его молекул с пылинками служит излучение в линиях. Таким образом, в уравнении баланса (14) необходимо учесть дополнительные источники и стоки энергии. Мы знаем, что в отсутствие источника в центре плотного облака, его типичные температуры - 10 ф 20 К [1]. Тогда тепловое излучение облака в отсутствие прямого источника излучения и будет определяться этой типичной температурой. В этом случае уравнение баланса можно переписать следующим образом:
к( *)•/ + к ^7= *г 4 d4c (15)
Решая уравнение переноса излучения и рассчитывая среднюю интенсивность, а также опираясь на наблюдательные данные о межзвездной пыли и рассчитывая коэффициенты поглощения ультрафиолетового и инфракрасного излучения, из условия (15) найдем равновесную температуру пыли.
1.4. Радиационная гидродинамика
Характерные масштабы L межзвездного облака - парсеки и десятки парсек. Длину свободного пробега электронов, атомов и молекул I в межзвездном облаке можно оценить как I р 1/тн, где н - сечение взаимодействия атомов или молекул водорода. Беря для оценки снизу атомар-
ный водород, получаем н = !Ян2 : 3 X 1012 ф 3 X 1016 см (здесь Ян - радиус атома водорода), что на два-шесть порядков меньше минимального размера облака L. Это означает, что в одной элементарной ячейке размером г3 (I = г = L) содержится большое число частиц. Следовательно, вещество в облаке на масштабах задачи можно рассматривать как сплошную среду.
Гравитационное взаимодействие между разными частями системы сжимает вещество, оно стремится упасть само на себя. Такой процесс приводит к перераспределению массы - в одной области плотность растет, в другой уменьшается. Газ представляет собой сплошную среду, в нем присутствуют тепловые движения, которые могут препятствовать сжатию под действием самогравитации.
Мощные источники излучения, присутствующие в толще облака, также могут сильно влиять на приток (и отток) энергии и импульса к веществу, а потому и на его движение в целом. При описании динамики такого течения можно считать, что жидкость состоит как из частиц вещества, так и из фотонов [15]. А значит, и необходимо учитывать вклад обоих типов частиц в уравнение движения (или в уравнение гидростатического равновесия) и в закон сохранения энергии. Взаимодействие излучения и вещества описывается при помощи уравнений радиационной гидродинамики. Излучение будем описывать, опираясь на приближение о локальном термодинамическом равновесии (ЛТР).
Если система находится в гидростатическом равновесии, то в неподвижном облаке весу вышележащих слоев противостоит давление. Балансом этих двух сил и определяется равновесное распределение плотности в среде. Полное давление в случае ЛТР есть сумма давления совершенного газа и давления излучения.
Тогда уравнение, определяющее условие баланса сил, запишется в виде:
^ d' (г) 4п “
■*г=-рг)?(2> 1 )н(2’'м’ <|6>
где Н(2,\) - полный поток излучения.
Совместно с этим уравнением необходимо решить уравнение Пуассона
d2-<Z> = 4жGp( 2) (17)
аг
для гравитационного потенциала ' (г) и уравнение переноса излучения, определяющее изменение интенсивности излучения при прохождении его через поглощающую, рассеивающую и излучающую среду,
И~0к =' I + ! I (^,^£ч, (18)
где I - оптическая толща, а
а
т =----Г - (19)
а + к у ’
альбедо однократного рассеяния, определяющее вероятность того, что фотон рассеется, а не
поглотится,
&( (2) = к*42)В(Та) +1 Foё(ju - р0)ел«0 (20)
а (г) 2 - (20)
коэффициент излучения в предположении об ЛТР;
В(Та) - функция Планка;
^ - поток излучения от внешнего источника;
Г - косинус угла между нормалью к поверхности и направлением падения фотонов от источника.
Зная интенсивность излучения в среде, можем рассчитать среднюю интенсивность 4(2) и поток излучения НДг). Равновесную температуру пылевых частиц, как уже упоминалось выше, можно рассчитать из условия лучистого равновесия:
к (V.) 4 + /0= к,а8вТ4. (21)
Здесь /0 - слагаемое, описывающее вклад проникающего фонового излучения (космические лучи, рентген).
Вообще говоря, распределение температуры, обеспечивающее равновесие, не известно. Поэтому, воспользовавшись некоторым начальным приближением распределения температуры, необходимо путем итераций перестроить зависимость Т(г) так, чтобы поле излучения стало удовлетворять уравнению энергетического баланса (21).
Также хотелось бы отметить, что атомарный и молекулярный газ взаимодействуют с ультрафиолетовым излучением по-разному. Молекулярный водород поглощает в ультрафиолете (0,1-0,11 мкм), а среда из атомарного водорода для ультрафиолетового излучения прозрачна [1]. Для описания взаимодействия излучения и молекул необходимо решать уравнение переноса излучения в линиях [3; 23]. В нашей задаче мы не учитываем взаимодействие излучения и молекул.
2. Математическая модель
2.1. Решение уравнения переноса излучения методом ё-Эддингтона
На больших оптических глубинах с большой степенью точности выполняется соотношение между нулевым и вторым моментами интенсивности [8; 15]
/ _ К (2) _ 1
4 (г) 3 . (22)
Основываясь на этом результате, Эддингтон предположил, что соотношение 4 =3 К выполняется во всей среде. Тогда интенсивность излучения можно задать выражением [14]:
I = /0 + М, (23)
где /0 = 4 -средняя интенсивность;
I = 3Н - утроенный средний поток излучения в среде.
Рассмотрим один однородный слой с заданной оптической толщиной, известными альбедо однократного рассеяния и индикатрисой. Полученные для него соотношения потом можно обобщить на случай неоднородной многослойной среды.
Определим показатель асимметрии [6], то есть среднее значение косинуса угла рассеяния, и показатель, определяющий долю излучения, рассеянного вперед [14], следующим образом:
1+1
g = ~^Cp(9 (C°S9 )dcos9 ’ (24)
1
f = 2Qp(e) P2 (c°s9 )dcos9, (25)
здесь P1(cosz) = cos/ и P1(cosz) = (1/2)(3cos2/ + 1) - полиномы Лежандра первого и второго порядка; / - угол между падающим излучением в направлении а = (t, ё) и рассеянным в направлении аи = (ш,ёи). Как уже отмечалось ранее, рассеяние ультрафиолетового излучения происходит преимущественно вперед, то есть индикатриса имеет острый пик. С большой степенью точнос-
ти [6] этот пик можно представить в виде дельта-функции [16]. Усредняя по азимутальному углу ё, получаем индикатрису рассеяния в приближении ^-Эддингтона:
Pt- Ed (л) = 2 Ж л - + (1 - f )(1 + 3g * р$), (26)
где g* определяется следующим образом:
* g - f
g =~^j ■ (27)
Такую индикатрису p-_Ed и будем использовать при решении уравнения переноса излучения.
Разделим суммарное излучение в среде на однонаправленный прямой поток TF от источника, падающий под углом, для которого t = cos ”0, и диффузный поток
ш 1 * ~
Ы + — F0S(л - л0)ехР(- т /л0), 0 < л Ь 1
1 л v 2 (28) I, - 1Ьл<0
Тогда, подставив в уравнение (18) для данной частоты вместо индикатрисыр(Г) выражение в виде (26) и учтя разделение излучения на рассеянное и прямое, получим уравнение переноса в следующем виде [14]:
А—г +1 = — с (1 + 3g>AUI(A,AVT*)dAu+—т°(1 + 3g*W>)e"T Л0, (29)
dT
*
где т =(1 - of )т, о =■
1 - а/
Воспользовавшись приближением Эддингтона (23) и взяв нулевой и первый моменты по углам, получим разделение интегро-дифференциального уравнения переноса на систему
и 2п
двух ОДУ с искомыми величинами G = п/0 и Н —— /1 [14]. Решив эту систему, найдем диффузный поток и среднюю интенсивность. В итоге, получаем полный поток излучения F(т*) = 2Н + ц0жF0exp<- т*/и0).
Все предыдущие рассуждения данного пункта касались однородного слоя, в реальности среда может оказаться существенно неоднородной. Для описания многослойной среды разделим ее на слои таким образом, чтобы в каждом слое величины Т, г\ и g* были примерно постоянными.
Тогда всю среду можно будет смоделировать как последовательность слоев каждый со своими оптическими характеристиками.
Для /-го слоя (/' = 1, ... , N - 1), где I. - перенормированная оптическая толщина от верхней границы до уровня /, имеем уравнения, применимые в пределах I. и IV г.+1:
СХ-т - Х-т - т /ип
, = х2,_ е' + х2,.е ‘ - ае М0, (30)
Н, = - Рх21. /'Т + р.х21е-Х/Т - ~е■ ^, (31)
где х2/-1 и х2/ - неизвестные константы для /-го слоя, которые находятся из граничных условий.
Условие непрерывности внутреннего потока дополняется условиями на верхней и нижней границах. В итоге получаем систему из 2^2 линейных уравнений, решив которую найдем неизвестные константы х.
Будем считать, что внешняя по отношения к источнику граница облака - полностью поглощающая, это позволит говорить о том, что в итоге все излучение уходит из среды, а диффузный поток = 0, то есть граница не отражающая. Так как источник находится в центре облака, то
кроме прямого излучения от источника на вещество будет оказывать влияние отраженное диффузное излучение от другой половины облака. Поэтому следует полагать F0“ = F0“(f = 0).
Таким образом, полный поток и средняя интенсивность на каждом слое, которые необходимы для решения задачи, равны
Fi = 2Hi + nF0 №'
1 *
J G. + F0e'1 л 0
і ___ і 0 5
ЇЇ
(32)
(33)
здесь F0 - поток радиации, идущей от источника.
2.2. Процедура расчета и управляющие параметры
Для того чтобы построить распределение температуры Т^) и плотности г(2), нужно решить систему уравнений (5), (7), (8), (15)—( 17), (29), которая фактически разбивается на две связанные между собой подсистемы. Назовем их гидростатической - это уравнение гидростатического равновесия и уравнение Пуассона для гравитационного потенциала системы и оптической - уравнение переноса излучения и уравнение баланса лучистой энергии.
В расчетной модели учитываются два источника прогревающего излучения: точеный источник, расположенный в начале координат, отвечающий излучению звезды, и распределенный по z источник фонового излучения, поддерживающий в отсутствии точечного источника равновесную температуру, во всем облаке равную Т0. В расчетах принималось Т0= 10 К.
Вообще говоря, равновесное распределение температуры нам не известно. Поэтому для некоторого начального приближения распределения температуры, в отсутствии источника излучения, решаем гидростатическую подсистему, из которой находим некоторое распределение плотности. Используя полученное распределение плотности, находим оптические параметры задачи для неоднородной среды: оптическую толщину I , альбедо однократного рассеяния (индикатрису рассеяния выбираем какую-либо модельную). Для полученных параметров решаем уравнение переноса излучения и из уравнения лучистого равновесия находим температуру.
Теперь сравним полученное температурное распределение с учетом излучения от источника с предыдущим, для которого было получено распределение плотности. Для этого определим невязку как
Є \т-С Zi) - T (z. )\
Є \T (Zi )\
(34)
Рис. 2. Итерационная процедура расчета распределения температуры и плотности
/і 0
___ i 1
i=1
Если ц больше некоторой заданной точности ё то по новому найденному значению температуры Т^) = Тг"^) находим распределение плотности и т. д. (см. рис. 2). И так продолжаем итерации до тех пор, пока не выполнится условие ц V ё . Полученное итерациями температурное распределение и, соответственно, плотность и будут интересующими нас равновесными распределениями. В задаче есть два естественных пространственных масштаба: гравитационный масштаб
,____ _ [ЛР, _ ,
I а = смЦGi 0, где С*с = -у ^=0 - скорость звука в центре облака; г0 - плотность в центре
облака и длина свободного пробега фотона 1рк = а01, определяемая через коэффициент ослабления а0 в центре облака. Отношение этих масштабов определяет первый безразмерный управляющий параметр задачи:
D = $0? /, (35)
имеющей смысл характерной оптической толщины среды.
Вторым безразмерным параметром можно было бы взять отношение характерных скоростей в задаче - скорости света с и скорости звука сх0. Однако вместо этого удобнее использовать другую безразмерную конструкцию:
^^ 4-1^- (36)
I с 15П Й3 (36)
1рк ° ипн0й
= г 0
параметр, который мы называем числом Соболева 2 Здесь пн 0 = т - концентрация газа в центре.
тн
А именно, число Соболева будет определять относительный вклад энергии излучения только от точечного источника в общий силовой баланс, причем вклад, естественно, зависит от оптической толщины среды. Фактически, сравнивается энергия излучения в цилиндре, высотой, равной длине свободного пробега фотона 1рк, и тепловая энергия в цилиндре, высота которого равна джинсковской длине Сила лучистого давления зависит от длины свободного пробега. Если 1рН > zh (где zh - полутолщина облака, определяемая джинсовским масштабом), то излучение свободно покидает среду, не оказывая силового воздействия на вещество облака. При 1рк, меньших полутолщины облака, излучение будет
оказывать силовое давление на вещество, определяющееся столкновениями фотонов и частиц среды.
Оценим диапазон значений, которые может принимать число Соболева для характерных параметров нашей задачи. Типичная концентрация газа в плотных облаках порядка 102-105 см-3. Характерную температуру выбираем 100 К. Такая температура соответствует температуре холодной фазы двухфазной модели среды. Отвечающая такой температуре звуковая скорость ем р 1 км/с. Тогда для таких параметров число Соболева будет лежать в пределах от 3 ф 102 до 6 ф 104. Такие большие числа Соболева обусловлены выбором нормировки на 100 К, а само число Sv(T) ф Т3,5.
Третьим безразмерным параметром является поток излучения точечного источника, нормированный на поток тепловой энергии, излучаемой с некоторой поверхности, обладающей температурой Т0:
£ = ^/г 8Х. (37)
Первые два безразмерных параметра характеризуют среду, а третий параметр характеризует источник.
В расчетах число узлов расчетной сетки по координате z принималось равным 300.
3. Результаты численного моделирования
Приведем рассчитанные распределения плотности и температуры газа для различных параметров задач.
Вначале рассмотрим случай отсутствия точечного источника. Этому отвечает параметр Sv = 0. При наличии только фонового прогревающего излучения температура остается постоян-
ной и равной Т=10 К. Варьируя характерную оптическую толщину D, мы фактически рассматриваем равновесные распределения с различной массой облака (рис. 3). Из рисунка 3 видно, что при увеличении плотности газ уплотняется под действием собственной гравитации.
Теперь построим распределения концентрации и температуры газа при включенном точечном источнике с заданным потоком F0 = 0,567 эрг/с ф см2 и различных числах Соболева Sv1 = 3 ф 105 и Sv1 = 6 ф 105. При этом параметр D всегда варьируется в диапазоне от 1 до 10.
Увеличение числа Соболева отвечает уменьшению длины свободного пробега излучения. В свою очередь, это приводит к более эффективному взаимодействию излучения со средой. Как видно из рисунков 4 и 5, температура среды в центре облака выше и спадает с увеличением расстояния от источника.
Из сравнения рисунков 6 и 7 видно, что при увеличении числа Соболева происходит качественная перестройка в распределении вещества в центральной части облака. В случае большего Sv газ эффективнее прогревается в центре, что заставляет вещество под действием его давления уйти на большие высоты и сформировать на высотах порядка z р 0,1-0,3 пк плотную оболочку, поддерживаемую как высоким давлением излучения, так и высоким термодинамическим давлением газа в прогретой центральной каверне. С другой стороны, внешние слои эффективно охлаждаются из-за собственного излучения и, в то же время, излучение, идущее от источника, сильно ослабляется внутренними слоями. Следовательно, не происходит эффективного разогрева внешних слоев. А значит, формируется инверсное распределение температуры, что в динамическом рассмотрении неизбежно приведет к конвективной неустойчивости.
Рис. 3. Распределение концентрации при отсутствии точечного источника излучения
Рис. 4. Температурное распределение. D = 1 ы 10, Sv =3 ф 10
п,см- п,см
х = 13.0029 ' I ^ 12.2415 ьМ11.4515
тЦю.бгбв ^=атвб8 -
5 = 8^1
х = 7.82Ж
ъ = 67118\Дт,
X =_5.413 \\\\\\ \^
1: = 3.721 \ \ \ \ \ '
I ' -
г,пк
Рис. 5. Температурное распределение. D = 1 ы 10, Sv = 6 ф 104
Рис. 6. Распределение концентрации. D = 1 ы 10, Sv = 3 ф 104
Рис. 7. Распределение концентрации. D = 1 ы 10, Sv = 6 ф 104 ISSN 2222-8896. Вестн. Волгогр. гос. ун-та. Сер. 1, Мат. Физ. 2012. № 1 (16) 69
ъ= 7.0046
О 0.05 0.1 0.15 0.2 0.25 0.3
г,пк
Рис. 8. Температурное распределение. Поток падающего излучения F, Sv = 6 ф 102
Я = 14,5793 1 1 1
13.5606
*М2.5259
т =Цт\4707
Ш^3881
г = 8.08£^^к
х = 6.8193
г = 5.3928 \\\ \
т = 3.5887 I
г,пк
Рис. 9. Температурное распределение. Поток падающего излучения 10 F1, Sv = 6 ф 102
(О о! II I I
X = 6.3453
т = 5.6855
Т = 5.0238
т = а.я^й4
~ т= 3.6864 ~
т = я ппяа
_ г = 2.3022 _
х = 1.5702
ъ = 0.77196
°0 0.05 0.1 0.15 0.2 0.25 0.3
2, ПК
Рис. 10. Распределение концентрации. Поток падающего излучения F1, Sv = 6 ф 102
Выявим влияние интенсивности излучения источника на распределение плотности и температуры в облаке. Для этого рассмотрим случаи с различающимися потоками излучения от источника F0 = 5,67 эрг/с ф см2 и F0 = 56,7 эрг/с ф см2. При этом число Соболева в обоих случаях совпадает - Sv = 6 ф 105, параметр D изменяется в одном диапазоне от 1 до 10. Как и следовало ожидать, увеличение потока излучения приведет к увеличению температуры в слоях облака, расположенных вблизи к источнику.
Из расчетов видно (рис. 8 и 9), что на малых высотах существуют области с высокой температурой, по сравнению с периферией облака, и более низкой плотностью, чем в вышележащих слоях (рис. 10 и 11). При увеличении высоты облака температура спадает, а концентрация газа увеличивается. Такое инверсное распределение в динамике может приводить к возникновению конвективной неустойчивости.
1 х = 14.5793 1 1 1 1
/ ч= 13.5606 N.
- / / т= 12.5259 хХ, -
- ' ’ / / Т= 11.4707 -
- / / / 1= 10.3881 -
/ / -с = 9.268
'п = 8 086^___ ^ = 6;8193__ —— . ц = 5.3928 ” Т^Т£яя7 1 1 1 1 1 -
О1 1 1 1 1 1
а 0.05 0.1 0.15 0.2 0.25 0.3
г,пк
Рис. 11. Распределение концентрации. Поток падающего излучения 10 Г, Sv = 6 ф 102
5. Выводы
1. Разработана гидростатическая модель газопылевого межзвездного облака, находящегося в равновесии с собственной тяжестью и внешним излучением. Модель включает в себя следующие составляющие: описание взаимодействия пылинок с излучением, учитывающее полидисперсность пыли; расчет переноса излучения в приближении ^-Эддингтона; расчет гидростатического баланса вещества.
2. Выполненные численные расчеты стационарных распределений в облаке показывают, что наиболее благоприятные условия для развития конвективной неустойчивости возникают при наличии мощного источника излучения, окруженного оптически плотным для его излучения газопылевым облаком.
ПРИМЕЧАНИЯ
1 Работа выполнена при финансовой поддержке грантов РФФИ N° 11-02-01332а, РФФИ-Поволжье № 11-02-97124р_поволжье_а и фонда некоммерческих программ «Династия».
2 В. В. Соболев - российский советский астрофизик, академик АН СССР.
СПИСОК ЛИТЕРА ТУРЫ
1. Бочкарев, Н. Г Основы физики межзвездной среды / Н. Г. Бочкарев. - М. : Изд-во МГУ, 1991. - 352 с.
2. Ван де Хюлст, Г. Рассеяние света малыми частицами : пер. с англ. / ван де Г. Хюлст. - М. : Изд-во иностр. лит., 1961. - 536 с.
3. Каплан, С. А. Физика межзвездной среды / С. А. Каплан. - М. : Наука, 1979. - 592 с.
4. Колесник, И. Г. Образование гигантских молекулярных облаков в сверхоблаках и возникновение сверхзвуковой турбулентности / И. Г Колесник // Кинематика и физика небесных тел. - 1987. - Т. 3. - С. 50-61.
5. Лекции по теории переноса излучения : учеб. пособие / Д. И. Нагирнер. - СПб. : Изд-во С.- Петерб. ун-та, 2001. - 284 с.
6. Ленобль, Ж. Перенос радиации в рассеивающих и поглощающих атмосферах : пер. с англ. / Ж. Ле-нобль. - Л. : Гидрометеоиздат, 1990. - 263 с.
7. Минин, И. Н. Теория переноса излучения в атмосферах планет / И. Н. Минин. -М. : Наука, Гл. ред. физ.-мат. лит., 1988. - 264 с.
8. Михалас, Д. Звездные атмосферы : пер. с англ. / Д. Михалас. - М. : Мир, 1982. - 352 с.
9. Огульчанский, Я. Ю. Модель эволюции сверхзвуковых движений в молекулярных облаках и характеристики фрагментированной среды / Я. Ю. Огульчанский // Кинематика и физика небесных тел. - 1992. - Т. 8, № 3. - С. 3-13.
10. Огульчанский, Я. Ю. Начальная стадия развития сверхзвуковой турбулентности при образовании гигантских молекулярных облаков / Я. Ю. Огульчанский // Кинематика и физика небесных тел. - 1988. - Т. 4, № 5. - С. 3-12.
11. Соболев, В. В. Курс теоретической астрофизики / В. В. Соболев. - 3-е изд., перераб. - М. : Наука, Гл. ред. физ.-мат. лит., 1985. - 504 с.
12. Спитцер, Л.-мл. Физические процессы в межзвездной среде : пер. с англ. / Л.-мл. Спитцер. - М. : Мир, 1981. - 350 с.
13. Шифрин, К. С. Рассеяние света в мутной среде / К. С. Шифрин. - М. ; Л., 1951. - 289 с.
14. Briegleb, B. P. A delta-Eddington multiple scattering / B. P. Briegleb, B. Light // National Center for Atmospheric Research, 2007. - 100 с.
15. Computational Methods for Astrophysical Fluid Flow / R. J. LeVeque [et al.]. - N. Y. : Springer-Verlag Berlin Heidelberg, 1998. - 508 с.
16. Joseph, J. H. The delta-eddington approximation for Radiative Flux Transfer / J. H. Joseph, W. J. Wicombe // Journal oe the atmospheric sciences. - 1976. - Vol. 33 - P. 2452-2459.
17. Koyama, H. An Origin of Supersonic Motions in Interstellar Clouds / H. Koyama, S.-I. Inutsuka // Astrophys. J. - 2002. - V. 564. - P. L97-L100.
18. Larson, R. B. Turbulence and star formation in molecular clouds / R. B. Larson // Mon. Not. Roy. Astron. Soc. - 1981. - V 194. - P. 809-826.
19. Lazarian, A. Turbulence statistics from spectral line observations / A. Lazarian // Journal of The Korean Astronomical Society. - 2004. - V. 37. - P. 563-570.
20. Mac Low, M.-M. Control of star formation by supersonic turbulence / M.-M. Mac Low, R. S. Klessen // Rev. Modern Physics. - 2004. - V 76. - P. 125-194.
21. Mathis, J. The Size Distribution of Interstellar Grain / J. Mathis, W. Rumpl, K. Nordsieck // Astrophys. J. -1977. - V 217. - P. 425-433.
22. Mie, G. // Ann. de Phys. - 1908. - V. 25. - Р 377.
23. Mihalas, D. Foundation of Radiation Hydrodinamic / D. Mihalas, B.W. Mihalas. - N. Y. : Oxford Universit press, 1984. - 718 р.
24. Pan, L. The temperature of interstellar clouds from turbulent heating / L. Pan, P. Padoan // Astrophys. J. -2009. - V 692. - Р 594-607.
25. Tiliens, A. G. G. M. Physics and chemistry of the interstellar vedium / A. G. G. M. Tiliens. - Cambridge : Cambridge University Press, 2005. - 491 р.
26. Weizsacker, C. F. von. // Astrophys. J. - 1951. - V. 114. - P. 165.
27. Weizsacker, C. F. von. In IAU and International Union of Theoretical and Applied Mechanics Symposium / C. F. von Weizsacker. - 1949. - P. 158.
28. Wicombe, W. J. Delta-Eddington Approximation for a Vertically Inhomogenus Atmosphere / W. J. Wicombe // National Center for Atmospheric Research, 1977. - 66 р.
A HYDROSTATIC MODEL OF SELF-GRAVITATING OPTICALLY DENSE INTERSTELLAR CLOUD
E. V. Zhukova, A.M. Zankovich, I.G. Kovalenko, K.M. Firsov
We consider a hydrostatic model of an optically dense gas and dust of interstellar clouds, which is in equilibrium with its own weight and external radiation. Numerical calculations of stationary distributions of temperature and concentration of gas in the cloud. Identified conditions for the emergence of convective instability.
Key words: interstellar dust, radiative transfer, gas-dust clouds, interaction of radiation with matter, convective instability, Eddington approximation.
