О ротационном движении экранированного ИСЗ в нецентральном гравитационном поле Текст научной статьи по специальности «Физика»

А. А. Тихонов

О РОТАЦИОННОМ ДВИЖЕНИИ ЭКРАНИРОВАННОГО ИСЗ В НЕЦЕНТРАЛЬНОМ ГРАВИТАЦИОННОМ ПОЛЕ*

Рассматривается ИСЗ с электростатически заряженным экраном на орбите, регрессирующей вследствие сжатия Земли. Учитываются вековые возмущения орбиты, вызываемые 2-й зональной гармоникой геопотенциала, т. е. уходы долготы восходящего узла Q с угловой скоростью и аргумента перигея с угловой скоростью кш при неизменной форме орбиты (определяемой фокальным параметром p и эксцентриситетом е) и постоянном наклонении i [1]. Исследуется вращательное движение ИСЗ относительно его центра масс под действием возмущающего момента M, включающего момент Mr гравитационных сил и момент Мл сил Лоренца, возникающих при движении ИСЗ относительно магнитного поля Земли (МПЗ). МПЗ аппроксимируется квадрупольным приближением [2]. Изучается эволюция ротационного движения ИСЗ [3]. Аналогичная задача без учета влияния сжатия Земли рассматривалась в [4] для заряженного ИСЗ, находящегося на кеплеровой круговой экваториальной орбите.

Используются следующие правые ортогональные декартовы системы координат. O3X*Y*Z* (орты т*,т*,т*)—инерциальная система с началом в центре Земли, осью O3Y*, направленной по оси суточного вращения Земли (с угловой скоростью из = изТ*) и осями O3X* и O3Z* в плоскости экватора. CXYZ (орты Ti,Т2,Тз)— «пери-гейная» система, ось CX которой параллельна касательной к орбите в перигее, ось CY — по нормали к плоскости орбиты, ось CZ — параллельна радиус-вектору перигея. Вследствие регрессии орбиты «перигейная» система координат вращается в инер-циальной системе с угловой скоростью к = + кш = + lJvT2, где к^ = Q =

—LVo£(R3/p)2cosi, кш = lUk = 0.5u!oe(R3/p)‘2(5cos2 i — 1), lvо = \Jц/a3 — среднее движение ИСЗ, л — гравитационная постоянная Земли, a — большая полуось орбиты ИСЗ, Rэ — экваториальный радиус Земли, ё « 0.0016239 — безразмерная постоянная, определяемая величиной сжатия Земли [1]. Система CL1L2L3 (орты /ь/2,/3), ось CL3 которой направлена вдоль вектора L кинетического момента ИСЗ. Ориентация системы CL1L2L3 относительно CXYZ определяется углами р и а [3], где р — угол отклонения вектора L от оси CY, а а —угол между осью CZ и линией пересечения плоскостей (X,Z) и (Y, L3,Li). C^-цС (орты ei0, е20, ез0)—орбитальная система, ось C£ направлена по положительной трансверсали к орбите, ось CQ — вдоль радиус-вектора R = O3C = рез0/(1 + е cos v), где v = (Z, () —истинная аномалия точки C. Cxyz (орты ei, в2, ез) —главные центральные оси инерции ИСЗ.

Взаимная ориентация введенных систем координат определяется матрицами направляющих косинусов P = (Pij), K = (Kij), A = (Aij), U = (Uj), T = (Tj) :

з з з з з

Ti ,PijTj, T ''У^,Kijlj, li ''У^,Aijej, T ''У^,Uijej, ei ''У^,Tijej. (1)

j=i j=i j=i j=i j=i

Получим уравнения ротационного движения ИСЗ на регрессирующей орбите. В си-

* Работа выполнена при финансовой поддержке Министерства образования России (грант №Т02-14.0-1443).

© А. А. Тихонов, 2004

лу равенства (Lx, Ly, Lz)T = K(0, 0, L)t, вытекающего из (1), локальная производная L, вычисляемая в системе координат CXYZ, равна K(0, 0, L)t + K(0, 0, L)t и, следовательно, теорему об изменении кинетического момента Lx * = Mx *, L y * = My *, L z* = Mz* можно представить в виде

K(0, 0, L)t + K(0, 0, L)t + к x L = (Mx* , My *, Mz* )T. (2)

Обозначая проекции вектора угловой скорости Шь системы координат CL1L2L3 относительно системы CXYZ на оси Li,L2,L3 соответственно через Шьі, шь2, Шьз, а проекции вектора к на оси Li, L2 —через кьі, кь2, после эквивалентных преобразований приведем векторное уравнение (2) к скалярным уравнениям

L = M3, шь2 = Mi/L — кь2, шьі = —M2/L — кьі. (3)

Поскольку Шь = p + а,

шьі = — àsin р, шь2 = p, шьз = (à cos p (4)

и на основании (3) получаем уравнения для переменных L, р, £ = а + — п/2:

M1 • M2

L = М3, ¿=------b fe sin * sin £, E=-------------1- feísin * eos E ctgp — eos i). (5)

L L sin p

Угол —£ представляет собой дополнение до п/2 угла между линией узлов и проекцией вектора L на плоскость орбиты ИСЗ. Определим ориентацию ИСЗ относительно системы координат CL1L2L3 с помощью s-параметров [5, 6]. Направляющие косинусы вектора L относительно осей X, y, z будут равны соответственно A31 = 4(2sis3+s2 (|s|2 — 1))/w, A32 = 4(2S2S3 — si(|s|2 — 1))/w, A33 = 1 — 8(s1 + s\)/w, где |s|2 = s2 + s2 + s3,

w = (|s|2 + 1)2.

Введем абсолютную угловую скорость ИСЗ

Ш = к + Шь + Шг, (6)

обозначив через Or угловую скорость ИСЗ относительно системы координат CL1L2L3. Кинематические уравнения в s-параметрах имеют вид [5, 6]:

(S1, S2, s3) = 'B(шrx, Шry, urz) ,

1 /si — s2 — s2 + 1 2(sis2 + S3) 2(sis3 — S2) \ (7)

® = — T I 2(sis2 — S3) —s2 + s\ — S3 + 1 2(S2S3 + si) I

\ 2(s1s3 + s2) 2(s2s3 — si) —s1 — s2 + s3 + 1/

Поскольку Lx = Aux = A31 L, Ly = Buy = A32L, Lz = Cuz = A33L, где A, B, C — главные центральные моменты инерции ИСЗ, на основании (6), (7) и равенства (шьі,шь2,шьз)т = A(ubx,uby,ubz)T получаем

Ат(кьі, кь2, кьз)Т + Ат(шьі, Шь2, шьз)Т + B 1(si, s2, s3)T = L(A3i/A, A32/В, A33/C)т.

Разрешая эту систему относительно si, s2, s3 с учетом легко проверяемого равенства BAT = BT, а также принимая во внимание уравнения (4) и (5), получим следующие кинематические уравнения для переменных si, s2, s3:

(¿i, ¿2, S3)T — LB(A3i/A, A32/B, A33/C)T —

— L_ iBT(—М2, M1, M2 ctgp + LkQ sin i cosS/ sinp)T. (8)

Входящие в эти уравнения элементы матриц A и B — известные функции переменных si, ¿2, ¿3• Система уравнений (5), (8) является обобщением соответствующих уравнений, выведенных в [6], и позволяет учесть регрессию орбиты ИСЗ.

Переходя к вычислению моментов Мл и Mdг, будем рассматривать динамически симметричный (A — B) заряженный ИСЗ (с суммарным зарядом Q), центр заряда которого (точка О) находится на оси динамической симметрии Cz и смещен относительно центра масс на величину zo. Если zo — 0, то можно пренебречь фактором градиентности МПЗ (как показано в [4]) и аппроксимировать момент Мл выражением

Мл — QP0 х TT(vc х Be), (9)

где р 0 — CО — zoe3, ve — R + к х R — 1J3 х R — скорость точки C относительно МПЗ, Вc — Bi + B2 —магнитная индукция МПЗ, представляемая в виде суммы дипольной и квадрупольной составляющих. Явные выражения для векторов Bi и 132 приведены в работах [2] и [4], однако вследствие регрессии орбиты в них следует заменить часовой угол на угол ф — (шз — ko)t, отличающийся учетом поправки, обусловленной движением линии узлов. Вводя обозначения pi — vcnBcz — vc^Bcrp Р2 — vc^Bc^ — vc^Bcz, P3 — vc^Bcv — vcvBcz, получи м

MЛj — Qz0 Pi (TiiAj2 — Ti2Aji) (j — 1, 2, 3). (10)

р i=i

Гравитационный момент Mdг с учетом сжатия Земли вычислен в работе [7]. Для динамически симметричного ИСЗ Mrz — 0,

з-54(7(Г-

T32 T33+

а проекция Mry может быть получена из Mrx путем умножения на —1, замены T32 на T31, в2 на ві- Здесь Y£/R = sin i sin и, ві, в2, вз —направляющие косинусы осей x,y, z относительно оси O3Y *, т. е. элементы второй строки матрицы PKA.

В проекциях на оси Li, L2, L3 имеем: Mrj = AjiMrx + Aj2Mry + Aj3Mrz- Таким

образом, найдены проекции Mj = Mлj + Mrj (j = 1, 2, 3).

Для выявления вековых эффектов в эволюции ротационного движения заряженного ИСЗ произведем усреднение системы (5), (8) по быстрым переменным Si,S2,S3-В работе [6] проанализировано конфигурационное s-пространство (обозначим его So) динамически симметричного твердого тела, невозмущенно вращающегося вокруг центра масс с кинетической энергией T и кинетическим моментом L. Показано, что So является правильным тором, форма и размеры которого характеризуются параметром

2 2T/L2 — 1/C

в0: sin 20о = 1 ————, и что усреднение любой непрерывной на So функции от

1 /A — 1 /C

s-параметров по невозмущенному движению можно заменить усреднением по So- В результате такого усреднения получены средние значения проекций момента Ліл:

{Mлl)s = Qzo cos2^o[—pi cos(a — v) + P3 sin(a — v)], {Mлз)s = 0,

{M^)s = Qzo cos 20o[p1 cospsin(<r — v) — p2 sinp + p3 cospsin(<r — v)]

Mrx

и средние значения проекций момента Air, в которых (Мгз)s = 0; (Mri)s и (Mr2)s здесь не приводятся ввиду своей громоздкости.

Подставим усредненный по s-параметрам возмущающий момент в дифференциальные уравнения (5), записав их в безразмерной форме путем перехода от независимой переменной t к новой независимой переменной — средней аномалии т = wot. В результате получим первый интеграл L = Lo = const и два дифференциальных уравнения для медленных переменных р и E:

dp (Mi)s kQ . , . dE {М2)s , kn . ,

— =-------------1-------sinismi, —— = —-—:--------------------------------------------------1-(sm г ctgpcos E — cos г). (11)

ат шoLo шо ат woLo sin р wo

С целью выявления вековых эффектов в эволюции ротационного движения заряженного ИСЗ произведем усреднение системы (11) по переменной т за период Ti = 2п/шо обращения центра масс ИСЗ и по переменной ф за период T2 = 2п/(шз — kn) обращения МПЗ относительно линии узлов. Поскольку (Mi)s и (M.2)s зависят явным образом от V, то усреднение по т заменим усреднением по v, пользуясь известным [1] ат (1 — e2)3/2

соотношением — = -----------------Г7Г. Поэтому, для усреднения функции т[и,ф) имеем

dv (1 + e cos v )2

формулу

1 íTl Г2 „и и (1-е2)3/2 Г2* Г2* f(vA)d.vd4 ^

/ ^,0 = 77777т/ f(i/^)dt1dt2 =--------—2--- , ,2 ■ 12

Í1Í2J0 Jo 4п2 Jo Jo (1 + e cos V )2

В результате, после усреднения на основании (12), дифференциальные уравнения (11) принимают следующий вид:

-ß- = ai cos(E — сОл-) + a,2{uJ-K) sin E + 03(0»^) cos E + r\ sin E+ ат

+ gi cos р sin E + g2 cos р sin(E — 2шп) + g3 sin р sin 2E+

+ g4 sin р sin 2(E — 2шп) + g5 sin р sin 2(E — шп),

dE (13)

— = —a 1 ctg /3SÍn(E — uJk) + a2(uJ1T) ctg pcos E — ctg psin E + aAw-n) +

ат

cos 2р cos 2р

+ r 1 ctg p cos E + Г2 + Об(^л-) cos p + o 1-cos E + 02------cos(E — 2ww) +

sin р sin р

+ g3 cosрcos2E + g4 cosрcos2(E — 2шп) + g5 cosрcos2(E — шп).

Здесь параметры aj (j = 1,4) обусловлены моментом Мл и содержат множителем ко-3Qzq R^

эффициент Kji = -—---------------------^(— 91) eos20o, параметры gj (j = 1,6) обусловлены моментом

2Lo p2

-> А — С I р , З.о \

Мр и содержат множителем коэффициент Кг = —-------, —^ (1-sin 2во), а парамет-

Lo \¡p 2

22 - э • • _ _ ~Лэ ,

ры г\ = —ё—тг sin i cos i и г 2=е—^ cos 8 обусловлены сжатием Земли. p2 p2

Рассмотрим вначале случай экваториальной орбиты (i = 0). На основании (13) получаем уравнения ротационного движения ИСЗ

ар

— = а\ sm <7 — 2<7б sm р sm а cos er,

ат (14)

2

—— = а\ ctg рcos а — 2д$ cos р cos а + (<75 + ge) eos р + ai,

ат

ъг л тс g 2 3

где a i = Kj\e, di = лл-п

9í P

„ _R|,-15ч 2 r. 3

Кгє—(—7— )e , дб = Кг-p2 8 2

R¿ W3 \ 3

p2 wo / 2

-R|

Є 2 ’ S'5 =

p2

ЯЭ 3

1 + Зє—^(l + -є2) p2 2

. Уравнения вида (14) были подробно

исследованы в работе [8], где показано, что они имеют первый интеграл

— <?5 sin2 рcos2 <т + ai sinрcosa — di cosp — — (<75 + ge) cos2 p = const,

позволяющий построить траектории апекса вектора L на сфере L = Lo в системе координат CXYZ. Установлено, что в зависимости от значений параметров вектор L может иметь 2, 4 или 6 стационарных значений. Получены критерии реализации каждого из возможных случаев. Не повторяя полученных в указанной работе результатов, отметим лишь, что если не учитывать сжатие Земли (е = 0) в выражении момента Mг, то g5 =0 и могут реализоваться лишь случаи двух или четырех полюсов, лежащих на меридианах а = 0 и а = п сферы L = Lo (см. рис. 3, 4 из [8]).

Учет сжатия Земли позволяет обнаружить дополнительно существование качественно иных траекторий апекса вектора L, а также появление двух новых полюсов, не лежащих на меридианах а = 0 и а = п (см. рис. 5 из [8], кроме случая 2в), который в данной задаче не реализуется вследствие того, что mi = — (g5 + ge)/(2g5) > 0).

Вместе с тем, малость постоянной е, входящей в параметр g5, приводит к тому, что упомянутые два новые полюса, угловые координаты которых удовлетворяют системе уравнений sin рcos а = ai/(2g5), cos р = —di/(g5 + ge), существуют лишь при условии |ai/g5| ^ 2, накладывающем жесткое ограничение на параметры ИСЗ и начальные условия его вращательного движения. Действительно, неравенство |ai/g5| ^ 2 удовлетворяется либо при |Мл| ^ |Mp|, либо при | cos2$o| ^ 1, что соответствует $о ~ п/2,

т. е. первоначальной закрутке ИСЗ вокруг оси, ортогональной к радиус-вектору 0.

Перейдем к исследованию эволюции ротационного движения заряженного ИСЗ на наклонной орбите. При этом будем предполагать, что невозмущенная орбита ИСЗ является круговой, а ее наклонение — мало и поэтому можно пренебречь малой величиной sin2 i по сравнению с sin i. При этом уравнения (13) будут иметь вид

dp ^ ^ dE cos 2р

— = bi sin Е + <71 cos p sin E, — = 61 ctg p cos E + c\ + <71-----cos E + g§ cos p, (15)

dr dr sin p

где bi = a2 + ri, ci = a4 + Г2, Г2 = sRl^/p2, ri = —Г2 sin i cos i,

к 9Í R3 . .

a-2 = К л -ñ — sin i

g° p

cos і (l + §r2) - -— 2 2w0

gi = Kr3r2 sin i cos i,

a4 = КЛЩ — (l + r2-------^-), <76 = (1 + 3r2).

g0 p uo 2

Отсюда следует, что учет сжатия Земли приводит не только к количественным изменениям в системе (15) (коэффициенты bi и ge), но и вызывает появление слагаемых (с множителем gi), качественно меняющих структуру уравнений (15).

Нетрудно проверить, что уравнения (15) допускают первый интеграл

bi sinрcos Е + <7i sinpcos peas Е — с\ cosр — -да cos2 р = h = const, (16)

с помощью которого можно выразить cos E через р, подставить в систему (15), ввести новую переменную x = cos р и получить для нее дифференциальное уравнение

(dp/ dr)2 = (1 — x2)(bi + <7ix)2 — (h + c\x + ^gex2)2.

Таким образом, задача интегрирования системы (15) сводится к квадратурам, что позволяет найти закон движения вектора L. Вместе с тем, для получения наглядных качественных результатов, характеризующих эволюцию ротационного движения ИСЗ, целесообразно, следуя [3], построить траектории вектора L на сфере L = Lo, вращающейся в абсолютной системе O3X*Y*Z* с угловой скоростью fe. Для построения этих траекторий введем угол а между вектором L и осью CXo, которая лежит в плоскости (X, Z) и направлена ортогонально к линии узлов в «наивысшую» точку орбиты. Поскольку cos а = sin р cosE, то интеграл (16) можно записать в виде

cos а = (h + ci cos р + ^дв cos2 р)/(Ъ\ + д\ cos р), (17)

позволяющем построить семейство кривых cos а(р) для различных значений h. Полюсы вектора L соответствуют особым точкам дифференциальных уравнений (15), определяемым из системы уравнений

cos 2 р

Ъ\ sin Е + oi cos р sin Е = О, Ъ\ ctg р cos Е + с\ + д\--cos Е + Об cos р = О

sin р

или равносильного ей объединения систем (18), (19):

cos 2 р

sin Е = 0, ±(&i ctg р + gi------) + дв cos р + с\ = 0, (18)

sin р

bi + gi cos р = 0, —gi sin р cosE + g6 cos р + c\ = 0. (19)

В результате исследования системы (18) установлено, что полюсы вектора L располагаются на меридианах E = 0 и E = п,а траектории вектора L подобны приведенным на рис. 3 и 4 из [8]. Анализ системы (19) позволяет дополнительно обнаружить возможность существования других двух полюсов вектора L, не лежащих на меридианах E = 0 и E = п, и соответствующих им новых траекторий вектора L. Эта возможность реализуется при выполнении системы неравенств

\bi/gi\ < 1, (gici — bige)2 < gl(gl — bl). (20)

Показано, что система (20) имеет решения при физически реальных и конструктивно допустимых значениях параметров. Угловые координаты указанных полюсов определяются равенствами

g1 c1 — b 1 ge

cos р =-Ъ1/д1, cosE = ±-— , (21)

\ g1\ g12 — b2

что свидетельствует об их симметричном расположении относительно меридианов E = 0 и E = п. Типичный вид траекторий вектора L приведен на рисунке.

Такой характер расположения полюсов и соответствующие им траектории вектора L обусловлены совместным влиянием трех факторов: сжатием Земли, моментом гравитационных сил и моментом лоренцевых сил, причем в последнем моменте решающую роль играет квадрупольная составляющая МПЗ. Действительно, если положить g° = 0, то а2 = 0 и, следовательно, bi/gi = —1/(3ÄT). Полученная величина bi/gi заведомо не удовлетворяет ограничениям (20) в силу малости параметра Кг в условиях ротационного движения ИСЗ и поэтому решения (21) перестают существовать.

Устойчивые полюсы 1, 2, 3, 4 типа «центр» определяются из системы (18), а неустойчивые полюсы 5 и 6 (последний из них расположен симметрично полюсу 5 и поэтому не виден) типа «седло» имеют координаты (21). Анализ интеграла (17) позволил установить, что сепаратрисам, разделяющим зоны влияния полюсов 1, 2 и 3, 4, соответствует значение постоянной h, равное h* = bi(ci — bige/(2gi))/gi. Если h < h*, то изображающая точка движется по траекториям, находящимся в зоне влияния полюсов 1 и 2, а если h > h*, то — в зоне влияния полюсов 3 и 4.

Результаты, полученные в статье, доложены на международной научной конференции по механике «Третьи Поляховские чтения» (СПб., 4-6 февраля 2003 г.).

Summary

A. A. Tikhonov. On the rotary motion of a shielded artificial Earth satellite in a noncentral gravitational field.

The artificial Earth satellite equipped with a charged shield is presented. Its orbit is osculating due to the Earth compression. The secular disturbances, caused by the 2-nd zonal harmonic of geopotential are taken into account. The satellite attitude motion under the action of disturbing moments of gravitational and Lorentz forces is investigated. The Earth magnetic field is taken in quadrupole approximation. The evolution of rotary motion of the satellite is studied on the basis of new differential equations constructed with the use of s-parameters. The basic characteristics of the secular evolution of the shielded satellite rotary motion are revealed.

Литература

1. Справочное руководство по небесной механике и астродинамике / Под ред. Г. Н.Ду-бошина. М.: Наука, 1976. Изд. 2-е.

2. Тихонов А. А., Петров К. Г. Мультипольные модели магнитного поля Земли // Космич. исслед. 2002. Т. 40. №3. С. 219-229.

3. Белецкий В. В. Движение искусственного спутника относительно центра масс. М.: Наука, 1965. 416 с.

4. Тихонов А. А. Уточнение модели «наклонный диполь» в задаче об эволюции ротационного движения заряженного тела в геомагнитном поле // Космич. исслед. 2002. Т. 40, №2. С. 171-177.

5. Marandi S. R., Modi V. J. A preferred coordinate system and the associated orientation representation in attitude dynamics // Acta Astronaut. 1987. Vol. 15. №11. P. 833-843.

6. Петров К. Г., Тихонов А. А. Уравнения ротационного движения твердого тела, основанные на использовании кватернионных параметров // Изв. РАН. Мех. тверд. тела. 2002. №3. С. 3-16.

7. Сарычев В. А. Влияние сжатия Земли на вращательное движение спутника Земли // Искусственные спутники Земли. Изд-во АН СССР. 1961. Вып. 6. С. 3-11.

8. Тихонов А. А. Влияние неоднородности геомагнитного поля на эволюцию ротационного движения заряженного твердого тела // Вестн. Ленингр. ун-та. Сер. 1, 1991. Вып. 2 (№8). С. 90-99.

Статья поступила в редакцию 23 декабря 2003 г.