АСТРОНОМИЧЕСКИЙ ЖУРНАЛ, 2019, том 96, № 5, с. 418-430

УДК 521.1

ОСРЕДНЕННЫЕ УРАВНЕНИЯ ДВИЖЕНИЯ

ПРИ ВОЗМУЩАЮЩЕМ УСКОРЕНИИ, МЕНЯЮЩЕМСЯ

ПО ЗАКОНУ ОБРАТНЫХ КВАДРАТОВ

¹Санкт-Петербургский государственный университет, Санкт-Петербург, Россия

²Институт прикладной астрономии РАН, Санкт-Петербург, Россия

Поступила в редакцию 01.10.2018 г.; после доработки 16.12.2018 г.; принята к публикации 17.12.2018 г.

Рассмотрено движение точки нулевой массы под действием притяжения к центральному телу S и

возмущающего ускорения P^′, величина которого обратно пропорциональна квадрату расстояния до

S. Направление P^′ постоянно в одной из трех наиболее употребительных в астрономии системах

отсчета: основная инерциальная O и две орбитальные O_s с осью x по радиусу-вектору при s = 1 и

по вектору скорости при s = 2. Отношение |P^′| к основному ускорению, вызванному притяжением

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

осредняющее преобразование в первом приближении по малому параметру. Получены замкнутые

выражения для правых частей уравнений движения в средних элементах. В системах O, O₁ они

выражены через элементарные функции; в системе O₂ появляются полные эллиптические интегралы.

Получены замкнутые выражения функций замены переменных. В системах O, O₁ все встречающиеся

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

интегралом от элементарной функции, а также рядом по степеням эксцентриситета, сходящимся

абсолютно и равномерно при 0 ≤ e ≤ 1. В системе O₂ все функции, кроме тех, что определя-

ют вариацию средней аномалии, выражены через неполные эллиптические интегралы. Вариация

средней аномалии вычисляется с помощью ряда Фурье по средней аномалии. Интегрирование

осредненных уравнений движения будет выполнено в последующих работах. Возможные приложения

рассмотренной модельной задачи: движение астероида с учетом эффекта Ярковского-Радзиевского

и движение космического аппарата с солнечным парусом, когда возмущающее воздействие обратно

пропорционально квадрату расстояния от Солнца. Разумеется, определение компонентов вектора P^′

требует знания теплофизических характеристик тела и параметров его вращательного движения в

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

DOI: 10.1134/S0004629919050050

1. ВВЕДЕНИЕ

по вектору скорости, главной нормали к оскули-

рующей орбите и бинормали. В качестве вспомо-

Пусть в R³ имеется система двух точечных масс:

гательной понадобится также система O₃ с ор-

неподвижное тело S (например, Солнце) массой

тами, направленными в перицентр оскулирующей

m₀ и малое тело A (напр., астероид), движущееся

орбиты, по нормали к i₃ в плоскости оскулирую-

под действием силы притяжения к точке S и возму-

щей орбиты в сторону движения и бинормали [1].

В [2, 3] рассмотрены задачи о движении A при

щающего ускорения P^′. Введем три координатные

системы с общим началом в S, но с разными

возмущающем ускорении P^′, постоянном в одной

направлениями осей: основная инерциальная O с

из описанных систем отсчета. Там же указаны аст-

ортами i, j, k и две сопутствующие O_s с орта-

рономические приложения этих модельных задач.

ми i_s, j_s, k_s. Орты системы O₁ направлены по

Теперь мы рассмотрим аналогичные задачи при

радиусу-вектору, трансверсали — перпендикуляр-

P′ = P/r², r = SA, где постоянным в одной из

но к радиусу-вектору в плоскости оскулирующей

систем отсчета является вектор P. Отношение

орбиты в сторону движения, и бинормали — по

модулей возмущающего и основного ускорения

вектору площадей. Орты системы O₂ направлены

(

)(

)-1

P /r²

κ²/r²

= P/κ² = μ постоянно и счита-

^*E-mail: TNSannikova@gmail.com

ется малым. Здесь κ² — произведение постоян-

^**E-mail: kvk@astro.spbu.ru

ной тяготения на массу S. В настоящей работе

418

ОСРЕДНЕННЫЕ УРАВНЕНИЯ ДВИЖЕНИЯ

419

мы получим уравнения движения типа Лагранжа

В основной системе O имеем:

в оскулирующих элементах, явный вид осредня-

]

ющей замены переменных и уравнения движения

^[a²

a²

f₁ =

Φ₁ sin θ -

(e + cos θ)Φ₂

(1)

в средних элементах. В следующей работе будет

a³η r²

r²

рассмотрено интегрирование последних уравнений.

[

Задачи при P^′ = const и P = const во мно-

f₂ =

(e + cos θ)Φ₁ sin θ +

κa3/2η

гом схожи, что позволяет часто ограничиваться

]

ссылками на [3]. Основное отличие заключается в

(3 + 4e cos θ + cos 2θ)Φ₂ ,

использовании истинной аномалии наряду с экс-

центрической, поскольку в этом случае множитель

r-2 становится тригонометрическим многочленом.

f₃ =

Φ₃ cos w,

κa3/2rη

Можно указать по меньшей мере два приложе-

ния рассматриваемой задачи. Это движение асте-

f₄ =

Φ₃ sin w,

роида с учетом эффекта Ярковского-Радзиевского

κa3/2rη sin i

и движение космического аппарата с солнечным

^[a

парусом. В обоих случаях возмущающее воздей-

f₅ = -

(3 + 2e cos θ - cos 2θ)Φ₁ -

2κa3/2eη r

ствие обратно пропорционально квадрату рассто-

]

яния от Солнца.

Φ₂ sin 2θ

- f4 cosi,

2. УРАВНЕНИЯ ДВИЖЕНИЯ

(3 - 2e cos θ - cos 2θ)Φ₁ -

κa3/2e

В качестве переменных выберем оскулирующие

]

элементы ω, e, i, Ω, σ, M — среднее движение,

(cos θ + 2e)Φ₂ sin θ .

эксцентриситет, наклон, долгота восходящего узла,

аргумент перицентра и средняя аномалия соот-

ветственно. Первые пять из них образуют вектор

Здесь и ниже a = κ2/3ω-2/3 — большая полуось,

√

медленных переменных x = (x₁, . . . , x₅), а послед-

p = aη² — фокальный параметр, η =

1-e², r=

ний — быструю переменную y. Во всех использу-

= p/(1 + ecos θ) = a(1 - ecos E), где θ — истин-

емых системах отсчета уравнения движения типа

ная и E — эксцентрическая аномалии, w = σ +

Эйлера [4] имеют форму

+ θ — аргумент широты. Компоненты P во вспо-

x = μf(x,y),

могательной системе O₃ обозначены через Φ_s. Они

y = x₁ + μg(x,y)

связаны с компонентами P_s в основной системе

соотношением

с малыми (порядка μ) правыми частями, различа-

ясь лишь видом функций f = (f₁, . . . , f₅), g. По-

Φ_s = b1sP₁ + b2sP₂ + b3sP₃,

следние фактически приведены в [1]: надо лишь

заменить f_s, g на f_s/r², g/r².

где b_ks — элементы матрицы вращения

⎛

⎞

cos σ cos Ω - cos i sin σ sin Ω

- sin σ cos Ω - cos i cos σ sin Ω

sin i sin Ω

⎜

⎟

⎜

⎟

⎜

os σ sin Ω + cos i sin σ cos Ω

- sin σ sin Ω + cos i cos σ cos Ω

- sin i cos Ω⎟

⎝^c

⎠

sin i sin σ

sin i cos σ

cos i

Таким образом, функции Φ_s зависят только от мед-

В сопровождающей системе O₁

ленных переменных i, σ, Ω. Тригонометрические

функции от аргумента широты равны

(

)

a³

cos w = cos(σ + θ) = cos σ cos θ - sin σ sin θ,

f₁ = -

S sin θ + η²

(2)

a³η

r²

r³

sin w = sin(σ + θ) = sin σ cos θ + cos σ sin θ.

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

420

САННИКОВА, ХОЛШЕВНИКОВ

[

(

) ]

a²

a³

Уравнения движения в новых (средних) элемен-

f₂ =

S sin θ + η²

κa3/2e

r²

r³

тах:

X=F(X),

(6)

f₃ =

W cos w,

κa3/2rη

= X₁ + G(X),

f₄ =

W sin w,

κa3/2rη sin i

где

η a

F(X) = Ef(X, Y ), G(X) = Eg(X, Y ).

f₅ = -

S cos θ +

κa3/2e r²

)

Здесь и ниже мы пользуемся интегральными опе-

sin θ

⁽a

раторами

+η²

T - f4 cosi,

κa3/2ηe r

r²

∫

(-3e + 2 cos θ + e cos 2θ)

S-

Ef(X,Y ) =

f (X, Y ) dY,

(7)

2κa3/2e

2π

)

−π

sin θ

⁽a

∫

+η²

κa3/2e r

r²

If(X,Y ) =

[f(X, Y ) - Ef(X, Y )] dy.

Постоянные S, T , W — компоненты P в системе

O₁. Заметим, что W = Φ₃.

Последняя первообразная определяется однознач-

но условием нулевого среднего

В сопровождающей системе O₂

EIf(X,Y ) = 0.

(8)

f₁ = -

(3)

aη

r²

Во всех используемых системах отсчета функции

2√p(e + cos θ) T

r√p sin θ N

f,g выражены явно через истинную аномалию.

f₂ =

κϑ

r²

κaϑ r²

При необходимости их легко выразить через экс-

2√p sin θ T

центрическую аномалию. Поэтому интегралы (7)

f₅ =

вычисляются переходом к θ или E:

κeϑ r²

(

)

2e +

1+e²

cos θ N

r²

- f4 cosi,

dy = dM =

dθ =

dE,

(9)

√peϑ

a²η

2sin θ(1 + e² + ecos θ) T

при этом пределы интегрирования -π, π не ме-

g=-

няются. Напомним, что среднее значение нечетной

κe√aϑ

функции равно нулю, а среднее значение четной

(1 - e²) cos θ N

функции можно найти интегралом от 0 до π.

κe√aϑ r,

Определяющие замену переменных величины

где

√

u, v даются соотношениями

ϑ(θ, e) =

1 + e² + 2ecosθ.

(4)

u(x, y) =

If(x,y),

(10)

Постоянные T, N, W — компоненты P в системе

X₁

O₂. Мы опустили функции f₃, f₄, поскольку они

совпадают с приведенными в (2).

v(x, y) =

I [u₁(x,y) + g(x,y)] .

Обратим внимание, что u₁, . . . , u₅ находятся одной

3. ПРОЦЕДУРА ОСРЕДНЕНИЯ

квадратурой, тогда как v — двумя. Вторая — это

Выполним, ограничиваясь первым порядком

интеграл от u₁.

малости, близкую к тождественной осредняющую

Замечание 1. Функции u,v имеют первый по-

замену переменных (x, y) ←→ (X, Y ), устраняю-

щую быструю угловую переменную Y из уравнений

рядок малости по μ. Разности u(x, y) - u(X, Y )

в средних элементах (X, Y ). Подробности см. в [2,

периодичны как по y, так и по Y и имеют второй

3, 5]. Заметим лишь, что ни резонансов, ни малых

порядок. Поэтому в (5) безразлично, считать ли ар-

знаменателей мы не встретим, поскольку налицо

гументами u, v средние или оскулирующие элемен-

лишь одна быстрая переменная Y .

ты. Напротив, в (6) мы имеем дело именно со сред-

ними элементами. Однако для средних элементов,

Замена переменных:

отвечающих оскулирующим ω, e, . . . мы не вводим

x = X + u(X,Y ), X = x - u(x,y),

(5)

новых обозначений типа ω^′, e^′, . . ., поскольку в (6)

y = Y + v(X,Y ), Y = y - v(x,y).

оскулирующие элементы не встречаются.

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

ОСРЕДНЕННЫЕ УРАВНЕНИЯ ДВИЖЕНИЯ

421

2ω

4. ОСРЕДНЕННЫЕ УРАВНЕНИЯ

F₅ =

K(e)N - F₄ cos i,

ДВИЖЕНИЯ

πκ²

2ωη

Величины Φ₁, Φ₂ зависят лишь от медленных

K(e)N.

переменных, а Φ₃, S, T , W , T, N постоянны. По-

πκ²

этому вычисление F, G в системах O, O₁ сводится

Здесь и ниже использованы стандартные обозна-

к нахождению нескольких коэффициентов Ганзена

чения для полных и неполных эллиптических инте-

с нулевым нижним индексом (см. Приложение 1).

гралов в форме Лежандра:

Приведем правые части уравнений (6) в систе-

ме O:

∫

3ω²e

K(k) =

F₁ = -

Φ₂,

h(x, k)

κ²η²

ω(1 + 2η)

∫

F₂ =

Φ₂,

κ²(1 + η)

E(k) =

h(x, k) dx,

ωecos σ

F₃ = -

Φ₃,

κ²η(1 + η)

∫

ωesin σ

sin² x dx

K(k) - E(k)

F₄ = -

Φ₃,

D(k) =

κ²η(1 + η)sin i

h(x, k)

k²

ω(2 + η)

F₅ = -

Φ₁ - F₄ cos i,

∫

κ²e(1 + η)

F₁(ϕ,k) =

h(x, k)

ω(1 + 2η + e²)

Φ₁.

κ²e(1 + η)

∫

В знаменателях G и первого слагаемого F₅ присут-

F₂(ϕ,k) = h(x,k) dx,

ствует e, а в знаменателях F₄ и тем самым второго

слагаемого F₅ присутствует sin i.

∫

Правые части (6) в системе O₁ еще проще:

sin² x dx

F₁(ϕ,k) - F₂(ϕ,k)

F₃(ϕ,k) =

3ω

h(x, k)

k²

F₁ = -

(11)

κ²η²

ωe

F₂ =

√

κ²(1 + η)

1 - k2 sin2 x. В формулах (12)

где h(x, k) =

ωecos σ

F₃ = -

2√e

κ²η(1 + η)

ωesin σ

1+e

F₄ = -

κ²η(1 + η)sin i

С ростом e от 0 до 1 параметр k также возрастает от

0 до 1, причем e ≤ k; равенство достигается только

F₅ = -F₄ cos i

при e = k = 0 и e = k = 1.

2ω

G=-

В уравнениях (12) мы опустили F₃, F₄, т.к. они

κ²

совпадают с приведенными в (11). Как и в системе

По-прежнему в знаменателях F₄ и F₅ присутствует

O₁, в знаменателях F₄ и F₅ присутствует sin i, а

sin i, однако эксцентриситет в знаменателях не

эксцентриситет в знаменателях не появляется.

появляется. Интересно, что радиальный компонент

Обратим также внимание на присутствие η или

S возмущающего ускорения не влияет на орбиту A,

1 - e в знаменателях F₁, F₃, F₄, F₅ во всех си-

а лишь на положение A на орбите.

стемах отсчета. Функции F₂, G свободны от этого

Перейдем к системе O₂, где коэффициенты Ган-

недостатка. Однако в системе O₂ они пропор-

зена не работают, и появляются эллиптические ин-

циональны комбинациям эллиптических интегра-

тегралы. Они собраны в Приложении 3. Приведем

лов, стремящихся к бесконечности при e → 1 (см.

окончательный результат:

формулы (12)). Эта сингулярность несущественна,

6ω

поскольку

F₁ = -

E(k)T,

(12)

[

]

πκ²(1 - e)

2D(k)

[

]

lim

K(k) -

= 1, lim ηK(e) = 0.

4ω

2D(k)

e→1

1+e

e→1

F₂ =

K(k) -

πκ²

1+e

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

422

САННИКОВА, ХОЛШЕВНИКОВ

[

Замечание 2. Эллиптические интегралы вы-

1+e²

u₅ = -

2η

(θ - M) -

числяются тем быстрее, чем меньше их модуль.

2κ²ηe

e²

Поскольку e < k при 0 < e, k < 1, целесообразно

]

2η²

2η

перейти от модуля k к модулю e по приведенным

−

(E - M) -

sin θ Φ₁ -

e²

в [6, 7] формулам¹ :

[

K(k) = (1 + e)K(e),

(13)

−

e cos θ - ln(1 + e cos θ) +

κ²e³

2E(e) - (1 - e²)K(e)

]

E(k) =

2η²

1+e

+ ln

+ η - η² Φ₂ - u4 cosi.

1+η

Отсюда следует

E(e) - (1 - e²)K(e)

Для определения v необходимо, кроме интеграла от

K(k) -

D(k) =

(14)

1+e

[

Формулы (13), (14) можно использовать при вы-

Ig=

η³(θ - M) +

(16)

числении функций F₁, F₂, определенных соотноше-

κ²e³

]

ниями (12).

+ (3e² - 1)(E - M) - ηe sin θ) Φ₁ +

[

(

)

5. ЗАМЕНА ПЕРЕМЕННЫХ

cos θ +

κ²e e

e²

Вычислим определяющие замену переменных

(

)

функции u_s, v по формулам (10).

2η²

× ln(1 + e cos θ) - ln

1+η

5.1. Основная система отсчета O

(

η)]

+ (1 - η)

Φ₂,

Все необходимые интегралы от функций (1) вы-

e²

числены в Приложениях 1 и 2. Приведем оконча-

тельный результат для u_s.

вычислить интеграл от u₁. Интегралы от cos θ + e

{

и sin θ вычислены в Приложении 2 (п. 2 и п. 3).

3ω

u₁ = -

(cos θ + e)Φ₁ +

(15)

Интеграл от разности аномалий разбивается на два

κ²η²

}

интеграла (см. Приложение 2, п. 4):

+ [sin θ + e(θ - M)]Φ₂

I(θ - M) = I(θ - E) + I(E - M) =

^[1

(

e⁾

e²

u₂ =

cos θ -

= I(θ - E) - e cosE +

cos 2E.

κ²

(

)

2η²

Поэтому

ln(1 + e cos θ) - ln

e²

1+η

{

(

)]

[

Iu₁ =

- η²Φ1 sin E +

+ (1 - η)

Φ₁ +

η(E - M) +

κ²η²

e²

κ²e²

]

^[ e(η + e²)

(

)

+ (η + e²) cos E -

2e² - 1

(θ - M) + e sin θ Φ₂,

{

]

}

e³

cos 2E - eI(θ - E) Φ₂

u₃ =

cos σ [η(θ - M) - (E - M)] +

κ²ηe

[

]}

что вместе с (16) дает

2η²

+ ηsinσ ln(1 + ecosθ) + 1 - η - ln

Φ₃,

[

]

1+η

η³(θ - M) - (E - M) - ηesin θ

Φ₁ +

{

κ²e³

[

(

)

u₄ =

sin σ [η(θ - M) - (E - M)] -

κ²ηesin i

cos θ +

κ²

e²

[

]}

(

)

2η²

- ηcosσ ln(1 + ecosθ) + 1 - η - ln

Φ₃,

× ln(1 + e cos θ) - ln

1+η

(

η-η²

η)

3e(η + e²)

η+e²

¹ Вторая из формул (13) приведена в [8, раздел 9C] с

опечаткой.

e²

2η²

η²

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

ОСРЕДНЕННЫЕ УРАВНЕНИЯ ДВИЖЕНИЯ

423

]

5.3. Сопровождающая система отсчета O₂

× cos E -

cos 2E -

I(θ - E) Φ₂.

4η²

η²

Все необходимые интегралы от функций (3) вы-

числены в Приложении 4.

Функцию I(θ - E) можно вычислять по форму-

Выпишем u_s:

лам (27) или (28) Приложения 2.

[

(

)

]

6ω

u₁ = -

F₂

E(k)M T,

κ²(1 - e)

{

5.2. Сопровождающая система отсчета O₁

(

)

u₂ =

F₁

K (k) M -

Действуя, как и в предыдущем разделе, получим

κ²

}

3ω

[

(

)

]

u₁ =

[e(cos θ + e)S - (e sin θ + θ - M)T ] ,

κ²η²

−

F₃

D (k) M T +

(1 + e)

u₂ = -

(cos θ + e)S +

[

]

κ²

2η

arctg

[η²K(e) - E(e)] N,

κ²e

[(θ - M) - η(E - M) + e sin θ] T,

[

]

κ²e

2η

[

u₅ = -

ϑ-

E(e) T +

κ²e²

u₅ = -

S sin θ -

e cos θ + e² +

[

κ²e

κ²e²

(

))

(

)

]

2η

F₁

E+

- K(e)

+ ln(1 + e cos θ) + 1 - η - ln

T - u4 cosi.

κ²

1+η

√

]

e sin E +

1 - e2 cos2 E

Выражения для u₃, u₄ мы опустили, поскольку они

N - u4 cosi,

совпадают с приведенными в (15) с учетом Φ₃ = W .

e²

[

Как и выше, находим:

η²

Ig=

arctg

K(e) +

Ig=

[-2e(E - M) + η sin θ] S +

(17)

κ²

κ²e

[

(

)]

e cos θ + e² + ln(1 + e cos θ) +

ϑ-

E(e) T +

κ²e²

e²

]

2η

[

(

)

+ 1 - η - ln

(

)

1+η

F₁

E+

- K(e)

M -

κ²

]

√

Iu₁ =

S sin E +

(18)

e sin E +

1 - e2 cos2 E

κ²

[

(

−

e⁾

e²

e(1 + η) cos E +

κ²η²

Выражения для u₃, u₄ мы опустили, поскольку они

]

совпадают с приведенными в (15) с учетом Φ₃ = W .

−

cos 2E - I(θ - E) T.

Для определения v необходимо, кроме приве-

денного выше интеграла от g, вычислить интеграл

В итоге, суммируя (17) и (18), получим

[

]

от u₁. Обозначим

(

)

(E - M) +

sin θ S +

E(k)

κ²

F₂

M =

(θ - M) + H(θ, k),

{

(

e⁾

3e²

(

)

3e(1 + η) cos E +

cos 2E +

E(k)

κ²η²

H(θ, k) = F₂

θ.

[

]

η³

2η²

cos θ +

ln(1 + e cos θ) - ln

Обе переменные θ - M и H нечетны и

2π-

e²

1+η

периодичны как функции от θ. Интеграл I(θ - M)

}

приведен в Приложении 2 (п. 4). Воспользуемся

(2 + η)η

- 3I(θ - E) T.

рядом Фурье для H:

1+η

∑

Функцию I(θ - E) можно вычислять по форму-

H(θ, k) =

(-1)m-1B_m(k)k2m sin mθ.

(19)

лам (27) или (28) Приложения 2.

m=1

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

424

САННИКОВА, ХОЛШЕВНИКОВ

Коэффициенты B_m(k) можно представить ря-

Интегрирование элементарно:

дом [6]:

∑

C_n

B_m(k) =

(20)

IH =

cos nM.

(24)

n=1

∑

(s + 1) · · · (s + m)

Учитывая, что e = O(k²), заключаем, что C_n =

(s + m + 1) · · · (s + 2m)

s=0

= O(k2n), и ряды (23), (24) быстро сходятся при

малых и умеренных эксцентриситетах.

^[ (2s + 2m - 1)!!]2

(2s + 2m)!!

2s + 2m - 1

Окончательно,

{

[

(

E(k)

e⁾

По уточненной формуле Валлиса [9]

Iu₁ =

e cos E +

κ²(1 - e) π

^[ (2s + 2m - 1)!!]2

]

}

E(k)

(2s + 2m)!!

π(s + m)

−

cos 2E

I(θ - E) - IH T,

{

Далее, при m ≥ 1

[

(s + 1) · · · (s + m)

2s + 2m + 1

2(1 - e) arctg

< 1.

κ²(1 - e)

(s + m + 1) · · · (s + 2m)

2s + 2m

(

)]

Эти оценки позволяют вывести из (20) неравенство

−

K(e) +

ϑ-

E(e)

e²

[

]

(

B_m(k) <

(21)

3E(k)

e⁾

e²

πm

e cos E +

cos 2E

∑

k2s

}

3E(k)

(2s + 2m + 1)(2s + 2m - 1)

s=0

I(θ - E) - 3IH T +

∑

[

(

)

(

)

πm

(2s + 2m + 1)(2s + 2m - 1)

s=0

F₁

E+

- K(e)

M -

κ²

Сумма в последнем выражении равна 1/2(2m - 1).

]

√

В итоге получаем оценку B_m(k):

e sin E +

1 - e2 cos2 E

−

e²

B_m(k) <

m ≥ 1.

πm(2m - 1)

Функцию I(θ - E) можно вычислять по форму-

Если оценить отдельно слагаемое в (20), отвеча-

лам (27) или (28) Приложения 2, а функцию IH —

ющее s = 0, и суммировать в (21) от единицы до

по формуле (24).

бесконечности, то последняя оценка при m > 1

несколько улучшается:

6. ЗАКЛЮЧЕНИЕ

B_m(k) <

m ≥ 2.

В задаче о движении тела в трех наиболее упо-

2πm²

требительных системах отсчета с возмущающим

Перейдем к средней аномалии

ускорением, обратно пропорциональным квадрату

расстояния r = SA, получены замкнутые формулы

∑

для уравнений движения в средних элементах в

sin mθ =

S0mn(e)sin nM,

(22)

первом порядке по малому параметру. В систе-

n=1

мах O и O₁ правые части уравнений движения —

где S0mn(e) = X0mn(e) - X0m-n(e) можно найти в [4,

элементарные функции элементов. В системе O₂

правые части уравнений движения выражаются

10-12]. Как известно, S0mn(e) имеет порядок

через полные эллиптические интегралы первого и

e^|n-m|. Подставляя (22) в (19), получим

второго рода в форме Лежандра.

∑

В системах O и O₁ функции замены переменных

H(θ, k) = - C_n sin nM,

(23)

элементарны для всех элементов, кроме средней

n=1

аномалии. Выражение для средней аномалии со-

∑

держит интеграл от сложной (хотя элементарной)

C_n =

(-1)^mB_m(k)S0mn(e)k2m.

функции, или бесконечный ряд, сходящийся абсо-

m=1

лютно и равномерно при всех 0 ≤ e ≤ 1.

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

ОСРЕДНЕННЫЕ УРАВНЕНИЯ ДВИЖЕНИЯ

425

В системе O₂ функции замены переменных эле-

Приложение 2

ментарны для i, Ω, а для a, e, σ содержат полные и

неполные эллиптические интегралы первого и вто-

ДЕЙСТВИЕ ОПЕРАТОРА I НА ФУНКЦИИ,

рого рода в форме Лежандра. Для средней анома-

ВХОДЯЩИЕ В СООТНОШЕНИЯ (1) И (2)

лии появляется бесконечный ряд, коэффициенты

которого сами являются рядами.

Согласно (7), (9)

∫

В следующей работе будет рассмотрено инте-

If = (f - Ef) dM =

грирование осредненных уравнений.

∫

r²

(f - Ef) dθ =

(f - Ef) dE.

БЛАГОДАРНОСТИ

a²η

Авторы благодарны анонимному рецензенту за

Постоянная интегрирования определяется из усло-

ценные замечания, учтенные в окончательной ре-

вия нулевого среднего (8). Полезной часто оказы-

дакции рукописи. Работа выполнена при финансо-

вается формула

вой поддержке Российского научного фонда (грант

∫

18-12-00050).

If = f dM - (Ef)M.

Приложение 1

Ниже под знаком интеграла оказываются ли-

бо четные, либо нечетные периодические функции

быстрой переменной. В первом случае необходи-

ДЕЙСТВИЕ ОПЕРАТОРА E НА ФУНКЦИИ,

ВХОДЯЩИЕ В СООТНОШЕНИЯ (1) И (2)

мо сначала вычислить Ef, что сделано в При-

ложении 1. Зато в результате получаем нечетную

Рассматриваемые функции имеют вид

функцию с нулевым средним. Во втором случае

(r/a)ⁿ cos kθ, их средние значения равны коэф-

автоматически Ef = 0. Зато после интегрирования

фициентам Ганзена с нулевым нижним индексом

получаем четную функцию, и нужно еще опреде-

Xnk0(e). Свойства коэффициентов Ганзена и таб-

лить ее среднее значение.

лица их первых значений приведены в [4, 13, 14].

1. Вспомогательные формулы

Выпишем нужные нам коэффициенты Ганзена

d(E - M)

X-3,00 = η-3, X-2,00 = η-1, X-2,10 = 0,

- 1,

(25)

)

X-1,00 = 1, X-1,10 = -e

X-1,20 =1-η

d(θ - M)

⁽a²

1+η

=η

r²

и их комбинации

d sin θ

ae cos 2θ

a cos θ

eX-2,00 + X-2,10 =^e

являются простыми следствиями (9). Из (25) нахо-

дим нужные нам линейные комбинации:

3X-1,00 + 4eX-1,10 + X-1,20 =2η(1+2η),

1+η

ae cos θ

e²

(26)

3X-1,00 + 2eX-1,10 - X-1,20 =2η(2+η),

1+η

d(θ - M)

d(E - M)

=η

3X-1,00 - 2eX-1,10 - X-1,20 =2(1+2η+e2),

1+η

ae cos 2θ

d sin θ

2a cos θ

= 2η

(1 - e²)X-3,00 - X-1,00 =^e

η(1 + η)

Используя (26), получим:

−3eX-1,00 + 2X-1,10 + eX-1,20 = -4e.

ae² cos 2θ

1-η

d sin θ

-e²

= 2eη

1+η

Понадобятся также соотношения

d(θ - M)

d(E - M)

E (cos θ) = X0,10 = -e,

− 2η

+ (2 - e²)

∫

E (cos E) =

cos E(1 - e cos E) dE = -

2. Теперь легко определяются интегралы от

четных функций :

)

∫ (

)

⁽a

E (cos nE) = 0, n > 1.

-1

dM = E - M,

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

426

САННИКОВА, ХОЛШЕВНИКОВ

)

∫ (

)

⁽a

a²

θ-M

Чтобы определить постоянные интегрирования, мы

dM =

r²

воспользовались формулой

)

∫ (

)

∫

⁽a

a³

dM =

ln(1 + e cos θ) dM =

r³

η³

∫

(1 + e cos θ) dθ -

(1 - e cos E)[ln η² - ln(1 - e cos E)] dE =

η³

∫

θ - M + esinθ

(1 - e cos E) ×

η³

∫ (

)

[

]

(

)

2η

cos θ

cos θ +

dM =

× ln

- ln(1 - 2β cos E + β²) dE,

1+η

[η(θ - M) - (E - M)] ,

где

∫ (

)

(

)

2β

1-η

β=

cos 2θ

cos 2θ -

dM =

1+η

1+β²

1+η

[

]

Поскольку при β² < 1 [7, п. 4.224.14, п. 4.397.6]

2η(e sin θ - θ + M) + (2 - e²)(E - M)

e²

справедливы выражения

(

)

∫

a²

∫

cos θ

cos θ dM =

r²

ln(1 - 2β cos E + β²) dE = 0,

∫

2π

cos θ dθ =

sin θ.

−π

∫

Понадобится также соотношение

cos E ln(1 - 2β cos E + β²) dE = -β,

2π

∫

-π

I(cos θ + e) = (cos θ + e) dM =

приходим к равенству

∫

=η²

cos E dE = η² sin E.

2η²

E ln(1 + e cos θ) = -1 + η + ln

1+η

3. Интегралы от нечетных функций берутся

непосредственно.

4. Интегралы от разности аномалий

∫

(

)

sin θdθ

sin θ

=η

I(E - M) = e sin E(1 - e cos E) dE =

1 + ecosθ

[

]

(

2η²

e⁾

e²

ln(1 + e cos θ) + 1 - η - ln

= -e cos E +

cos 2E.

1+η

∫

(

)

sin 2θdθ

С истинной аномалией положение сложнее. Соот-

sin 2θ

=η

1 + ecosθ

ношение [13, 15]

[

2η

ln(1 + e cos θ) - e cos θ -

β sin E

θ - E = 2arctg

e²

1 - β cosE

]

2η

- η(1 - η) - ln

влечет

1+η

∫

(

)

∫

cos θ + e

I(θ - E) = 2 ψ(E) dE,

sin θ

sin θ dθ = -

r²

β sin E

ψ(E) = (1 - e cos E) arctg

dE.

Понадобится также соотношение

1 - β cosE

∫

(

e⁾

Взять интеграл аналитически не удалось. Можно

I(sin θ) = η sin E dE = -η cos E +

предложить два варианта действий.

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

ОСРЕДНЕННЫЕ УРАВНЕНИЯ ДВИЖЕНИЯ

427

4a. Численное интегрирование

Приложение 3

I(θ - E) = I(E) - I₀,

ДЕЙСТВИЕ ОПЕРАТОРА E НА ФУНКЦИИ,

где

ВХОДЯЩИЕ В СООТНОШЕНИЯ (3)

∫E

Ключевую роль играет величина (4). Ее целесо-

I(E) = 2

ψ(t) dt, I₀ = EI(E).

образно представить в разных формах:

-π

√

ϑ=

1 + e² + 2ecosθ =

Вычисление I₀ можно свести к однократному инте-

√

1 + ecosE

грированию [16, п. 597, задача 11]:

= (1 + e)

1 - k2 sin2 ϕ = η

1 - ecosE

∫

где

I₀ =

dE ψ(t) dt =

(π - t)ψ(t) dt.

2√e

−π

-π

ϕ=

1+e

Функция ψ(t) нечетна, а tψ(t) четна, поэтому

С помощью (9) найдем средние значения трех

∫

четных функций:

I₀ = -

tψ(t) dt.

∫

a²ϑ

2(1 + e)

ϑ dθ =

r²

πη

Окончательно,

∫

√

∫E

∫

2(1 + e)

1 - k2 sin2 ϕ dϕ =

E(k),

I(θ - E) = 2 ψ(t) dt +

tψ(t) dt.

(27)

πη

−π

∫

a²

dθ

E r2ϑ ⁼

πη

πη(1 + e)

4b. Ряд по степеням β [4, 13] есть

∑

∫

βⁿ

θ-E=2

sin nE.

dϕ

√

K(k),

n=1

πη(1 + e)

1 - k2 sin2 ϕ

Отсюда

∫

a² cos θ

(1 - 2 sin² ϕ) dϕ

(1 - e cos E)(θ - E) =

√

r²ϑ

πη(1 + e)

1+β²

1 - k2 sin2 ϕ

∑

βⁿ

(1 - 2β cos E + β²) sin nE =

[K(k) - 2D(k)] .

πη(1 + e)

n=1

β(2 + β²)

Среднее значение трех других функций находим

sin E -

переходом к эксцентрической аномалии:

1+β²

∫

∑

n + 1 - (n - 1)β²

1 - ecosE

βⁿ sin nE.

√

dE =

(29)

n(n² - 1)

rϑ

πη

1 - e2 cos2 E

n=2

Интегрируя, получим

∫

)

√

K(e),

β(2 + β²)⁽e

πη

1 - e2 cos2 E

πη

I(θ - E) = -

+ cos E

(28)

1+β²

∫

∑

n + 1 - (n - 1)β²

a cos θ

cos E - e

βⁿ cos nE.

√

dE =

1+β²

n²(n² - 1)

rϑ

πη

1 - e2 cos2 E

n=2

Ряд сходится абсолютно и равномерно при 0 ≤ e ≤

≤ 1.

=-πηK(e),

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

428

САННИКОВА, ХОЛШЕВНИКОВ

∫

√

Интегралы от двух оставшихся функций из (3)

2η

Eϑ=

1 - e2 cos2 E dE =

E(e).

находим переходом к эксцентрической аномалии:

∫

1 - ecosx

Здесь использованы очевидные соотношения

√

dx -

Irϑ⁼

1 - e2 cos2 x

∫π

h₁(cos x) dx = 0,

K(e)M + const.

πη

∫

Разобьем интеграл на два. Первый (содержащий

единицу в числителе) подстановкой x = x^′ - π/2

h₂(cos x) dx = 2 h₂(cos x) dx,

сводится к стандартному эллиптическому, а второй

(содержащий e cos x в числителе) элементарен:

если h₁(-y) = -h₁(y), h₂(-y) = h₂(y).

∫

(

)

1 - ecosx

√

dx = F₁

E+

1 - e2 cos2 x

Приложение 4

(

√

)

- ln e sin E +

1 - e2 cos2 E

+ const.

ДЕЙСТВИЕ ОПЕРАТОРА I НА ФУНКЦИИ,

Постоянное слагаемое определим условием нечет-

ВХОДЯЩИЕ В СООТНОШЕНИЯ (3)

ности I(a/rϑ):

[

]

(

)

1. Интегралы от четных функций

F₁

E+

- K(e) -

K(e)M

rϑ

С учетом найденных в Приложении 3 средних

√

e sin E +

1 - e2 cos2 E

значений, получим

−

∫

√

a²ϑ

2(1 + e)

Мы использовали нечетность эллиптических инте-

1 - k2 sin2 ϕ dϕ -

r²

гралов F_s(ϕ, k) по первой переменной и свойство

F₁(ϕ + π,k) = F₁(ϕ,k) + 2K(k).

2(1 + e)

E(k)M =

πη

Аналогично получим

[

(

)

]

2(1 + e)

∫

F₂

E(k)M

a cos θ

cos x - e

√

dx +

rϑ

1 - e2 cos2 x

∫

dϕ

√

K(e)M + const,

r²ϑ

η(1 + e)

1 - k2 sin2 ϕ

πη

a cos θ

K(k)M =

rϑ

πη(1 + e)

[

]

(

)

[

(

)

]

F₁

E+

- K(e) -

K(e)M

F₁

K(k)M

η(1 + e)

√

e sin E +

1 - e2 cos2 E

∫

a² cos θ

(1 - 2 sin² ϕ) dϕ

eη

√

r²ϑ

η(1 + e)

1 - k2 sin2 ϕ

2. Интегралы от нечетных функций

[K(k) - 2D(k)] M =

Используя переход к истинной аномалии, выра-

πη(1 + e)

зим интегралы через элементарные функции и E(e)

[

(

)

(

)

∫

a² sin θ

sin θ dθ

F₁

- 2F₃

η(1 + e)

r²ϑ

]

(

)

2η

−

K(k)M +

D(k)M .

ϑ-

E(e)

ηe

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

ОСРЕДНЕННЫЕ УРАВНЕНИЯ ДВИЖЕНИЯ

429

∫

a sin θ

sin θ dθ

найдем

=η

rϑ

(1 + e cos θ)ϑ

2η²

π arcsin e

K(e) -

E(e) +

arctg

+ const,

∫

так что

a sin θ cos θ

sin θ cos θ dθ

=η

√

rϑ

(1 + e cos θ)ϑ

1-e

(

)

Ef = arctg

arcsin e +

2η

1+e

arctg

ϑ-

E(e)

+ const.

e²

[η²K(e) - E(e)].

Осталось найти среднее значение арктангенса.

Обозначим временно

Упростим последнее выражение. Положим

√

1 + ecosE

√

f (E) = arctg

= arctg

1-e

1 - ecosE

α = arcsine, β = arctg

1+e

e sin E

f^′(E) = -

√

Простая тригонометрия показывает, что

1 - e2 cos2 E

и проинтегрируем интеграл для среднего значения

⁽π

α⁾

по частям

tg β = tg

β =

∫π

то есть

πEf = f(E) d(E - esin E) =

Ef =

[η²K(e) - E(e)].

_π



= (E - e sin E)f(E)

Окончательно,

∫π

a sin θ

- (E - e sin E)f^′(E) dE =

(30)

rϑ

[

]

√

1-e

arctg

[η²K(e) - E(e)] ,

= π arctg

1+e

a sin θ cos θ

∫π

sin E(E - e sin E)

rϑ

√

dE.

[

]

1 - e2 cos2 E

arctg

[η²K(e) - E(e)] -

e²

Разобьем отрезок интегрирования для A на два

(

)

2η

точкой E = π/2 и во втором интеграле сделаем

−

ϑ-

E(e)

e²

подстановку E → π - E:

∫

sin E(E - e sin E)

√

dE +

1 - e2 cos2 E

СПИСОК ЛИТЕРАТУРЫ

1. Т. Н. Санникова, К. В. Холшевников, Вестн.

∫

sin E(π - E - e sin E)

СПбГУ. Сер. 1: Математика. Механика. Астроно-

√

dE =

мия № 4, 134 (2013).

1 - e2 cos2 E

2. Т. Н. Санникова, Вестн. СПбГУ. Сер. 1: Математи-

ка. Механика. Астрономия №1, 171 (2014).

∫

sin E(π - 2e sin E)

√

dE.

3. Т. Н. Санникова, К. В. Холшевников, Астрон. журн.

1 - e2 cos2 E

91(12), 1060 (2014).

4. М. Ф. Субботин, Введение в теоретическую

Представляя числитель последнего интеграла в

астрономию (М.: Наука, 1968).

виде

5. Н. Н. Боголюбов, Ю. А. Митропольский, Асимп-

тотические методы в теории нелинейных ко-

-2e +

(1 - e² cos² E) + π sin E,

лебаний (М.: ФМ, 1963).

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019

430

САННИКОВА, ХОЛШЕВНИКОВ

6. А. М. Журавский, Справочник по эллиптиче-

2-е. (М.: Наука, 1976).

ским функциям. (М.-Л.: Изд. АН СССР, 1941).

11. A. Cayley, Mem. Roy. Astron. Soc. 29, 191 (1861).

7. И. С. Градштейн, И. М. Рыжик, Таблицы интегра-

12. M. P. Jarnagin, Astron. Paper 18, 2 (1965).

лов, сумм, рядов и произведений (М.: Физматгиз,

13. К. В. Холшевников, В. Б. Титов, Задача двух тел

1963).

(СПб.: Изд. СПбГУ, 2007).

8. Е. Янке, Ф. Эмде, Ф. Лёш, Специальные функции

14. А. В. Грибанов, Труды АО ЛГУ 38, 165 (1983).

(М.: Наука, 1964).

15. А. Уинтнер, Аналитические основы небесной ме-

9. В. А. Антонов, Е. И. Тимошкова, К. В. Холшевни-

ханики (М.: Наука, 1967).

ков, Введение в теорию ньютоновского потен-

16. Г. М. Фихтенгольц, Курс дифференциального и

циала (М.: Наука, 1988).

10. Г. Н. Дубошин (ред.), Справочное руководство

интегрального исчисления. Т. 3 (М.-Л.: Физмат-

по небесной механике и астродинамике. Изд-е

лит, 1960).

АСТРОНОМИЧЕСКИЙ ЖУРНАЛ том 96

№5

2019