Научный журнал
Фундаментальные исследования
ISSN 1812-7339
"Перечень" ВАК
ИФ РИНЦ = 1,074

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ТУРБУЛЕНТНЫХ ТЕЧЕНИЙ В ПРОТОПЛАНЕТНОМ ДИСКЕ

Кукса М.М. 1 Маров М.Я. 2
1 Обнинский институт атомной энергетики
2 ФГБУН «Институт геохимии и аналитической химии им. В.И. Вернадского» Российской академии наук
Рассматривается замкнутая система уравнений сжимаемой магнитной гидродинамики (МГД) масштаба среднего движения. Система уравнений положена в основу двумерной модели осесимметричного протопланетного диска в цилиндрических координатах. При построении модели наряду с традиционным теоретико-вероятностным осреднением МГД-уравнений используется весовое осреднение Фавра. Применен модифицированный коэффициент турбулентной кинематической вязкости, учитывающий влияние магнитного поля на турбулентное течение через посредство пути смешения. Разработан параллельный код для численного решения системы осредненных МГД-уравнений. В процессе вычислительных экспериментов подобраны «свободные» граничные условия, позволяющие веществу свободно проникать в область моделирования и покидать ее. При этом установившийся темп аккреции газа на звезду согласуется с наблюдательными данными о классических звездах на стадии Т Тельца. В результате численного моделирования получены распределения средней плотности, гидродинамической скорости и магнитной индукции в турбулизованном протопланетном диске на расстоянии 1 а.е. от звезды.
протопланетный диск
турбулентная магнитная гидродинамика
численное моделирование
1. Колесниченко А.В., Маров М.Я. К магнитогидродинамическому моделированию протопланетного диска Солнца // Астрономический вестник. – 2009. – Т. 43. – № 5. – С. 424–448.
2. Колесниченко А.В., Маров М.Я. О влиянии спиральности на эволюцию турбулентности в солнечном протопланетном облаке // Астрономический вестник. – 2007. – Т. 41. – № 1. – С. 3–23.
3. Колесниченко А.В., Маров М.Я. Термодинамическая модель МГД-турбулентности и некоторые ее приложения к аккреционным дискам // Астрономический вестник. – 2008. – Т. 42. – № 3. – С. 1–50.
4. Balbus S.A., Hawley J.F. A powerful local shear instability in weakly magnetized disks. I. Linear analysis // The Astrophysical Journal. – 1991. – Vol. 376. – P. 214–222.
5. Brandenburg A., Nordlund A., Stein R.F., Torkelsson U. Dynamo-generated turbulence and large-scale magnetic fields in a keplerian shear flow // The Astrophysical Journal. – 1995. – Vol. 446. – P. 741–754.
6. Campbell C.G. Magnetohydrodynamics in binary stars. – Dordrecht: Kluwer Academic Publishers, 1997. – 306 p.
7. Kuksa M.M. The magnetic field impact in accretion rate in a protoplanetary disk // Odessa Astronomical Publications. – 2012. – Vol. 25. – Issue 2. – P. 104–106.
8. Lin D.N.C., Papaloizou J. On the structure and evolution of the primordial solar nebula // Monthly Notices of the Royal Astronomical Society. – 1980. – Vol. 191. – P. 37–48.
9. Lin D.N.C., Pringle J.E. A viscosity prescription for a self-gravitating accretion disc // Monthly Notices of the Royal Astronomical Society. – 1987. – Vol. 225. – P. 607–613.
10. Lynden-Bell D., Pringle J.E. The evolution of viscous discs and the origin of the nebular variables // Monthly Notices of the Royal Astronomical Society. – 1974. – Vol. 168. – P. 603–637.
11. Shakura N.I., Sunyaev R.A. Black holes in binary systems. Observational appearance // Astronomy and Astrophysics. – 1973. – Vol. 24. – P. 337–355.
12. Shakura N.I., Sunyaev R.A., Zilinkevich S.S. On the turbulent energy transport in accretion discs // Astronomy and Astrophysics. – 1978. – Vol. 62. – P. 179–187.
13. Velikhov E.P. Stability of an ideally conducting liquid flowing between cylinders rotating in a magnetic fluid // Sov. Phys. JETP. – 1959. – Vol. 9. – P. 995–998.

Реконструирование аккреционного протопланетного диска является одной из фундаментальных проблем космогонии Солнечной системы. До сих пор открытым остается вопрос о механизме переноса углового момента в аккреционных дисках, поскольку молекулярная вязкость не может обеспечить темп аккреции, наблюдаемый в дисках вокруг молодых звезд солнечного типа. Н.И. Шакура и Р.А. Сюняев [11], Д. Линден-Белл и Дж. Прингл [10] предположили, что повышенная вязкость может быть вызвана турбулентностью. В качестве источников турбулентности предлагались конвекция [8], нелинейная гидродинамическая неустойчивость [12], гравитационная неустойчивость [9] и возмущения, вызванные внешними воздействиями, однако ни один из них не мог обеспечить перенос углового момента за требуемое время. Прогресс был достигнут с привлечением магниторотационной неустойчивости, открытой Е.П. Велиховым [13] и развитой Ст. Бальбусом и Дж. Хоули [4]. Численные расчеты [5] показали, что магниторотационная неустойчивость порождает турбулентность, которая генерирует и поддерживает магнитное поле в присутствии диссипации.

Существование даже слабого магнитного поля существенно усложняет гидродинамические течения в протопланетном диске. В работах А.В. Колесниченко и М.Я. Марова [1, 3] в приближении одножидкостной магнитной гидродинамики получена замкнутая система магнитогидродинамических уравнений масштаба среднего движения, предназначенная для моделирования сдвиговых и конвективных турбулентных течений слабо ионизованной дисковой среды в присутствии магнитного поля. Данная модель послужила отправной точкой для настоящего исследования.

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

Математическая модель

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

Во-первых, теоретико-вероятностное осреднение по Рейнольдсу

Eqn69.wmf

где мгновенное значение величины A представляется в виде суммы осредненной Eqn70.wmf и пульсационной A′ составляющих

Eqn71.wmf

Во-вторых, средневзвешенное осреднение по Фавру

Eqn72.wmf

где мгновенное значение величины A есть сумма средневзвешенного значения ⟨A⟩ и турбулентной флуктуации A′′

Eqn73.wmf

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

Eqn74.wmf (1)

где t – время; Eqn75.wmf и Eqn76.wmf – осредненная плотность и ее начальное значение соответственно; Eqn77.wmf – средневзвешенная скорость; Eqn78.wmf – начальная скорость звука; γ = 5/3 – показатель адиабаты; Φgrav – гравитационный потенциал, создаваемый протосолнцем (самогравитация диска не учитывается); Eqn79.wmf – осредненный симметричный тензор деформаций; νturb – коэффициент турбулентной кинематической вязкости; ηturb – коэффициент турбулентной магнитной диффузии; Eqn80.wmf – осредненное собственное поле диска; Bext – внешнее постоянное поле (например, создаваемое протосолнцем или остаточное поле межзвездной среды) Eqn81.wmf – результирующее осредненное магнитное поле.

Применительно к модели протопланетного диска дифференциальные операторы записываются в цилиндрической системе координат (r, φ, z) с учетом ∂/∂φ = 0, поскольку диск считается симметричным относительно оси вращения Oz.

Для моделирования коэффициентов турбулентного переноса используются два подхода. Первый подход основан на (стандартной) α-модели турбулентной вязкости [11]

Eqn82.wmf (2)

где Eqn83.wmf – средневзвешенная угловая скорость жидкой частицы в дифференциально вращающемся диске.

Второй подход основан на модифицированной модели [2], учитывающей влияние магнитного поля на турбулентное течение через посредство пути смешения

Eqn84.wmf (3)

где путь смешения

Eqn85.wmf

где в свою очередь магнитогидродинамическое число Ричардсона

Eqn86.wmf

Коэффициент турбулентной магнитной диффузии считается постоянным и, согласно [6], определяется по аналогии со стандартной моделью

Eqn87.wmf (4)

где ΩK,1 au – кеплеровская угловая скорость на расстоянии 1 а.е. от звезды.

Эмпирический параметр α принимается равным 0,01, что является общепринятым значением в моделях газовых протопланетных дисков.

Численный метод, начальные и граничные условия

Автором разработаны программные модули для среды Matlab, в которых реализована апробированная ранее явная численная схема третьего порядка аппроксимации [7] и средства визуализации результатов моделирования.

Решение системы дифференциальных уравнений (1) ищется в узлах регулярной ортогональной сетки с числом узлов N = Nr×Nz, где Nr – число узлов по оси Or, Nz – число узлов по оси Oz. С каждой стороны сетки дополняется тремя слоями граничных узлов, в которых задаются граничные условия.

Пространственные производные приближаются конечными разностями шестого порядка аппроксимации

Eqn88.wmf(5)

Eqn89.wmf(6)

где Δr – шаг сетки по радиусу. Производные по z вычисляются аналогичным образом.

Производная по времени вычисляется с помощью явной трехэтапной схемы Рунге‒Кутта третьего порядка аппроксимации. Шаг по времени Δt на каждом временном слое определяется из условия Куранта‒Фридрихса‒Леви, что обеспечивает сходимость явной численной схемы.

В начальный момент времени плазма равномерно распределена в области моделирования с плотностью Eqn90.wmf, что при начальной скорости звука Eqn91.wmf и толщине диска h = 0,05 а.е. дает начальную поверхностную плотность Eqn92.wmf. Начальная радиальная и вертикальная скорости равны нулю Eqn93.wmf. Начальная азимутальная скорость кеплеровская Eqn94.wmf. Собственное магнитное поле диска в начальный момент времени отсутствует Eqn95.wmf. Внешнее магнитное поле имеет небольшую вертикальную составляющую Eqn96.wmf.

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

Результаты исследования и их обсуждение

Результаты численного моделирования с модифицированным коэффициентом турбулентной вязкости (3) представлены на рис. 1, 2.

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

Eqn97.wmf (7)

Результаты численного моделирования с коэффициентом (2) и с модифицированным коэффициентом (3) имеют следующие общие черты:

Параметр

Граничные условия

в радиальном направлении

в вертикальном направлении

слева

справа

сверху

снизу

Плотность

односторонняя производная

односторонняя производная

фиксированное значение

фиксированное значение

Радиальная скорость

копирование крайнего значения

копирование крайнего значения

асимметричная экстраполяция

асимметричная экстраполяция

Азимутальная скорость

асимметричная экстраполяция

асимметричная экстраполяция

асимметричная экстраполяция

асимметричная экстраполяция

Вертикальная скорость

асимметричная экстраполяция

асимметричная экстраполяция

копирование крайнего значения

копирование крайнего значения

Радиальная индукция

копирование крайнего значения

копирование крайнего значения

копирование крайнего значения

копирование крайнего значения

Азимутальная индукция

копирование крайнего значения

копирование крайнего значения

копирование крайнего значения

копирование крайнего значения

Вертикальная индукция

копирование крайнего значения

копирование крайнего значения

копирование крайнего значения

копирование крайнего значения

pic_23.tif

Рис. 1. Осредненные плотность и гидродинамическая скорость в диске

– Радиальная и азимутальная компоненты магнитного поля меняют знак при переходе через экваториальную плоскость диска.

– Собственное магнитное поле диска развивается (растет по модулю), однако для выяснения причины этого явления требуются дальнейшие исследования.

– В обоих случаях подтвердилась гипотеза, высказанная в [1] о том, что азимутальная компонента магнитной индукции Eqn98.wmf гораздо сильнее меняется по высоте, чем по радиусу, а значит, дает основание учитывать ее градиент в моделях турбулентного переноса.

Результаты, полученные с модифицированным коэффициентом, отличаются тем, что:

– Вертикальная компонента Eqn99.wmf имеет противоположный знак;

– Магнитное поле развивается медленнее, чем в первом случае.

Выводы

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

pic_24.tif

Рис. 2. Осредненная собственная магнитная индукция в диске

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

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

Рецензенты:

Камаев Д.А., д.т.н., заведующий лабораторией, ФГБУ Нпо «Тайфун», г. Обнинск;

Перегуда А.И., д.т.н., профессор, Калужский филиал Московского финансово-юридического университета, г. Малоярославец.

Работа поступила в редакцию 01.08.2013.


Библиографическая ссылка

Кукса М.М., Маров М.Я. ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ТУРБУЛЕНТНЫХ ТЕЧЕНИЙ В ПРОТОПЛАНЕТНОМ ДИСКЕ // Фундаментальные исследования. – 2013. – № 10-1. – С. 35-39;
URL: http://www.fundamental-research.ru/ru/article/view?id=32210 (дата обращения: 13.11.2019).

Предлагаем вашему вниманию журналы, издающиеся в издательстве «Академия Естествознания»
(Высокий импакт-фактор РИНЦ, тематика журналов охватывает все научные направления)

«Фундаментальные исследования» список ВАК ИФ РИНЦ = 1.074