Символьные и численные методы в математическом моделировании гравитирующей быстровращающейся сверхплотной конфигурации с реалистическими уравнениями состояния Беспалько Евгений Валерьевич

Диссертация - 480 руб., доставка 10 минут, круглосуточно, без выходных и праздников

Автореферат - бесплатно, доставка 10 минут, круглосуточно, без выходных и праздников

Беспалько Евгений Валерьевич. Символьные и численные методы в математическом моделировании гравитирующей быстровращающейся сверхплотной конфигурации с реалистическими уравнениями состояния : Дис. ... канд. физ.-мат. наук : 05.13.18 Тверь, 2005 110 с. РГБ ОД, 61:06-1/271

Содержание к диссертации

Введение

1 Уравнение гидростатического равновесия для гравитирующей стационарно вращающейся, замагниченной конфигурации 13

1.1 Основное уравнение математической модели быстро-вращающейся намагниченной конфигурации 13

1.2 Аналитическое представление вклада давления в уравнение гидростатического равновесия 18

2 Математическая модель гравитирующей быстровра-щающейся сверхплотной конфигурации с реалистиче скими уравнениями состояния 22

2.1 Основные предположения и постановка задачи 22

2.2 Схема расчетов характеристик модели 26

2.3 Обсуждение результатов 34

2.4 Регуляризация метода Ньютона и выбор оптимального итерационного параметра 48

3 Вычисление ньютоновского потенциала гравити рующей конфигурации с поверхностью близкой к сфероиду 53

3.1 Постановка задачи 53

3.2 Метод рядов Бурмана-Лагранжа при выборе параметров возмущенной эллипсоидаль ной поверхности 54

3.3 Ньютоновский гравитационный потенциал на внутреннюю точку Фт(Р, L) 59

3.4 Ньютоновский гравитационный потенциал на внешнюю точку <$>out(N,P,L) 65

3.5 Обсуждение результатов 69

Оценка погрешности решений уравнений описываю щих быстр о вращающуюся гравитирующую конфигу рацию 71

4.1 Оценка погрешности метода 71

4.2 Погрешность аппроксимации поверхности 79

4.3 Метод оценки погрешности в линейном приближении 81

Заключение 83

Приложения 85

Литература

Аналитическое представление вклада давления в уравнение гидростатического равновесия
Схема расчетов характеристик модели
Метод рядов Бурмана-Лагранжа при выборе параметров возмущенной эллипсоидаль ной поверхности
Погрешность аппроксимации поверхности

Введение к работе

Математическое моделирование самогравитирующих систем вещества стало в последнее десятилетие одним из приоритетных направлений в астрофизике. В первую очередь этому способствовали открытия наблюдательной астрономии таких космических объектов как пульсары, квазары и компактные рентгеновские источники (рентгеновские пульсары).

В настоящие время наибольший интерес вызывают пульсары. Пульсар-это вращающаяся намагниченная нейтронная звезда, обладающая осью симметрии, наклонной к оси вращения [1-5]. Пульсары были открыты А. Хъюишом и другими в 1968 году [6]. Проблема фигур равновесия быстровращающейся самогравитирующей жидкой капли долгое время рассматривалась как чисто математическая задача ньютоновской теории гравитации. Изучением данной проблемы занимались такие выдающиеся математики, как К.Г. Яко-би, A.M. Ляпунов и многие другие.

В конце двадцатого века С. Чендрасекхар и другие впервые исследовали постньютоновские поправки и реакцию гравитационного излучения на фигуру равновесия [7, 8]. Открытие пульсаров перевело эту проблему в практическую плоскость. Доминирующими факторами в динамике и эволюции пульсаров является сильное собственное гравитационное поле, быстрое вращение, частота которого достигает нескольких сотен оборотов в секунду, а также большое давление в центре, поэтому любая реалистическая модель таких объектов должна основываться на общей теории относительности

[9, 10, 11, 12].

Актуальность математического моделирования сверхплотных конфигураций в значительной степени обусловлена проблемой регистрации гравитационных волн, поскольку отклонения фигуры звезды от осевой симметрии, вызываемые, например, магнитным полем в сочетании с эффектами вековой неустойчивости, делают ее источником монохроматического гравитационного излучения. Интерес к этой проблеме резко возрос после открытия в 80-х и 90-х годах большого числа миллисекундных пульсаров [13, 14, 15, 16], так как интенсивность гравитационного излучения пропорциональна шестой степени частоты вращения. Важным является то, что уверенный прием гравитационных волн от пульсаров будет, по сути, первым шагом на пути к созданию гравитационно-волновой астрономии, открывая новый канал получения информации из космоса, наряду с электромагнитным и нейтринным каналами, действующими в настоящее время [17, 18].

Отметим еще одно важнейшее направление в исследовании пульсаров, имеющее конкретные перспективы технического применения. Речь идет о создании уникального стандарта частоты и шкалы времени между импульсами радиоизлучения пульсаров [19, 20, 21]. Подобные "часы^,,лишены основного недостатка квантовых часов любого типа, который заключается в непредсказуемом накоплении ошибки при измерении больших промежутков времени, что имеет принципиальное значение для повышения астрономических наблюдений и космической навигации даже в пределах солнечной системы; проект создания пульсарной шкалы времени, предусматривающий, в частности, вывод на орбиту Земли радиолокационной антенны для приема импульсов, начал прорабатываться в нашей стране и ряде других стран в 1989 году, но был отложен по соображениям ненаучного характера на неопределенный срок. Совершенно очевидно,

что в научном обосновании этого проекта проблема конфигурации пульсаров играет очень важную роль.

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

Примерами таких методов являются разложение Клеро-Лежандра, разложение Чандрасекара-Милна и метод квазисферической аппроксимации. Последний еще называют методом двойной аппроксимации: 1) во внутреннем ядре, где центробежная сила всегда мала, применяется разложение первого порядка по параметру v 2) пренебрегается влиянием массы внешних слоев, считая, что сила тяготения порождается только веществом слегка сплюснутого ядра. Главное преимущество этого метода состоит в том, что вращение без особых затруднений удается включить в обычные программы расчета эволюции звезд.

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

Метод заключается в том, что при помощи уравнения:

ф = ф(р;Д (1)

где j-заданная функция от ш, а ш- доля массы заключенного в цилиндре, задав подходящее пробное распределение плотности

po(j,z), найдем функцию Фо(,.г).

Исходя из этой функции и учитывая уравнение:

Р = Р(Ф), (2)

можно в свою очередь получить уточненное распределение плотности pi(Q_yz). Подставляя эту плотность в уравнение (1), получим уточненный потенциал Фі(,,г) и т.д. Таким образом, попеременно решая уравнения (2) и (1), мы придем к согласованному решению.

Следует отметить, что для уравнения поверхности, сильно отличающегося от сферы, есть и другие подходы - чисто разностная схема, вариационные методы и.д. [22, 23, 24].

В работе [25] была сделана попытка исследовать структуру газовых масс на примере политроп и случая белых карликов численными методами с использованием компьютерных методов.

Политропному случаю соответствует уравнение состояния:

Р = Кр¹⁺і,

где n-политропный индекс, ЙТ-константа пропорциональности зависит от величины энтропии, приходящей на нуклон, и от химического состава, но не зависит от г и ^(центральная плотность).

В самом общем случае электронное давление р_е в белом карлике зависит от плотности р, температуры Т и химического состава. Однако электронный газ в основном объеме белого карлика столь сильно вырожден, что даже при довольно высоких температурах (скажем, Т « 10⁷ К) в большинстве случаев прекрасным приближением является полное вырождение (Т = О К) по крайней мере в том, что касается глобального внутреннего строения звезды. Другими словами, в первом приближении белый карлик можно рассматривать как баротропу, следовательно:

P = af(x), *=()*,

где f(x) = х(2х³ - 3)(х² + 1)з + ШпН-^гх, 6 = 9- 10%, a \i_e-молекулярный вес электрона. Необходимо отметить, что холодный полностью вырожденный белый карлик можно рассматривать как политропную конфигурацию в предельных случаях низкой плотности (га = 1, 5) и высокой (га = 3) [22, 23, 26].

Структура конфигурации в [25] определяется теоремой Гаусса:

<9Ф др,

_{+ 4(1}

= -4жвр,

(1 - М²)

11 h) ₊1 д_

г² дг \ дг ) г² dfx вместе с уравнениями гидростатического равновесия:

_=р_ ₊ ^г(1-л³)

ЭР дЧ ₂₂

дР__ «ЭФ

здесь Ф-гравитационный потенциал, ^-расстояние от центра масс конфигурации. То есть, чтобы найти строение звезды около центра, плотность р и гравитационный потенциал Ф разлагаются в степенные ряды по радиальной переменной. Коэффициенты в этих разложениях в свою очередь разлагаются по многочленам Лежандра Pi(cos0). Аналитическое продолжение, а затем пошаговое интегрирование описывают строение внешних областей. Физические параметры твердотельно вращающихся политроп найдены в диапазоне 0 < п < 3. При п > 3 метод Джеймса принципиально не применим [22, 23, 25], т.к. при га > 3 становится все труднее хоть с какой-нибудь разумной степенью точности определить внешнюю границу.

Джеймс показал, при п < 0,808 на каждой последовательности осесимметричных твердотельно вращающихся политроп имеется точка бифуркации, в которой ответвляется неосесимметрич-ные фигуры равновесия. Если га > 0,808, то на последовательности

твердотельно вращающихся политроп бифуркации нет. Реальные конфигурации имеют реалистические уравнения состояния: Бете-Джонсона, Рейда.

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

Первая область ранр < Р < Рпис & 2.8 Ю¹⁴?^³" ^сР^авнительно хорошо изучена. Равновесная материя состоит из обогащенных нейтронами ядер, образующих кулоновскую решетку, электронов и свободных нейтронов. При возрастании плотности свободные нейтроны обеспечивают все большую долю полного давления. При р

*~*^ Рпис

начинается деформация и разрушение ядер, т.е. ядра начинают распадаться и сливаться.

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

При сверхвысоких плотностях, р > Ю¹²^з в материи появляется заметное количество гиперонов, и взаимодействие между нуклонами должно рассматриваться с учетом релятивистских эффектов. Надо отметить, что уравнения состояния, полученные до настоящего времени, содержат множество неопределенностей.

В данной работе будут использованы три уравнения состояния:

1.Уравнение состояния Бете-Джонсона(ВЛ) описывает состояния конденсированного вещества при 1.7-10¹⁴р < 3.2-10¹⁶, т.е. при сверхвысоких плотностях. Предполагается, что материя содержит нейтроны, протоны и гипероны с массами, не превышающими массу А-резонанса, взаимодействие между которыми описывается моди-

фицированным потенциалом Рейда.

2.Уравнение состояния Рейда(Ы) описывает состояния конденсированного вещества при р > 7 10¹⁴. Причем основной компонент вещества - нейтроны, взаимодействие между которыми описывается потенциалом Рейда с мягким кором, приспособленным к ядерной материи.

Для сравнения будет рассмотрено уравнение состояния Оппенге-ймера-Волкова (OV), описывающее состояния конденсированного вещества при 0 < р < оо, причем, основной компонент вещества -нейтроны в этом случае не взаимодействуют между собой, т.е. представляют идеальный невырожденный ферми газ.

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

Центральная задачи при расчете конфигурации - это задача вычисления ньютоновского гравитационного потенциала на внутреннюю точку.

Эта задача представляет также самостоятельный интерес, а не только в связи с расчетом конфигурации.

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

потенциал входит в уравнения, определяющие конфигурацию, и в уравнения, описывающие динамику гравитирующих масс.

В теории ньютоновского гравитационного потенциала возникает три типа задач по трем типам симметрии. Наиболее простой случай, когда конфигурация не вращается; тогда она будет иметь сферически симметричную форму. Этот случай наиболее хорошо изучен [27]. Но часто необходимо учитывать вращение, тогда конфигурация будет иметь сферически несимметричную форму [23]. В этом случае мы приходим к необходимости использовать более сложные поверхности, в частности эллипсоидальные, граница которых представляет эллипсоид [28, 29]. Изучением эллипсоидальных фигур равновесия занимались: Маклорен, Якоби, Ляпунов и многие другие.

Разработаны различные методы для нахождения потенциала как на внутреннюю точку, так и на внешнюю. Наиболее известный метод вычисления потенциала эллипсоида - это метод Лагранжа. Но вычисление внешнего потенциала эллипсоида методом Лагранжа представляет известные трудности, главная из которых состоит в том, что область интегрирования зависит от положения притягиваемой точки в пространстве; когда притягиваемая точка находится внутри притягивающего объема, область интегрирования охватывает всю поверхность притягивающего объема.

Для вычисления внешнего потенциала эллипсоида можно воспользоваться теоремой Айвори, в которой говорится о том, что проекции сил притяжения, действующие на одну и ту же точку со стороны двух софокусных эллипсоидов (Е) и (") (причем точка лежит вне эллипсоида (Е), а эллипсоид (") проходит через эту точку) на одну и ту же главную ось соотносятся между собой, как площади сечений эллипсоидов (Е) и (Е') плоскостью, перпендикулярной к выбранной главной оси. Таким образом, данная теорема позволяет

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

Однако в реальных случаях для неоднородных конфигураций возникает задача, в которой поверхность конфигурации представляет из себя более сложную структуру. В работе [30] предложен метод аппроксимации поверхности псевдоповерхностью, а именно возмущенной эллипсоидальной поверхностью, параметры которой определяются из условия минимума квадрата плотности на этой поверхности. Для этих целей эффективно может быть использован метод разложения в ряд Бурмана-Лагранжа по малому параметру.

Надо отметить, что теория потенциала является бурно развивающимся разделом математической физики. И ее результаты используются далеко за пределами породивших ее астрономии и геофизики - в частности во многих разделах чистой математики. И опубликовано много работ, трактующих теорию потенциала с чисто математических позиций [31, 32, 33].

Ввиду чрезвычайной сложности как аналитических, так и численных расчетов возникает необходимость использования компьютерных методов. Мы воспользовались пакетом символьной и численной математики MAPLE.

Система компьютерной математики MAPLE была выбрана не случайно. В рамках этой системы можно быстро и эффективно выполнять не только символьные, но и численные расчеты, при-

чем это сочетается с превосходным средством графической визуализации и подготовки электронных документов. MAPLE является одним из лидеров среди подобных себе систем, не зря ядро MAPLE используется в ряде других математических систем, например, MATLAB и Mathcard, для реализации в них символьных вычислений [34, 35, 36, 37].

MAPLE- типичная интегрированная система. Она объединяет в себе:

-мощный язык программирования

-редактор для подготовки и редактирования программ

-ядро алгоритмов и правил преобразования математических выражений

-численный и символьный процессор

-систему диагностики

-библиотеки встроенных и дополнительных функций

-пакеты функций сторонних производителей и поддержку некоторых других языков программирования и программ.

Необходимо отметить, что последние реализации MAPLE являются одними из самих надежных систем компьютерной математики. Надежных прежде всего в смысле высокой достоверности получения правильных результатов при сложных символьных вычислениях. Это первая система компьютерной математики, успешно прошедшая тестирование на задачах повышенной сложности, предлагаемых для оценки качества подобных систем.

Целью исследования данной диссертации является построение математической модели сверхплотной конфигурации с реалистическими уравнениями состояния на основе символьных и численных методов.

Аналитическое представление вклада давления в уравнение гидростатического равновесия

Как известно, уравнение состояния ядерного вещества связывает давление Р, плотность, температуру Т и химический состав звезды. Символически это можно записать следующим образом: Р = Р(р,Г,АьА2,...).

Однако в случае холодных белых карликов или нейтронных звезд температура фактически всюду равна нулю и следовательно всюду s = 0. В силу этого уравнение состояния можно взять в виде р = Р(р).

В данной работе будут использованы три уравнения состояния ядерной материи Бете-Джонсона (BJ), Оппенгеймера-Волкова (OV) и Рейда (R) [22, 51].

Уравнения состояния ядерного вещества очень часто в литературе представлены в виде графиков. Как следствие, нам понадобится процедура оцифровки. Поэтому возникает задача об аналитическом представлении уравнений состояния и как следствие р интеграла J — с последующим использованием его в уравнении о гидростатического равновесия. Будем считать, что Р(р) - функция гладкая и монотонно возрастающая. Интеграл из уравнения гидростатического равновесия удобно искать в виде Ф(у), где р р Є [О, А)], У Є [-1,0], р = ро(1 + У)- Причем у = Y, РаЬс Х&з, abc а + Ъ + Р, более подробно мы разберем это ниже. Таким образом, получаем: (У pj р Pj Уро(1 + у) PjV dP(y) (1 + 2/) 0 -1 -1 где PQ = Р(ро) Продифференцировав Ф(у) по у получаем: (1.2.1) №_ _ 1 1 dP dy Pl + y dy откуда, у = (1 + у)Щу)- JvV{y)dy. -і Надо отметить тот факт, что, поделив Р(р) на Р0, а р на ро, мы ввели безразмерные величины и -. Обозначим— через Ру. Для нахождения Ф(у) необходимо решить задачу N J2(Pi - Hto)? - rnin, (1.2.2) г=0 где Pi = - j У І = 15 причем, Pi и pi получены посредствам оцифровки. Будем искать Ф(у) в виде многочлена: л. %) = EV- (1-2.3) Тогда, проинтегрировав (1.2.3), получаем: НУ) = (у + Л 4- (1-2-4) Продифференцировав (1.2.2) по 5k,к = 0..п и учитывая (1.2.4), получим систему уравнений: ( + пйг )=. fc = 0-"- (L2-5)

При А; = 0 получим тождество 0 — 0. Поэтому имеем систему из п уравнений, из которой найдем 5k- Для этого систему уравнений (1.2.5) перепишем в виде: п J2Ab& = Bk, к = 1..п, (1.2.6) 3=0 где

Коэффициент 50 находится из условия Ф(—1) = 0. Учитывая вышесказанное, построен алгоритм для аналитического представления вклада давления в уравнение гидростатического равновесия (1.1.14).

1. Воспользовавшись программой оцифровки Grafula, создаем таблицу значений pi и PPi-,i = 0..iV. Или таблицу логарифмов pi и PPi, если график представлен в дважды логарифмических осях.

2. Если график дан в дважды логарифмических осях, полученную таблицу логарифмов pi и PPi переводим в таблицу значений pi иРРі) і = 0..N и получаем безразмерные величины / — -р -, f- и вычисляем УІ = Я1 — 1, і = 0..N. Уї РО 3. Воспользовавшись системой символьной математики MAPLE [52-57] для решения системы (1.2.6), находим &. 4. Вычисляем погрешность аппроксимации по формуле d = 5Z»=o- (- " (&))2- В силу технических причин пришлось ограни читься случаем N = 2, погрешность в этом случае составила 5-Ю-3. лл

В предыдущей главе мы показали, что уравнение гидростатического равновесия для гравитирующеи, стационарно вращающейся, замагниченной конфигурации в постньютоновском приближении удобно записать в виде [58-67,23]: Ф + (/ ) - у ( 2 + У2) = П(РЛг) + IU (2.1.1) р ад=/ - (2.1.1-) где Ф-гравитационный потенциал конфигурации, P(/J)-давление, П(р )-вклад ностньютоновских поправок в уравнение гидростатического равновесия, Пт-вклад магнитных напряжений- Конкретный вид функции Р(р) зависит от уравнения состояния ядерного вещества конфигурации. Напомним, что для расчетов мы будем использовать численные данные для уравнения состояния ядерной материи Бете-Джонсона, Оппенгеймера-Волкова, Рейда [22]- Таким образом, уравнение (2-1-1) представляет интегральное уравнение с подвижной границей в R3 относительно / .

Как будет показано в данной работе, параметры конфигурации зависят от выбора уравнения состояния, что является критерием для их отбора. В частности, при критическом значении е = е& происходит резкое увеличение интенсивности гравитационного излучения (ГИ) конфигурации, при регистрации частоты которого определяется б& и появляются аргументы в пользу выбора конкретного уравнения состояния.

Схема расчетов характеристик модели

Вычисление ФаЬс ДОЯ значений Р = 4, 6,8 и L = 2,4,6 нереально без использования системы символьных вычислений на компьютере. С этой целью нами была составлена программа в системе символьной математики MAPLE для нахождения Фаьс Поскольку программа оказалась большая и сложная, возникает необходимость ее тестирования- Нами были выполнены два теста. Проводился предельный переход е 1к сферически симметричному случаю, который легко рассчитывается по независимой программе- Кроме того ДФ = AnGpj т.е. действуя на Ф оператором Лапласа, мы получаем плотность умноженную на AnGp. Оба теста прошли успешно, для различных значений L и Р. Аналогичная программа в системе MAPLE была составлена и для функционала А.

Представим 0 и П(т) также в виде многочленов степени Р по координатам жі, x2j х$ р г а,Ь}с atb,c В2 SirGpQaf где Бо-характерная напряженность магнитного поля в центре конфигурации. Тогда система уравнений (2ЛЛ),(2Л-6) может быть записана в виде системы алгебраических уравнений относительно раьс и Z \ —ФаЬс + K0Qabc - c(62aSbO + 2b ao) c0 = Щт)аЬс, (2.2Л8) Фу = 0, Ф = о, где к Р ш"

Коэффициенты, определяющие структуру конфигурации раЬс и Zijfa разобьем на симметричные и антисимметричные части относительно оси вращения в виде разложения по степеням параметра X до квадратичного члена включительно: РоЬс = /а\( /М,АН-Ь,с + Pl{ab)cX + Р[аЬ]СХ (2-2Л9) Ziik — Т7ТТТТГ7 Ї-М.АІ + Z\U4\kX + Z[ij]fcX Здесь и далее a&c и ijk являются четными, а вводимые величины удовлетворяют соотношениям симметрии: Рі(аЬ)с = р1(Ьа)с, Р[аЬ}с = —р[Ьа]с (2.2,20) %l(ij)k = %l(ji)ki %[ij]k = —Z\ji]k Симметризируя и антисимметризируя по первым двум индексам систему уравнений (2.2.18), получаем новую систему уравнений для

Определения р{аЬ)с, Р[аЬ]с, Z(ab)c, Z[ab\c -Ф(аЬ)с + Ко(аЬ)с - «(UaUo + йЛо)йю = 0 (2.2,20a) — Ф[аЬ]с+ 0[а6]с = ктЩт)[аЬ]с (2.2.20&)

При выводе (2.2.20а) мы отбросили малый члеп Щт) {аЬ)сі кото-рый приводит только к малой поправке порядка fem.

Коэффициентами в уравнениях (2.2.20а), (2.2.20Ь) являются интегралы Jac% которые мы будем вычислять численно. Вычисление этих интегралов и все другие расчеты будем проводить в формате ю-30.

С физической точки зрения є свободный параметр, а е вычисляемый. Однако из удобства вычислений проще считать е свободным параметром, а є вычисляемым- В результате все параметры конфигурации будут зависеть от е,/?о,.Ро Подставляя (2.2,19) в (2.2.18), получаем систему уравнений представляющую многочлен по степеням малого параметра X. Поэтому возможно использовать метод разложения по малому параметру X. В первом приближении положим X = 0 и найдем значения рас и 2 , соответствующие фигуре вращения- В этом приближении положим 7 0 — #оо є — 6о- Двумерные массивы неизвестных в системе (20а) Zik, рас. Ко = -Кось с = б0 обозначим, как ут{т = 1,2,-iVi).

Легко найти формулу для определения Ni : N, = 1(Р + 2)(Р + 4) + l(L + 2){L + 4). (2.2.21) При Р = 6 и L = 2 имеем iVi — із. Связь между переменными будет выглядеть следующим образом: Ш = Р02, Уг = PQ4, УЗ = А», У4 - Р20і Ї/5 = Р22: г/6 = р24, г/7 = 0. 2/8 = 42, г/9 = РбСЬ У10 — 02 Ун = z20, г/і2 = о, г/із = е В этом случае уравнения (2.2,20а) представляют собой систему из 13 нелинейных алгебраических уравнений и могут быть записаны в векторном виде: Ду,е) = 0, У = (уі,г/2,.-.,ї/із). (2,2,22) Вследствие плохой обусловленности матрицы Якоби / (у, е) для численного решения (2-2.22) будем использовать регуляризованный аналог метода Ньютона [80] с параметром регуляризации а 10_6-Имеем следующую итерационную схему: У(п+1)(е) = УИ(е) - rn[a/2(y (e),e) + 7( ( ), /(у И.е)]- xT(yin\e),e)f(yW(e),e), (2.2.23) где тг-номер итерации rn-(G0 rn С 1) [80], / (у (еХе) -матрица Якоби, / (у (е)?е) " транспонированная матрица Якоби. Величина //2(у( )(е), є) — w(e) представляет собой невязку и определяет точность решения системы уравнений (2.2.22), которая в нашем случае составила 10 30.

Метод рядов Бурмана-Лагранжа при выборе параметров возмущенной эллипсоидаль ной поверхности

Как отмечалось во второй главе, граница конфигурации Е находится из условия равенства плотности на границе нулю: p(x,y,z) = 01 (x,y,z)eZ. (3.2.2) Потенциал Ф будет удовлетворять на границе условиям: Ф = Фіпв и (УФоиОв = ФІП)Е. Поверхность , как и ранее, выберем в виде возмущенной эллипсоидальной поверхности: SD: ЗЯ + ХІ+ХІ + 2 ад - 1. (3.2.3) Напомним, что параметр сплюснутости е — . Условия близости (3.2.2) и (3.2.3) можно сформулировать введением функционала Л: Ц /Vdfi. (3.2.4) 4тгр,

В отличие от ранее приведенного условия близости (2.1.5) здесь будет отсутствовать весовой множитель -\. В данной постановке задачи нам не нужно, чтобы интегралы Л и Ф имели одинаковую структуру, так как система, к которой приводит минимум Л, будет решаться самостоятельно.

Условие минимума Л приводит к системе алгебраических уравнений относительно аьаз и Z : « = 4- = 0 (3-2-5) Фі = 0! -А( = 0) = 0 ф2 = аз -Л( ь = 0) = 0. Представим плотность конфигурации /), как и ранее, в виде полинома степени Р: р а,Ь,с После перехода к сферическим координатам R, 9, ф: xt — Rotk i — sin в cos ф, 6І2 = sin в sin ф, a3 = cos в, выражение для возмущенной эллипсоидальной поверхности (3.2.3) примет вид: R = 1 Y, ZiJk Y&\4ak3 = 1 + Ф(Д, a , )- (3.2.7) Дальнейшие вычисления будут основаны на использовании варианта теоремы Лагранжа, из которой мы получаем: Rh+2 в 2 + (л + 2) ± [а ЩаУ] . (3.2.9а) IL LLCL 3=1 Более детально (3.2.9а) будет выглядеть следующим образом: RM=а +(л+2) J: L Ztjkm&i&i± s=l ijk (3.2.9b) где Zyt(s)-MHOr04JieH от степени, не превосходящей 5. Рассмотрим более подробно производные в выражении (3.2.9Ь): — _ __ = (произведем заменуа = у + 1, (у-близко к нулю) = _ 1 d -1 (i + y) +j+ +fc+i _ і 27 dy -1 (i + )s ;l = 5 s Для определения явного вида Фя понадобится соотношение, полученное в [68]: (1 + у)" 1{а) Ф8 = S{a, Ь) = (1 + !) +!_ = П(6 + 1-2г), 2/=0 г=1 где а = s — l,ft = f + j + fc + A + L Из него и (3.2.9Ь) следует: 5—эо sL Д = Е Е 2чк(а)с\с&с$Ф9(% +з + + Л - 1), (3.2.9с), 5! S=0 tjft где 0(s) — 0 при s = 0 и &(s) = 1 при s О. Ф3 — 1 при я Ои s = l. Квадрат плотности будет выглядеть следующим образом: р р А = Е Е /WaW . 3-2.10) обе а 6 с Так как мы перешли к сферическим координатам, то, учитывая (4,2.9с), выражение для квадрата плотности примет вид: Z,a+a/+i b+br+j c+c +k 5! s=0 abc a b d ijk со Р Р sL ґ- y is) Рр = 2 2 S S Г PabcPalVc Zi3k{s)& ФДЛі + ї + j + A-l), (3.2.11) где /z-i — a + a! + 6 + br + с + /. Учитывая (4-2.4) и (4.2.11), Л можно представить в явном виде, удерживая в разложении Л члены Z k степени не выше sm, где $т-максимально учитываемая степень по з : А, , Л (-i) fteM AlsmJ = 2 L Z t 2 1 РаЬсРа Ъ с Лф abc,a b c s—Q ijk Ф (Лі + і + j + fc - 1)QAI,4,H3 (3.2.12) QA, ,A3 = Ja d a dQ, (3.2.13) где Ai = a + a + t, A2 = b + У + j, Ai = с + d + /с. Интегралы (3.2.13) могут быть упрощены: Qa+a +i,b+b +j,c+ +k = (3.2.14) (a + a + z -l)!!(b + y + j-1)11(0 + (/ +fc-1)!! (Лі + г+і + А:)!! Трехмерные массивы неизвестных в системе (3-2.5) %цъ обозначим через уш{та = 1, 2, .,iVi), При Р — 6, L — 2 имеем Ni — 5. Связь между переменными будет выглядеть следующим образом: Уі = у, V2 = Є, Уз = 002, УА = 020, Уъ — 200, где I- характерный масштаб распределения масс гравитирующей системы. В этом случае систему уравнений (3.2.5) удобно записать в векторном виде: /(у,е)=0, (уъУ2,...у,) (3.2.15)

Для численного решения системы (3.2.15) использован регуляризо-ванный аналог метода Ньютона [80] с параметром регуляризации а = 1(Г3: +,)(,, = - ( ) ) + 7( 4 ),6)/ (6),6)1- x7(y{n)(e),e)f(yW(e),e), (3.2.16) где тп-(во тп 1)-оптимальный шаг на п-ой итерации, причем из соображений сходимости во ОЛ.

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

Погрешность аппроксимации поверхности

Как говорилось ранее, точную границу конфигурации мы заменили псевдограпицей 8D9 а условия близости точной границы и псевдограницы сформулировали введением функционала Л (2-1.5). Псевдограницу мы представили в виде возмущенной эллипсоидальной поверхности, т.е. с добавкой многочлена по координатам степени

В данном разделе будут приведены графики, иллюстрирующие отклонение плотности от 0 на псевдоповерхности конфигурации SD при различных L для полученных в диссертации значений р.

Из графиков можно сделать вывод, что во всех трех случаях уравнения состояния при увеличении L на 2, то есть при L — 4, точность аппроксимации точной поверхности конфигурации увели-чивается на несколько порядков. Эти результаты дают нам возможность предположить, что при увеличении L в несколько раз мы можем добиться аппроксимации поверхности с высокой степенью точности.

В диссертационной работе проведено исследование построенной математической модели гравитирующей, быстровращающейся сверхплотной конфигурации с использованием уравнений состояния ядерной материи: Бете-Джонсона, Оппенгеймер-Волкова и Рейда- Разработанная математическая модель позволяет исследовать эволюцию пульсаров, их гравитационное излучение, а также, что представляет наибольший интерес, гравитационное излучение вблизи критических точек. Основной целью данной работы было численное доказательство существования критических точек(точек бифуркации) по параметрам є и є решений уравнения гидростатического равновесия рассматриваемой конфигурации. В результате проведенных исследований были получены следующие результаты:

1. Разработана и реализована программа в системе символьной математики MAPLE для аналитического представления ньютоновского гравитационного потенциала неоднородной возмущенной эллипсоидальной конфигурации на внутреннюю точку в виде полинома координат для заданных значений Р и L.

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

3. Уравнение гидростатического равновесия, описывающее бы-стровращающуюся сверхплотную замагниченную конфигурацию, сведено к системе нелинейных алгебраических уравнений для коэффициентов разложения плотности ра&с, полуосей а Оз и коэффициентов возмущенной эллипсоидальной поверхности Zijk. Полученная система представлена в виде разложения по степеням малого параметра X = р[щ0 = (/ 200 _ Аш)- В нулевом по X приближении симметричная система решена регуляризованным аналогом метода Ньютона, Составлена соответствующая программа в пакете MAPLE- Найдены численные значения коэффициентов разложения плотности и параметры, определяющие форму аппроксимирующей псевдоповерхности, как функции параметра е. Вблизи критических точек исследуемая система уравнений найдена в приближении с точностью до членов X3 и сведена к решению соответствующего кубического уравнения.

4. В линейном по параметру X приближении с помощью символьных и численных расчетов доказано существование точек бифуркации по параметрам е и є решений уравнения гидростатического равновесия конфигурации, в которых происходит ответвление асимметричных решений относительно оси вращения для распределения плотности при Г)т — 0.

5. Проведено исследование этих решений при г)т — 0, которое показало, что при е е (є Єь) кубическое уравнение имеет одно вещественное решение X — 0, что соответствует аксиально-симметричной конфигурации, Но при е е&( 6k) уже будут два корпя X = 0 (симметричное решение) и X = \/%{бк — е)(асимметричное решение).

б- Проведено исследование полученного кубического уравнения при rjm Ф 0, которое показало, что при е& — е 3 т/m кубическое уравнение имеет три решения, а при efc — е С r#i оно имеет одно решение. В точке бифуркации е = е&, Хк = X{ek) = (1) э те-асимметрия распределения плотности на много порядков больше, чем в линейном приближении.