Предвычисление орбитальных параметров искусственного спутника Земли

14

2. ПРЕДВЫЧИСЛЕНИЕ ОРБИТАЛЬНЫХ ПАРАМЕТРОВ ИСЗ

2.1.Разновидности уравнений движения ИСЗ .

2.1.1. Сначала вспомним вывод дифференциального уравнения второго порядка для описания движения твердого тела в центральноv гравитационноv поле, образованном моделью  однородного земного шара. Здесь можно полагать, что масса М  Земли сосредоточена в начальной точке геоцентрической системы координат рис. 2.1а, а масса m спутника сосредоточена в точке, имеющей искомые координаты x=x(t), y=y(t), z=z(t) конца радиуса-вектора r (рисунок построен для случая  z=z(t)=0) . (Далее вектора  обозначаются жирным курсивом).

По второму закону Ньютона ускорение ИСЗ a=d²r/dt²=dv/dt определяется действующей силой F из равенства F=am. В нашей модели действует сила тяготения F=fMm(-r)/r³=mm(-r)/r³, где f - постоянная всемирнjго тяготения, а m = 3,986008.10¹⁴м³/с²  - гравитационная постоянная Земли. Приравнивание двух выражений для силы F (после сокращения на массу ИСЗ) и приводит к  первому уравнению движения (единицы массы ИСЗ) в векторном виде

a= -mr/r³  или d²r/dt²=dv/dt = -mr/r³                                            (2.1)

2.1.2.Понижение порядка дифференциального уравнения достигается, если левую и правую части (2.1) умножить векторно на r. Справа получим ноль, так как векторное произведение одинаковых векторов [r´r]=0. Оставшееся равенство [a´r]=0 сопоставим с выражением для производной от векторного произведения   d[v´r]/dt=[ar] + [vv] =[a´r]. Можно видеть,что эта производная равна нулю, то есть векторное произведение [v´r] равно некоторому постоянному вектору:

[v´r] =С,          причем модуль                 С=rvsin{rÙv}=2 DS,                        (2.2) где DS – неизменное приращение  площади, «заметаемой» радиусом-вектором в единицу времени – эта площадьтреугольника со сторонами r  и vимеет размерность м²/с. (Заметим, что если единицу времени устремить к нулю, то  обозначение DS необходимо заменить на производную dS/dt.)

Неизменность «заметаемой радиусом вектором  за единицу времени» площади DS была выявлена около 400 лет назад немецким астрономом Иоганом Кеплером, но не аналитически, а в процессе обработки многочисленных наблюдений астронома Тихо-Браге за движением планет вокруг Солнца.

2.1.3. Поскольку вектора rи vперпендикулярны вектору С, не меняющему во времени своего направления, то движение ИСЗ происходит только в одной плоскости. Известно, что траектория такого невозмущенного движения имеет вид эллипса с большой  а  малой  b  полуосями и площадью S_э=pab, причем если величина DS определена, то легко находится и период обращения по эллиптической траектории Т=S_э/DS и средняя угловая скорость (среднее движение) ИСЗ n₀=2p/T=2pDS/S_э, Целесообразно определить формулы связи «заметаемой за секунду» площади DS  с параметрами эллипса.

На рис.2.2 фокус F эллипса совмещен с центром Земли и использована перигейная декартова система координат, в которой ось Ох, совмещенная  с большой полуосью эллипса, направлена в точку перигея, которая находится от центра Земли  на минимальном расстоянии r_п=a-(a²-b²)^0,5 , тогда как точке апогея соответствует расстояние

15

r_A=a+(a²-b²)^0/5. В этих точках векторы rиv взаимно перпендикулярны. Что приводит к полезным далее соотношениям

2DS =v_Пr_Ф = v_Аr_А ,      v_А =v_П (r_П/r_А) < v_П,        r_П +r_A =2a,      r_Пr_А=b²=S_Э²/(pa)²,         (2.3)

где v  и v   - модули векторов скорости ИСЗ.

Теперь применим закон постоянства суммы Э кинетической  Э_к=0.5v²  и потенциальной Э_п=-m/r энергий спутника, движущегося в центральном поле. Из равенства сумм энергий в точках апогея и перигея вытекает, что 0,5v_П²-m/r_П=0.5v_A²-m/r_П. Подставляя сюда v_A из (2.3) получим v_П=(1/r_П)[mb²/a]^0,5 или 2DS=dS/dt=[mb²/a]^0,5, что с учетом выражения для b² из (2.3) приводит к формуле для периода обращения Т=2p(а³/m)^0/5. Заметим. что период обращения при 0<b<a, зависит лишь от большой полуоси эллиптической орбиты .

		Рис. 2.2.  Геометрические  параметры орбиты ИСЗ

 

2.1.3. Перейдем от обобщенных параметров орбиты к определению текущих координат ИСЗ. В декартовой перигейной системе  справедливо известное уравнение эллипса [x+(a²-b²)^0,5]²/a²+y²/b²=1. Более широко используется  полярная система координат (r,q) с центром в фокусе F и полярной осью, направленной на перигей. Здесь q=q(t) – угол между полярной осью (радиусом-вектором перигея) и радиусом-вектором ИСЗ. Этот угол носит веками наименование истинная аномапия. Декартовы и полярные координаты связаны равенствами x=rcosq, y=rsinq. Подставляя их в предыдущее уравнение после громоздких выкладок получают известное уравнение эллипса в полярной системе

16

r=r(t)= (b²/a)/(1+ecosq)=(1-e²)a/(1+ecosq),                                   (2.4)

где e=(1-b²/a²)^0,5 – эксцентриситет, который в СРНС невелик – около 1%. 

2.2. Зависимость истинной аномалии от времени. Формулы Кеплера.

2.2.1. Для предвычисления текущих координат необходимо получить выражения зависимости истинной аномалии q=q(t) от времени, отсчитываемого обычно с момента t_П прохождения перигея. Традиционный для последних веков путь поиска подобной закономерности – составление дифференциального уравнения для искомой величины. Это уравнение для истинной аномалии q=q(t) можно /6/ получить, обратив внимание, что площадь вытянутого треугольника DS на рис. 2.1 выражается (при стремящейся к нулю единице времени) с использованием его высоты h=rDq=r(dq/dt) равеством DS=0,5r²dq/dt. Воспользовавшись полученным в п.2.2 выражением для DS  и формулой (2.4) для r получим искомое нелинейное уравнение

dq/dt = (ma/b²)(1+ecosq)²                                       (2.5)

2.2.2 Лишь  для окружности (е=0, a=b) истинная аномалия является линейной функцией времени q=m/а=2p(t-t_П)/Т. Правая часть этой формулы и в случае эллиптической траектории характеризует отношение «заметаемой» радиусом-вектором ИСЗза интервал t-t_п площади S(t-t_п)=(t-t_п)×DS к площади S_Э эллипса: это отношение также линейно зависит от времени: S(t-t_п)/S_Э=(t-t_П)/Т. Ввиду важности этой дроби для описания  движения ИСЗ и необходимости подчеркнуть ее периодичность вводится умножение на 2p, а величина

М=2p(t-t_П)/Т=2p(t-t_П)(DS/S_Э)                                          (2.6)  

называется «средняя аномалия».

2.2.3. Достойно /6/ удивления и восхищения, что почти 400 лет назад – когда не существовало дифференциального исчисления – искомое аналитическое выражение для q=q(t) – т.е. решение уравнения (2.5) -  было выявлено Иоганом Кеплером и обосновано с помощью простых    геометрических соотношений, соответствующих рис. 2.3, на котором  орбита-эллипс описана окружностью радиусом а, равным большой полуоси эллипса. правый фокус которого совпадает с центром Земли. 

Эллипс можно рассматривать как проекцию круга  радиуса a, повернутого вокруг большой оси эллипса на угол, косинус которого равен отношению b/a полуосей эллипса (см. рис. 2.3). Следовательно, площадь фигуры S_FЭП, «заметаемой» радиусом-вектором за время t - t_п, и площадь фигуры S_FКП соотносятся, как малая  и большая полуоси эллипса, то есть S_FЭП =(b/a)×S_FКП=(b/a)(S_ОКП-S_ОКF).

Центральный угол круга КОП  называется «эксцентрическая  аномалия» и обозначается Е. Из рассмотрения рис. (2.3 ) видно,  что площадь  кругового сектора S_ОКП=pa²(E/2p)=(a²E)/2, а площадь треугльника ОКF, (основание которого (a²-b² )^0.5 равно половине расстояния между фокусами) составит S_ОКF=0.5(a²-b²)^0.5aSinE. Таким образом

                S_FКП=(a²/2){E-[(a²-b²)^0.5/a]sinE} = (a²/ 2)(E-esinE),              S_FЭП=(ab/2)(E-esinE).

17

		Рис.  2. 3.  Взаимосвязь истинной и эксцентрической аномалий

 

Так как (ab)/2=S_э /2p, а 2p (S_FЭП/S_э)=2p(t-t_п)¤T=M, то получаем знаменитое  уравнение 

Кеплера для определения Е по значению М из (2.6)