ГАЗОДИНАМИЧЕСКАЯ СТРУКТУРА ТЕЧЕНИЯ В СИСТЕМЕ SS CYG Жилкин А.Г. 1,2, Бисикало Д.В. 1 1 Институт астрономии РАН 2 Челябинский государственный университет Всероссийская астрономическая конференция "Космические рубежи XXI века" (ВАК ) сентября 2007 г., Казань
Полуразделенная двойная система Звезда- донор Звезда- аккретор Полость Роша Аккреционный диск Точка Лагранжа L 1
Постановка задачи Внутренняя критическая поверхность отождествляется с полостью Роша в ограниченной задаче трех тел. В системе отсчета, вращающейся с угловой скоростью системы вокруг общего центра масс, потенциал Роша определяется выражением: где G – гравитационная постоянная, r a и r d – радиус-векторы центров аккретора и донора, а M a и M d – их массы, соответственно, r c – радиус-вектор центра масс системы. Модуль вектора угловой скорости где A – расстояние между компонентами системы.
Постановка задачи Обмен веществом в двойной системе происходит через внутреннюю точку Лагранжа L 1, поскольку в этой точке градиент давления не уравновешивается силой тяготения. Скорость вытекания вещества равна местной скорости звука. Для описания течения в двойной системе используем декартову систему координат (x, y, z), заданную следующим образом. Начало координат находится в центре звезды-аккретора r a = (0,0,0). Центр звезды-донора лежит на расстоянии A на оси x от центра аккретора: r d = (-A,0,0). Ось z направлена вдоль оси вращения: = (0,0, ).
Магнитное поле аккретора
Звезда-аккретор обладает собственным дипольным магнитным полем, индукция которого определяется выражением: где - магнитный момент аккретора. Будем считать, что собственное вращение аккретора является синхронным. Поэтому в выбранной системе отсчета вектор магнитного момента не меняется с течением времени и имеет следующие компоненты: где – угол наклона вектора относительно оси z, а – угол между осью x и проекцией вектора на плоскость xy.
Основные уравнения
Функции нагрева и охлаждения Рис. 1. Зависимости функций нагрева и охлаждения от температуры T (сплошные линии). Пунктирными линиями показаны их линейные аппроксимации в окрестности равновесной температуры T=13600K.
Переход к безразмерным переменным В качестве масштабов длины и времени выберем следующие величины: В качестве масштаба плотности выберем значение плотности во внутренней точке Лагранжа Остальные масштабные соотношения можно записать в виде:
Исключение фонового поля Во многих практических задачах рассматриваемого типа фоновое магнитное поле аккретора (белый карлик или нейтронная звезда) во внутренних областях аккреционного диска оказывается очень сильным. Это может привести к накоплению ошибок в процессе численного моделирования при операциях с большими числами. Для уменьшения таких ошибок можно использовать расщепление полного магнитного поля B на фоновое поле B * звезды-аккретора и поле b, индуцированное токами в аккреционном диске и во внешней оболочке: Фоновое магнитное поле аккретора B * не зависит от времени и является потенциальным:
Адаптивная сетка Основана на преобразовании к произвольной нестационарной системе координат. Скорость сетки Рис. 2. Структура расчетной сетки в задаче о численном моделировании объемно распределенного взрыва в момент времени t = 0.26 для случая h =
Разностная схема 1. Аппроксимирует уравнения в произвольной нестационарной системе координат (адаптивная сетка). 2. Учитывает дополнительный член в уравнении индукции для очистки дивергенции магнитного поля (8-волновой метод). 3. Схема годуновского типа повышенного порядка точности. Конечно-разностная схема относится к классу схем, удовлетворяющих принципу неувеличения полной вариации решения (TVD), не требует вычисления собственных векторов матрицы гиперболичности, а использует только максимальные по модулю ее собственные значения (спектральный радиус). Она имеет третий порядок аппроксимации по пространственной переменной в областях гладкости решения и первый порядок аппроксимации по времени. На основе развитой конечно-разностной схемы разработан соответствующий параллельный численный код. Адаптация численного кода к многопроцессорной архитектуре проводилась путем разбиения расчетной области по трем осям системы координат в переменных,,. Для проверки производительности параллельного кода были проведены тестовые численные расчеты на вычислительном кластере Института астрономии РАН.
Тестовые расчеты Для проверки правильности работы кода были проведены тестовые численные расчеты следующих задач. 1. Задача Римана о распаде произвольного разрыва (газодинамический и МГД варианты). 2. Задача об объемно-распределенном взрыве в среде, пронизанной однородным магнитным полем. 3. Задача об аккреции изотермического ( = 1.001) и адиабатического ( = 5/3) газа на гравитирующую точечную массу. Следует учесть, что в первых двух тестовых задачах, получаемые численные решения отличаются от истинных обобщенных решений соответствующих газодинамических и МГД задач. Это связано с тем, что развитый численный код в качестве одной из консервативных переменных использует плотность энтропии, а не плотность энергии газа. В результате скорости сильных разрывов и значения скачков величин на них отличаются от истинных. Однако в задачах моделирования аккреции в полураздельных двойных системах сильных ударных волн, как правило,не возникает. Поэтому в этих случаях получаемые с помощью кода численные решения описывают картину возникающего течения с достаточно хорошей точностью.
Тестирование производительности Рис. 3. Эффективность (слева) и ускорение (справа) работы параллельной программы в зависимости от числа процессоров.
Тестирование производительности Рис. 4. Изолинии плотности и разбиение сетки по процессам в параллельном расчете задачи об объемно распределенном взрыве в момент времени t = 0.26.
Система SS Cyg Масса донора M d = 0.56 M, эффективная температура донора 4000 K Масса аккретора M a = 0.97 M Период обращения системы P orb = сек Большая полуось A = см. Внутренняя точка Лагранжа L 1 находится на расстоянии примерно 0.56 A от центра аккретора. Магнитное поле на поверхности аккретора, оцениваемое из наблюдений, составляет величину порядка B a = Гс. Отсюда можно оценить магнитный момент аккретора: = B a R a 3, где R a =0.68 R - радиус звезды-аккретора. Звезда-донор – красный карлик Звезда-аккретор – белый карлик Основные параметры системы:
Разностная сетка В расчетах структуры течения в системе SS Cyg использовалась геометрически адаптивная сетка (h = 0). Для более высокого разрешения вертикальной структуры аккреционного диска вначале производилось вертикальное сжатие сетки к плоскости симметрии. На втором этапе преобразования координаты центров ячеек сдвигались на некоторое расстояние к поверхности аккретора. Это позволило получить более высокое разрешение в этой области. Особенно это важно для адекватного разрешения структуры магнитного поля вблизи поверхности аккретора. В расчетах использовалась сетка с числом ячеек Все приводимые ниже расчеты проводились на вычислительных кластерах Межведомственного Суперкомпьютерного Центра РАН МВС-6000IM и МВС-15000BM с использованием 128 процессоров.
Разностная сетка Рис 5. Структура расчетной сетки в экваториальной (xy) плоскости.
Разностная сетка Рис 6. Структура расчетной сетки в вертикальной (xz) плоскости.
Структура газодинамического течения Рис. 7. Распределение плотности и скорости в экваториальной плоскости к моменту времени 7.4 P orb в немагнитном случае. Пунктирной линией показана граница полости Роша. Показана линия тока, начинающаяся из внутренней точки Лагранжа.
Структура газодинамического течения Рис. 8. Распределение изолиний плотности в экваториальной плоскости к моменту времени 7.4 P orb в немагнитном случае.
Структура газодинамического течения Рис. 9. Распределение температуры и скорости в экваториальной плоскости к моменту времени 7.4 P orb в немагнитном случае.
Структура газодинамического течения Рис. 10. Распределение плотности и скорости в вертикальной плоскости к моменту времени 7.4 P orb в немагнитном случае. Пунктирной линией показана граница полости Роша.
Структура газодинамического течения Рис. 11. Распределение температуры и скорости в вертикальной плоскости к моменту времени 7.4 P orb в немагнитном случае.
Синтетическая доплеровская томограмма Рис. 12. Распределение логарифма интенсивности рекомбинационного излучения в экваториальной плоскости (слева) и синтетическая доплеровская томограмма (справа). Маркерами обозначены основные элементы течения.
Заключение 1. Для численного моделирования МГД течений в полураздельных двойных системах развит трехмерный параллельный численный код, основанный на конечно-разностной схеме годуновского типа повышенного порядка точности для уравнений магнитной газодинамики в унифицированных переменных. 2. При построении математической модели предполагается, что собственное магнитное поле звезды-аккретора является дипольным. Для уменьшения ошибок при операциях с большими числами в разностной схеме вычисляется только магнитное поле, индуцированное токами в аккреционном диске и во внешней оболочке. 3. Для проверки вычислительных качеств и правильности работы кода проведены тестовые численные расчеты задачи Римана о распаде произвольного разрыва (газодинамический и МГД варианты), задачи об объемно-распределенном взрыве в среде, пронизанной однородным магнитным полем и задачи об аккреции изотермического и адиабатического газа на гравитирующую точечную массу. 4. Проведено численное моделирование структуры течения в системе SS Cyg. В газодинамическом случае за время порядка нескольких орбитальных периодов системы вокруг звезды-аккретора сформировался аккреционный диск с характерными особенностями структуры: ''горячей линией'', приливными ударными волнами, ''прецессионной'' спиральной волной плотности и др.