close

Вход

Забыли?

вход по аккаунту

?

Интегрирование уравнений движения малых тел Солнечной системы методом оскулирующих элементов.

код для вставкиСкачать
Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки. — 2009. — № 1 (18). — С. 222–227
Небесная механика и астрономия
УДК 523.642
ИНТЕГРИРОВАНИЕ УРАВНЕНИЙ ДВИЖЕНИЯ МАЛЫХ ТЕЛ
СОЛНЕЧНОЙ СИСТЕМЫ МЕТОДОМ ОСКУЛИРУЮЩИХ
ЭЛЕМЕНТОВ
А. Ф. Заусаев, Д. А. Заусаев
Самарский государственный технический университет,
443100, Самара, ул. Молодогвардейская, 244.
E-mails: Zausaev_AF@mail.ru
Предложен алгоритм метода оскулирующих элементов для численного интегрирования уравнений движения малых тел Солнечной системы.
Ключевые слова: интерполяция, численное интегрирование, дифференциальные
уравнения движения, метод оскулирующих элементов.
Дифференциальные уравнения движения. В связи с решением проблемы
астероидной опасности исследование эволюции орбит малых тел Солнечной
системы (астероидов, комет, крупных фрагментов в метеорных потоках) является актуальной задачей.
Большое значение имеет выбор математической модели, описывающей
движение небесных тел. Учет лишь гравитационных взаимодействий, который заложен в стандартных уравнениях движения задачи n тел, оказывается
недостаточным. Это приводит к более сложной форме дифференциальных
уравнений движения, где наряду с гравитационными эффектами учитываются релятивистские эффекты, фигура планет и Солнца, реактивные силы
и т. д.
В работе [1] предложен интерполяционный метод оскулирующих элементов, используемый для получения координат и элементов орбит больших планет и малых тел Солнечной системы. Этот метод предназначен для работы с
банком данных координат и элементов орбит малых тел Солнечной системы и
обеспечивает требуемую точность на ограниченном интервале времени — порядка 100 дней. Данное ограничение связано с тем, что смещение исследуемого объекта на каждом шаге интегрирования проводилось по невозмущенным
кепплеровским орбитам. Возмущение от больших планет учитывалось только в конце шага интегрирования. Однако при использовании данного метода
для численного интегрирования уравнений движения на интервале времени
порядка нескольких сотен лет следует учитывать изменение элементов орбит
Заусаев Анатолий Федорович — профессор кафедры прикладной математики и информатики; д.ф.-м.н., проф.
Заусаев Дмитрий Анатольевич — студент специальности «Прикладная математика и
информатика».
222
Интегрирование уравнений движения малых тел Солнечной системы . . .
на каждом шаге интегрирования. Для решения поставленной задачи необходимо из уравнений движения, приведенных в работе [1], получить дифференциальные уравнения движения для задачи двух тел и проинтегрировать
их.
Рассмотрим сначала абсолютное движение двух тел относительно произвольной инерционной системы отсчета. Обозначим через ρ0 и ρ векторы тел
S и M , определяющие их положения относительно начала координат O выбранной инерционной системы координат. Положение тела M относительно
S будет определяться вектором r = ρ − ρ0 . Тогда дифференциальные уравнения движения, приведенные в работе [1] в векторной форме, определяющие
движение тел S и M , запишутся так:
2
d2 ρ0
r
3aom rom
p
p
· ,
=
2
3
3
2
3
3
3
3
2
dt
r + r r − rom + (r − rom ) r
r
2
d2 ρ
3aos ros
p
p
·
−
,
=
3 + 3 (r 3 − r 3 )2
dt2
r
r 2 + r 3 r 3 − ros
os
(1)
где rom , ros — эффективные радиусы тел M и S; aom , aos — соответствующие
ускорения на расстоянии rom и ros , r — расстояние между центрами S и M .
Вычитая из второго уравнения (1) первое, получим уравнение, определяющее движение тела M относительно S:
d2 r
=−
dt2
2
3aos ros
p
p
3 + 3 (r 3 − r 3 )2
r 2 + r 3 r 3 − ros
os
2
3aom rom
p
p
+
3 + 3 (r 3 − r 3 )2
r 2 + r 3 r 3 − rom
om
!
r
· . (2)
r
Полагая, что rom и aom тела M пренебрежимо малы по сравнению с ros и
aos тела S, что справедливо, если в качестве тела S рассматривать Солнце
или большую планету, а в качестве тела M — астероид или комету, получим
уравнение движения тела M относительно S в следующей форме:
2
d2 r
r
3aos ros
p
p
· .
=
−
2
3
3
2
3
3
3
3
2
dt
r + r r − ros + (r − ros ) r
(3)
r3
Ограничиваясь членами первого порядка относительно r03 в разложениях
1 2
r3 3
r3 3
и 1 − r03
в бесконечный ряд, уравнения (3) можно представить
1 − r03
в виде
2 r 3 d2 r
ao ros
ros
=−
1+ 3 .
(4)
dt2
r3
3r
2 = µ;
Полагая, что ao ros
3
ros
3
= ε, запишем уравнение (4) следующим образом:
d2 r
µr ε
=
1
+
.
dt2
r3
r3
(5)
223
З а у с а е в А. Ф., З а у с а е в Д. А.
Как показано в работе [3], из уравнений (5) следует, что возмущения в
большой полуоси, эксцентриситете и в движении линии апсид можно вычислить по следующим формулам:
2eε
da
=− 2
sin ϕ(1 + e cos ϕ)3 ,
dϕ
a (1 − e2 )5
de
ε
=− 3
sin ϕ(1 + e cos ϕ)3 ,
dϕ
a (1 − e2 )3
ε
dω
= 3
cos ϕ(1 + e cos ϕ)3 ,
dϕ
ea (1 − e2 )3
(6)
где a — большая полуось, e — эксцентриситет, ϕ — истинная аномалия возмущаемого тела, ω — аргумент перигелия.
Изменение элементов орбит за время одного полного оборота находится
интегрированием по полярному углу ϕ (истинная аномалия).
Легко убедиться в том, что большая полуось и эксцентриситет не имеют
вековых изменений, в то время как линия апсид в течении каждого оборота
поворачивается в прямом направлении на угол
Z 2π
ε
∆ω = 3
cos ϕ(1 + e cos ϕ)3 dϕ,
(7)
ea (1 − e2 )3 0
величина которого с учётом членов первого порядка относительно эксцентриситета такова:
3πε
.
∆ω = 3
a (1 − e2 )3
При соответствующем выборе постоянной ε можно получить удовлетворительное количественное объяснение невязок в движении линии апсид планетных орбит.
Таким образом, эффект смещения линии апсид в прямом направлении,
как следует из (7), должен наблюдаться в невозмущенной задаче двух тел.
При численном интегрировании уравнений движения малых тел Солнечной
системы методом оскулирующих элементов начальными данными являются
координаты x, y, z и скорости ẋ, ẏ, ż исследуемого объекта, а также координаты и скорости больших планет xi , yi , zi ; ẋi , ẏi , żi в экваториальной системе
координат на момент времени t0 .
Нахождение элементов орбит по координатам и скоростям. По известным
координатам и скоростям находим элементы орбит по формулам [4]
1 2 v 2 (8)
=
− ,
a r
µ
p
где a — большая полуось орбиты; r = x2 + y 2 + z 2 , v 2 = ẋ2 + ẏ 2 + ż 2 , µ =
2 .
= a0 ros
Эксцентриситет e и эксцентрическая аномалия E на момент времени t0
вычисляются из соотношений
r ṙ
e sin E = √ ,
µ a
224
r
e cos E = 1 − ,
a
(9)
Интегрирование уравнений движения малых тел Солнечной системы . . .
где r ṙ = xẋ + y ẏ + z ż. Затем с помощью уравнения Кеплера
M0 = E − e sin E
вычисляется средняя аномалия в эпоху t0 .
Из соотношений
√
µp sin i sin Ω = y ż − z ẏ,
√
µp sin i sin cos Ω = xż − z ẋ,
√
µp cos i = xẏ − y ẋ,
(10)
(11)
находятся наклонение i и долгота восходящего узла Ω, отнесенные к экватору,
а также параметр p.
Аргумент перигелия ω находится по формулам
√
pr ṙ
z cos i
,
tg u =
,
tg v =
k(p − r)
x cos Ω + y sin Ω
ω = u − v,
ω = ω + ∆ω.
(12)
Вычисление координат и скоростей по элементам орбит. Экваториальные
координаты x, y, z и скорости ẋ, ẏ, ż на следующем шаге h = t − t0 вычисляются с помощью следующих формул [4]:
M = n(t − t0 ) + M0 ,
E − e sin E = N,
r
v
1+e E
tg =
tg ,
2
1−e
2
2
a(1 − e )
r=
,
1 + e cos v
ω = v + ω,
(13)
(14)
(15)
(16)
(17)
x = r(cos u cos Ω − sin u sin Ω cos i),
y = r(cos u sin Ω + sin u cos Ω cos i),
(18)
z = r sin u sin i.
q
Здесь M — средняя аномалия; n = aµ3 — среднесуточное движение, t0 — на-
чальный момент времени (эпоха), M0 — средняя аномалия в эпоху, E — эксцентрическая аномалия, v — истинная аномалия, u — аргумент широты.
Пусть V — скорость, Vr — радиальная скорость, Vk — трансверсальная скорость. Тогда
2 1
2
V =µ
−
,
(19)
r a
r
µ
Vr =
e sin v,
(20)
p
r
µ
Vn =
(1 + e cos v).
(21)
p
225
З а у с а е в А. Ф., З а у с а е в Д. А.
Для нахождения экваториальных координат скорости продифференцируем по времени формулы (18), получим
x
Vr + (− sin u cos Ω − cos u sin Ω cos i)Vn ,
r
y
ẏ = Vr + (− sin u sin Ω + cos u cos Ω cos i)Vn ,
r
z
ż = Vr + cos u sin iVn .
r
ẋ =
(22)
Для учёта возмущающего действия от больших планет уравнения движения возмущаемого тела следует представить в гелиоцентрической системе
координат. При переходе от барицентрической системы координат к гелиоцентрической уравнения движения возмущаемого тела, приведенные в работе
[1], с учётом формулы (3) запишутся в виде
2
d2 x
x
3aoo ros
p
p
· + X,
=
−
2
3
3
3 +
3 )2 r
dt
r 2 + r r 3 − ros
(r 3 − ros
2
y
3aoo ros
d2 y
p
p
·
=
−
+ Y,
2
3
3 + 3 (r 3 − r 3 )2 r
dt
r 2 + r r 3 − ros
os
где
X=
Y =
Z=
P
P
P
i
i
i
(23)
2
z
3aoo ros
d2 z
p
p
· + Z,
=
−
3
3
3
3
2
2
3
3
dt2
r + r r − ros + (r − ros ) r
xi −x
∆i
·
∆2i +∆i
yi −y
∆i
·
∆2i +∆i 3
yi −y
∆i
·
∆2i +∆i 3
√
3
3 +
∆3i −roi
√
3
3 )2
(∆3i −roi
2
3aoi roi
3
3
∆i −roi + 3
3 )2
(∆3i −roi
2
3aoi roi
3
3
∆i −roi + 3
3 )2
(∆3i −roi
√
√
2
3aoi roi
√
√
+
+
+
2
3aoi roi
xi
√
√
3
ri r 2 +ri r 3 −r 3 + 3 (r 3 −r 3 )2
i
i
oi
i
oi
2
yi
√ 3aoi roi√
ri r 2 +ri 3 r 3 −r 3 + 3 (r 3 −r 3 )2
i
i
oi
i
oi
2
zi
√ 3aoi roi√
ri r 2 +ri 3 r 3 −r 3 + 3 (r 3 −r 3 )2
i
i
oi
i
oi
,
,
,
∆2i = (xi − x)2 + (yi − y)2 + (zi − z)2 ; roi — эффективный радиус i-того тела;
aoi — соответствующее ускорение для i-того тела на расстоянии roi от центра
массы; X, Y , Z — возмущающие ускорения от больших планет.
Исправляя найденные координаты и скорости за счет возмущающего действия от больших планет, окончательно получим координаты и скорости на
следующем шаге интегрирования:
h2
h2
h2
′
′
X, yn+1 = yn+1
+ Y, zn+1 = zn+1
+ Z,
2
2
2
′
′
′
= ẋn+1 + hX, ẏn+1 = ẏn+1 + hY, żn+1 = żn+1 + hZ,
xn+1 = x′n+1 +
ẋn+1
(24)
(25)
′
′
′
′
— координаты и компоненты вектора
, żn+1
, ẋ′n+1 , ẏn+1
, zn+1
где x′n+1 , yn+1
скорости возмущаемого тела на n + 1 шаге, вычисленные без учёта возмущающего действия от больших планет; xn+1 , yn+1 , zn+1 , ẋn+1 , ẏn+1 , żn+1 —
координаты и компоненты вектора скорости возмущаемого тела на n+1 шаге,
полученные с учётом возмущающего действия от больших планет и Луны.
Таким образом, координаты радиус-вектора и скорости возмущаемого тела в точке tn+1 , если они известны в точке tn , методом оскулирующих элементов вычисляются следующим образом:
226
Интегрирование уравнений движения малых тел Солнечной системы . . .
a) по формулам (8)–(12) находятся элементы орбит исследуемого объекта
на момент времени tn ;
б) с использованием формул (13)–(22) вычисляются координаты и скорости данного объекта на момент времени tn+1 , при этом исправляется
величина аргументы перигелия по формулам (12);
в) с помощью формулы (24) определяются окончательные координаты и
компоненты скорости объекта на момент времени tn+1 с учётом возмущающего действия от больших планет и Луны.
Эффективность предложенного метода можно существенно повысить, если для координат и скоростей больших планет использовать заранее полученный банк данных.
Работа выполнена при финансовой поддержке Федерального агентства по образованию (РНП. 2.1.1/745).
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Заусаев А. Ф., Заусаев Д. А. Эволюция методы интерполяции, используемые для получения координат и элементов орбит больших планет и малых тел Солненой системы //
Вестн. Сам. гос. техн. ун-та. Сер. Физ.-мат. науки, 2008. — № 2(17). — C. 231–238.
2. Standish E. M. JPL. Planetary and Lunar Ephemerides, DE405/LE / In: Jet Lab Technical
Report. IOM 321. F-048, 1998. — P. 1–7.
3. Богородский А. Ф. Всемирное тяготение. — Киев: Наукова думка, 1971. — 353 с.
4. Заусаев А. Ф., Заусаев А. А. Математическое моделирование орбитальной эволюции малых тел Солнечной системы. — М.: Машиностроение-1, 2008. — 250 с.
Поступила в редакцию 08/XII/2008;
в окончательном варианте — 09/II/2009.
MSC: 85-08
INTEGRATION OF EQUATIONS OF SOLAR SYSTEM SMALL
BODIES MOTION WITH OSCULATING ELEMENTS METHOD
A. F. Zausaev, D. A. Zausaev
Samara State Technical University,
244, Molodogvardeyskaya str., Samara, 443100.
E-mails: Zausaev_AF@mail.ru
Algorithm of osculating elements method for numerical integration of equations of
Solar system small bodies’ motion is introduced.
Key words: interpolation, numerical integration, differential equations of motion,
method of osculating elements.
Original article submitted 08/XII/2008;
revision submitted 09/II/2009.
Zausaev Anatoliy Fedorovich, Ph. D. (Phys. & Math.), Prof., Dept. of Applied Mathematics and
Computer Science.
Zausaev Dmitriy Anatol’evich, Student, Dept. of Applied Mathematics and Computer Science.
227
Документ
Категория
Без категории
Просмотров
3
Размер файла
286 Кб
Теги
методов, уравнения, движение, солнечной, система, элементов, малыш, тел, оскулирующих, интегрированный
1/--страниц
Пожаловаться на содержимое документа