CC BY
7УДК 629.785:523.31
Влияние возмущающих гравитационных сил, связанных с нецентральностью гравитационного поля Земли, на эволюцию орбиты
космического аппарата
А.В. Ващенко
В работе исследуется возмущенное движение космического аппарата (КА). Описывается влияние возмущающих гравитационных сил, связанных с нецентральностью гравитационного поля Земли, на эволюцию орбиты КА. В работе также представлены результаты расчетов эволюции околокруговой орбиты КА.
При анализе эволюции орбиты необходимо учитывать большое число возмущающих факторов, действующих на КА во время полета. К их числу относятся возмущающие гравитационные силы, связанные с нецентральностью гравитационного поля Земли, и возмущающие аэродинамические силы от сопротивления земной атмосферы. Такие возмущающие силы, как силы тяготения Луны, Солнца и других небесных тел, силы светового давления, электромагнитные силы и реактивные, связанные с сублимацией материалов, также воздействуют на КА, однако они не учитываются. Это связанно с тем, что на рассматриваемых диапазонах высот орбит (до 600 км) их влияния на движение КА незначительно.
При баллистических расчетах исследуется эволюция околокруговых орбит, учитывается влияние возмущающих гравитационных сил, связанных с «неправильностью» формы Земли и возмущающие аэродинамические силы от сопротивления земной атмосферы. Исследуются орбиты высотой 600 км и наклонением орбиты i от 0 до 180°. Расчет производится с помощью программы баллистических расчетов (на базе математической системы компьютерной алгебры Maple 7).
В расчетах используется система дифференциальных уравнений в оскулирующих элементах с модификацией - вместо времени прохождения перицентра используется аргумент широты. Данная система может использоваться для анализа движения любого типа (эллипс, парабола, гипербола), но она имеет особенности при е=0 и i=0. Для исключения особенности е=0 используются две новых формальных переменных (составляющие вектора Лапласа): l(t)=e(t)sin(ro(t)) и q(t)=e(t)cos(ro(t)).
Метод оскулирующих элементов сводится к тому, что исследование возмущенной траектории КА может быть сведено к анализу совокупности невозмущенных траекторий, соответствующих каждому моменту времени. Т.е., при анализе возмущенного движения мы можем считать, что в любой момент времени КА находится на дуге конического сечения с определенными значениями элементов орбиты, а коническое сечение при движении КА изменяется. При этом изменяются элементы орбиты [1]:
р=рО); е=е(0; ю=ю(t); i=i(t); и=иО).
Элементы орбиты, рассматриваемые как функции времени, через которые координаты и составляющие скорости в возмущенном движении выражаются теми же формулами, что и в невозмущенном движении, называются оскулирующими элементами [1].
Обозначения: р(^ - параметр, е^) - эксцентриситет, - агрумент широты перицентра, О^) -долгота восходящего узла, ОД - наклонение, и(^ - аргумент широты, А^) - большая полуось, ОД -радиус, ОД - истинная аномалия, _S(t), _Т^), W(t) проекции возмущающих ускорений на оси орбитальной системы координат (по радиус-вектору, по нормали к нему в плоскости орбиты и по нормали к плоскости орбиты).
После всех преобразований система дифференциальных уравнений движения в оскулирующих элементах при исключении особенностей е^-0 принимает вид:
Рисунок1 Элементы возмущенного движения
st01 :
дt Р(t) 2 1 + q(t}ео8(и(t))+1(t) sin(u(t))
V
_Т(t) (-я(t) sin(u(t)) cos(u(t))- 2 1(t) - 2 sin(u(t)) + 1(t) cos(u(t))2)
1 + t) t)) + 1(t) sin(u(t))
1 + t) cos(u(t)) + 1(t) sin(u(t))
зЮ3 := д q(t) = _Б(t^рПГ t))
t) соз(и(t))2 + q(t) + 1(t) 81п(и(t)) соз(и(t)) + 2 соз(и(t))) _Т(t)
+ 1 + q(t) соз(и(t)) + 1(t) sin(u(t))
1(t) _W(t) cot(i(t)) sin(u(t)) д/ртр" + 1 + я(t)соз(и(t)) + 1(t) sin(u(t))
_ а -W(' )sin(u(г))
:_ г) = (1 + q(() соз(и(г)) + 1(г) sin(u(г))) sin(i(г))
8г05 :_ А,г) = ^/ртг -W(г)cos(u(г))
5 : дг1( ) 1 + q(0^(и(г)) + 1(г) sin(u(г))
а
згдб :_ — и( г) =
-W(г) р(г)2 со1(^г)) sin(u(г)) (1 + q(г) ^(и(г)) + 1(г) sin(u(г)))3т
1 , . ~. л . ч ,,■ , Л Л". I (1 + q( г) ^(и( г)) + 1( г) sin(u( г )))2
(з/2)
р( г)
д \/ А(г)2(l + q(г)«оКи(г)) + 1(г) sin(u(г)))_Т(г)
807 :_ — А( г) = 2 у рп-—-
дг р(г)
2 л/рп) А(г)2 (sin(u(г)) q(г)- cos(u(г)) 1(г)) г)
р(0
8г08 :_ дг Г8( г) = 0
:_ — £9( г) = 0 дг
8г010 :_ ддг А0( г) = 0
Данная система дифференциальных уравнений используется при исследовании эволюции орбиты КА в программе баллистических расчетов.
Возмущения от гравитационного поля Земли.
Если бы Земля имела форму шара со сферическим распределением плотности, то гравитационная сила, действующая на КА на орбите, была бы та же, что и в задаче двух тел. Однако форма Земли значительно сложнее. Прежде всего, Земля имеет полюсное сжатие, которое влияет на гравитационное ускорение. Коэффициент полюсного сжатия Земли имеет следующее значение:
а = (гс-Кз.п.)/Кз.э~0,0033523,
где гс - средний экваториальный радиус, равный 6378,245 км.;
Rз.п. - полярный радиус, равный 6356,863 км.
Потенциал гравитационного поля Земли как сжатого сфероида равен:
и = ^-Щ ¿(3зт2 р-1), г 3г
где Ат - гравитационный параметр Земли, г - расстояние от центра Земли до КА, р - широта точки, 8 = 66,07*10Л3 кмл2 [1].
Однако для точных расчетов орбиты спутника необходимо учитывать аномалии земного потенциала. Соответственно потенциал гравитационного поля Земли может быть записан, как [1]:
т г- Ыпк п г- А
и = ^ (1 - е /,0 (г-у Р(и яп(ф)) + е е 1п*(г-)п (1 - *т(ф)2) 2 diff (Р(п, яп(ф)), яп(ф)к) со8(к(1 - )
Г 1=2 Г п=2 к=1 Г
где гс - средний экваториальный радиус Земли, г,ф - геоцентрический радиус, геоцентрическая широта и долгота, РОЖф))- полином Лежандра 1-го порядка,
(к)
(1 -Б1п(ф)2)2 diff (Р(п,81п(ф)),81п(ф)к)соб(к(1-1п,к)) - присоединенная функция пк-го порядка,
выраженная через к-ю производную полинома Лежандра п-го порядка, 1,,о, 1п,к , 1п,к - постоянные, характеризующие фигуру Земли.
Зональные составляющие:
Потенциал гравитационного поля Земли для зональных составляющих принимается в форме
[1]:
п г-/т(У, I, о(—УР(1,8т(ф)))
Тессериальные и секториальные составляющие:
Потенциал гравитационного поля Земли для тессериальных и секториальных составляющих принимается в форме [1]:
Nnk n 1k
(S (ë (-)n (! - sin(^)2)2 diff (P(n, sin(f)), sin(f), k) cos(k(A - A))))
Unk(r, ф, A) = n=2 k=1-r-
r
значения тессеральной гармоники для n=10 и k=4.
значения секториальной гармоники для n=10.
Подставляя соотношение для долготы Х_1+й-ю31 в выражение для тессериальных и секториальных составляющих получим:
(S In.k (~)" (1 - sin(f)2)2 diff(F(n, sin(f)), sin(f), k) cos(kl) cos(k(Ank - Q + ®3t)) + sin(kl) sin(k(Ank - Q + ®3t)))))
^ ^ )n (1 - sin(f)2)2 Unk (r t) =-n=^_k=l r
г
Используя геометрические соотношения sin(ф)_sin(u)sin(i)
. sin(м)cos(/) cos(u)
Sin(/) = I C0S(/) = I 2=4
-у/1 - sm(u) sin(/) у/1 - sm(u) sin(/)
представим силовую функцию в виде:
иЕ(г,иДА1)_Шо(г,ид)+Шк(г,иДА0
И проекции ускорений на соответствующие направления через элементы орбиты:
— ™ • —ЦЕ(т, и, i,w, г) — иЕ(г, и,/, О, г) = —Ж(г,и,О,г), ЕТ = —и У 7 , ЕЖ = — У__
г г sin(u)
Постоянные силовой функции (потенциала) принимаются для йп_398603,2 гс_6378,165.
По данным из [1] и [2]:
Коэффициенты 1п0 зональных гармоник потенциала:
1п0
М=1/1000000
120 130 140 150 160 170
-1082,47т 2,57т 1,84т 0,06т 0,39т 0,47т
180 190 0,02т 0,11т
Коэффициенты 1пк тессериальных и зональных гармоник потенциала
1пк к=1 к=2 к=3 к=4
т=1/1000000 п=2 0,07т 1,77т 0 0
т=1/1000000 п=3 1,95т 0,38т 0
0,22т
т=1/1000000 п=4 0,66т 0,15т 0,06т 0,007т
Фазы Хпк тессериальных и зональных гармоник потенциала
Хпк к=1 к=2 к=3 К=4
п=2 1,086 -0,255 0 0
п=3 0,1327 -0,363 0,0398 0
п=4 3,7856 0,5027 -0,6981 0,424
Результаты баллистических расчетов
На рисунке 5 показана зависимость изменения долготы восходящего узла в зависимости от наклонения орбиты КА.
Характеристики космического аппарата следующие:
Масса КА на опорной орбите 1950 кг.
Наклонение орбиты от 0 до 180°
Высота орбиты 600 км
Характерная площадь аппарата 20 м2
80 60 Зависимость измен
Рисунок 5 Зависимость изменения долготы восходящего узла от наклонения орбиты КА.
Как видно из графика угол долготы восходящего узла изменяется, то есть происходит поворот плоскости орбиты вокруг оси Земли, называемый прецессией узла орбиты. С наклонением орбиты близким к 90° прецессия орбиты уменьшается. Этот свойство орбиты необходимо учитывать при анализе обеспечения КА энергией, получаемой от солнечных батарей. При больших изменениях угла долготы восходящего узла КА будет чаще оказываться в тени Земли. Соответственно можно выбрать орбиту так, чтобы КА всегда был освещен Солнцем. В таком случае орбита КА будет поворачиваться каждые сутки вместе с Землей на ~1°.
На рисунках 6 и 7 показаны зависимости изменения большой полуоси и долготы восходящего узла от массы аппарата и наклонения орбиты.
Зависимость
□
□
□
□
Рисунок 6 Зависимость изменения большой полуоси от массы КА и наклонения орбиты
Рисунок 7 Зависимость изменения долготы восходящего узла от массы КА и наклонения орбиты
Список использованной литературы:
1. Дубошин Т.Н. Справочное руководство по небесной механике. М.: Изд.«Наука», 1971.-431 с
2. Дубошин Т.Н. Небесная механика (основные задачи и методы). М.: Изд.«Наука», 1975.-385 с
3. Константинов М.С., Каменков Е.Ф., Перелыгин Б.П., Безвербый В.К. Механика космического полета. М.:Машиностроение, 1989.-408 с.
4. Гушин В.Н. Основы устройства космических аппаратов. М.:Машиностроение, 2003.-272 с.
Сведения об авторе
Ващенко Алексей Викторович, аспирант кафедры космические системы и ракетостроение Московского авиационного института (государственного технического университета); email: AlekseyVash@yandex.ru
