Космические исследования, 2022, T. 60, № 1, стр. 43-56

1. ВВЕДЕНИЕ

В данной работе рассчитываются низкочастотные (квазистатические) микроускорения на искусственном спутнике Земли типа спутников Бион М-1 и Фотон М-4 [1, 2]. Спутник совершает полет по низкой, почти круговой орбите в режиме солнечной ориентации. Его форма близка к цилиндру, на котором неподвижно закреплены две панели солнечных батарей, расположенные вдоль продольной оси симметрично относительно нее (рис. 1). В режиме солнечной ориентации нормаль к плоскости солнечных батарей направлена на Солнце, продольная ось спутника лежит вблизи плоскости орбиты, угловая скорость мала. Ориентация реализуется с помощью гиросистемы, образованной четырьмя двигателями-маховиками (ДМ). Указанный режим обеспечивает низкий уровень микроускорений и максимальный приток электроэнергии. Недостаток этого режима состоит в необходимости проведения разгрузок кинетического момента гиросистемы, ограничивающих время полета спутника без значительных возмущений. Ниже рассматривается вариант режима солнечной ориентации, при котором накопление суммарного кинетического момента системы ДМ ограничивается за счет управления углом поворота спутника вокруг нормали к светочувствительной стороне солнечных батарей. Исследовано движение спутника в режиме комбинации его гравитационной и солнечной ориентации. С использованием достаточно детальной математической модели движения спутника рассчитаны квазистатические микроускорения, возникающие на борту спутника в указанных режимах.

2. МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ МИКРОУСКОРЕНИЙ

Квазистатические микроускорения на низкоорбитальном спутнике Земли вызываются четырьмя причинами: 1) движением спутника относительно центра масс как твердого тела, 2) градиентом гравитационного поля, 3) аэродинамическим торможением, 4) действием силы, создаваемой органами управления. Если спутник совершает неуправляемое движение или для управления им используются ДМ, то последняя из перечисленных причин исчезает. В таком случае квазистатическое микроускорение в заданной фиксированной точке борта описывается простой формулой, причем чтобы воспользоваться ею, достаточно знать только орбиту и вращательное движение спутника.

Пусть спутник представляет собой твердое тело, и точка $P$ жестко связана с его корпусом. Микроускорением ${\mathbf{b}}$ в точке $P$ называется разность между напряженностью гравитационного поля в этой точке и абсолютным ускорением последней. Роль вектора ${\mathbf{b}}$ в орбитальных экспериментах аналогична роли ускорения свободного падения в экспериментах на поверхности Земли. В частности, если в точке $P$ закрепить пробное тело с исчезающе малой массой ${{m}_{p}},$ то сила реакции, действующая на это тело со стороны спутника, будет равна $ - {{m}_{p}}{\mathbf{b}}.$ Приближенная формула для расчета микроускорений имеет вид [3]

Здесь ${\mathbf{\rho }}$ – радиус-вектор точки $P$ относительно центра масс спутника – точки O, ${\mathbf{\omega }}$ – абсолютная угловая скорость спутника, точка над буквой означает дифференцирование по времени t, μ_E – гравитационный параметр Земли, ${\mathbf{r}}$ – геоцентрический радиус-вектор точки O, $r = \left| {\mathbf{r}} \right|,$ ${\mathbf{v}}$ – скорость этой точки относительно поверхности Земли, ρ_a – плотность атмосферы в точке O, $c$ – баллистический коэффициент спутника. Слагаемые в правой части формулы (1) отвечают первым трем указанным выше причинам возникновения микроускорений.

Формула (1) выведена для общего случая без каких-либо частотных ограничений. Однако если спутник имеет относительно большие инерционные характеристики и его вращательное движение рассчитывается как движение твердого тела (такое движение обычно очень медленное), то формула (1) дает именно квазистатическое микроускорение.

Формула (1) использовалась для расчета реальных квазистатических микроускорений, имевших место на летавших спутниках [3–8]. Ее можно использовать и для прогноза микроускорений [9, 10]. В этом случае составляются уравнения движения спутника, выбирается режим движения, вычисляется решение уравнений движения, моделирующее этот режим, и вдоль найденного решения микроускорение в заданной точке борта рассчитывается по формуле (1). Именно таким образом формула (1) применяется ниже.

3. УРАВНЕНИЯ ДВИЖЕНИЯ СПУТНИКА

Спутник считаем гиростатом, для описания его движения будем использовать три правых декартовых системы координат.

$O{{x}_{1}}{{x}_{2}}{{x}_{3}}$ – связанная система координат, образованная главными центральными осями инерции спутника. Полагаем, что спутник имеет форму прямого кругового цилиндра (рис. 1) радиуса ${{R}_{c}}$ и высоты ${{L}_{c}},$ с двумя прикрепленными к нему одинаковыми прямоугольными пластинами – солнечными батареями суммарной площади S_b. Ось цилиндра совпадает с осью $O{{x}_{1}}.$ Солнечные батареи расположены в плоскости $O{{x}_{1}}{{x}_{3}}$ симметрично относительно оси $O{{x}_{1}},$ стороны батарей параллельны осям $O{{x}_{1}}$ и $O{{x}_{3}},$ ось $O{{x}_{2}}$ перпендикулярна плоскости солнечных батарей. Светочувствительная сторона солнечных батарей обращена к полупространству х₂ > 0. Координаты геометрических центров цилиндра и пластин солнечных батарей обозначим (х_с, 0, 0) и (x_b, 0, 0) соответственно, базисные орты этой системы – e₁,e₂,e₃. Ниже, если не оговорено особо, компоненты векторов и координаты точек относятся к системе $O{{x}_{1}}{{x}_{2}}{{x}_{3}}.$

$C{{y}_{1}}{{y}_{2}}{{y}_{3}}$ – гринвичская система координат. Ее начало находится в центре Земли, плоскость $C{{y}_{1}}{{y}_{2}}$ совпадает с плоскостью экватора, ось $C{{y}_{1}}$ пересекает гринвичский меридиан, ось $C{{y}_{3}}$ направлена к Северному полюсу.

$C{{Z}_{1}}{{Z}_{2}}{{Z}_{3}}$ – квазиинерциальная система координат. Ось $C{{Z}_{2}}$ параллельна вектору кинетического момента орбитального движения спутника, ось $C{{Z}_{3}}$ лежит в плоскости экватора и направлена в восходящий узел оскулирующей орбиты спутника. Базис этой системы образован ортами ${{E}_{1}}{{E}_{2}}{{E}_{3}}.$

Матрицу перехода от системы $C{{Z}_{1}}{{Z}_{2}}{{Z}_{3}}$ к системе $C{{y}_{1}}{{y}_{2}}{{y}_{3}}$ обозначим $C = \left\| {{{c}_{{ij}}}} \right\|_{{i,j = 1}}^{3},$ где ${{c}_{{ij}}}$ – косинус угла между осями $C{{y}_{i}}$ и $C{{Z}_{j}}.$ Элементы этой матрицы выражаются через координаты и компоненты скорости центра масс спутника в гринвичской системе координат. Матрицы перехода от системы $O{{x}_{1}}{{x}_{2}}{{x}_{3}}$ к гринвичской системе координат и системе $C{{Z}_{1}}{{Z}_{2}}{{Z}_{3}}$ обозначим соответственно $\left\| {{{u}_{{ij}}}} \right\|_{{i,j = 1}}^{3}$ и $\left\| {{{a}_{{ij}}}} \right\|_{{i,j = 1}}^{3}.$ Здесь ${{u}_{{ij}}}$ и ${{a}_{{ij}}}$ – косинусы углов, которые образует ось $O{{x}_{j}}$ с осями $C{{y}_{i}}$ и $C{{Z}_{i}}.$ Справедливы соотношения ${{a}_{{ij}}} = {{{\mathbf{E}}}_{i}} \cdot {{{\mathbf{e}}}_{j}},$ $U = CA.$

Положение системы $O{{x}_{1}}{{x}_{2}}{{x}_{3}}$ относительно системы $C{{Z}_{1}}{{Z}_{2}}{{Z}_{3}}$ будем также задавать углами ψ, θ и φ, которые введем следующим образом. Если точку $C$ перенести в точку O, то систему $C{{Z}_{1}}{{Z}_{2}}{{Z}_{3}}$ можно перевести в систему $O{{x}_{1}}{{x}_{2}}{{x}_{3}}$ тремя последовательными поворотами: 1) на угол ψ вокруг оси $C{{Z}_{2}},$ 2) на угол π/2 – θ вокруг новой оси $C{{Z}_{1}},$ 3) на угол φ + π вокруг новой оси $C{{Z}_{2}},$ совпадающей с осью $O{{x}_{2}}.$ Элементы матрицы $A$ выражаются через эти углы с помощью следующих формул:

Уравнения движения спутника состоят из двух подсистем. Одна подсистема описывает движение центра масс спутника в гринвичской системе координат [11]. В ней учитываются нецентральность гравитационного поля Земли и сопротивление атмосферы. Нецентральность поля учитывается с точностью до членов порядка (16,16) включительно в разложении гравитационного потенциала Земли в ряд по шаровым функциям. Атмосфера считается вращающейся вместе с Землей, ее плотность рассчитывается согласно модели ГОСТ Р 25645.166-2004. Параметры атмосферы и баллистический коэффициент спутника считаются неизменными на всем интервале интегрирования уравнений движения.

Другая подсистема описывает движение спутника относительно центра масс (вращательное движение). Она образована уравнениями, выражающими теорему об изменении кинетического момента спутника в его движении относительно центра масс, кинематическими уравнениями Пуассона для элементов первых двух строк матрицы $U$ и уравнениями, описывающими изменение кинетического момента гиросистемы. В уравнениях, выражающих теорему об изменении кинетического момента спутника, учитываются гравитационный и аэродинамический моменты. Для гравитационного момента существует простое аналитическое выражение [12]. Формула аэродинамического момента имеет вид

Здесь ${{v}_{i}}$ – компоненты вектора ${\mathbf{v}}.$ При выводе последней формулы считалось, что молекулы атмосферы при столкновении с корпусом спутника испытывают абсолютно неупругий удар [12, 13], и не учитывалось взаимное затенение корпуса спутника и солнечных батарей от набегающего аэродинамического потока. Такое упрощение оправдано, поскольку для большинства движений спутника относительная продолжительность отрезков времени, на которых указанное затенение существенно, невелика.

Допущения, сделанные при выводе формулы аэродинамического момента, позволяют выписать явное выражение для входящего в формулу (1) баллистического коэффициента спутника. Этот коэффициент имеет вид $c = {{{{S}_{v}}} \mathord{\left/ {\vphantom {{{{S}_{v}}} m}} \right. \kern-0em} m},$ где $m$ – масса спутника, ${{S}_{v}}$ – площадь геометрической фигуры, являющейся проекцией внешней оболочки спутника на плоскость перпендикулярную вектору ${\mathbf{v}}.$ В данном случае

и аэродинамический член формулы (1) принимает вид

Все приводимые в данной работе расчеты микроускорений выполнены по формулам (1), (3).

При выводе выражений для аэродинамического момента и баллистического коэффициента не учитывалось возможное взаимное затенение корпуса спутника и солнечных батарей от набегающего аэродинамического потока. Такое упрощение оправдано, поскольку для большинства движений спутника относительная продолжительность отрезков времени, на которых указанное затенение существенно, невелика.

Заметим, что при расчете микроускорений баллистический коэффициент $c$ считался переменным, а в подсистеме уравнений движения центра масс спутника тот же коэффициент принимался постоянным. Сделанное упрощение оправдано тем, что влияние сопротивления атмосферы на движение центра масс весьма мало и достаточно точно описывается с использованием постоянного значения баллистического коэффициента, которое определяется в результате обработки траекторных измерений.

Подсистема уравнений вращательного движения спутника имеет вид

Здесь символом ${{\tilde {d}} \mathord{\left/ {\vphantom {{\tilde {d}} {dt}}} \right. \kern-0em} {dt}}$ обозначена локальная производная вектора в системе $O{{x}_{1}}{{x}_{2}}{{x}_{3}},$ $\hat {I} = {\text{diag}}({{I}_{1}},{{I}_{2}},{{I}_{3}})$ – тензор инерции спутника в этой системе, ${\mathbf{H}}$ – собственный кинетический момент гиросистемы, ${{{\mathbf{M}}}_{c}}$ – момент, действующий со стороны этой системы, на корпус спутника, ${{{{\omega }}}_{{\text{E}}}}$ – угловая скорость вращения Земли, ${{{\mathbf{u}}}_{1}}$ и ${{{\mathbf{u}}}_{2}}$ – первая и вторая строки матрицы перехода U. Третья строка этой матрицы вычисляется по формуле ${{{\mathbf{u}}}_{3}} = {{{\mathbf{u}}}_{1}} \times {{{\mathbf{u}}}_{2}}.$ Строки ${{{\mathbf{u}}}_{1}}$ и ${{{\mathbf{u}}}_{2}}$ связаны условиями ортогональности матрицы $U$ (${{{\mathbf{u}}}_{i}}$ – орты осей $C{{y}_{i}}$), которые учитываются при задании начальных условий для этих переменных.

Чтобы замкнуть подсистему уравнений вращательного движения, надо добавить к уравнениям (4) выражение для ${{{\mathbf{M}}}_{c}}.$ Оно приведено ниже.

В расчетах использовались следующие параметры описанной модели. Параметры спутника: $m$ = 6440 кг, ${{I}_{1}}$ = 2600 кг/м², ${{I}_{2}}$ = 11 100 кг/м², ${{I}_{3}}$ = = 10 900 кг/м², ${{R}_{c}}$ = 1.3 м, ${{L}_{c}}$ = 5.0 м, ${{S}_{b}}$ = 33 м², ${{x}_{b}}$ = = –1 м, ${{x}_{c}}$ = 0.3 м. Микроускорения рассчитывались в точке $P$ с координатами (–1, 0.7, 0.5 м) (рис. 1). Эта точка находится на внутренней стенке рабочего отсека спутника, примерно на ее середине. Вблизи этой точки возможна установка научной аппаратуры. Параметры модели атмосферы: ${{F}_{{10.7}}}$ = ${{F}_{{81}}}$ = 150, ${{A}_{p}}$ = 12.

Рассматривались два варианта начальных условий движения центра масс спутника. Соответствующие им решения уравнений орбитального движения назовем орбитами I и II. Начальные условия орбиты I задавались в восходящем узле орбиты в момент времени 07.13.07 UTC 5.V.2013. При этом элементы орбиты составляли: высота в апогее 575.2 км, высота в перигее 546.8 км, наклонение 64.87°, аргумент широты перигея – 124.65°, долгота восходящего узла в гринвичской системе координат 343.27°. Начальные элементы орбиты II имели такие же значении, но относились к моменту времени 07.13.07 UTC 21.XII.2013. В случае орбиты I в некоторой точке отрезка времени, на котором моделировалось движение спутника, максимальное значение угла $\vartheta $ между ортом направления “Земля–Солнце” и плоскостью орбиты спутника достигало 47°, для орбиты II максимальное значение угла $\vartheta $ достигало 88°. Таким образом, орбиты I и II по-разному расположены относительно Солнца.

Начальные условия уравнений (4) задавались в тот же момент времени, что и начальные условия орбитального движения. Этот момент служил началом отсчета времени – точкой $t = 0.$

4. СОБСТВЕННЫЙ КИНЕТИЧЕСКИЙ МОМЕНТ СИСТЕМЫ ДМ

Для определенности положим, что гиросистема спутника состоит из четырех ДМ. Векторы ${\mathbf{H}}$ и ${{{\mathbf{M}}}_{c}}$ реализуются в виде суммы четырех векторов, отвечающих соответствующим ДМ. Полагаем, что оси вращения ДМ расположены параллельно боковым ребрам правильной четырехгранной пирамиды. Высота пирамиды параллельна оси $O{{x}_{1}},$ линии пересечения граней пирамиды с плоскостью $O{{x}_{2}}{{x}_{3}}$ параллельны или перпендикулярны осям $O{{x}_{2}},$ $O{{x}_{3}}$ (рис. 2). Орты ${{{\mathbf{g}}}_{k}}$ осей вращения ДМ имеют компоненты

которые параметризируем углами α и β: ${{d}_{1}} = \cos {{\alpha ,}}$ ${{d}_{2}} = \sin {{\alpha }}\sin {{\beta ,}}$ ${{d}_{3}} = \sin {{\alpha }}\cos {{\beta }}{\text{.}}$

Суммарный кинетический момент системы ДМ выражается формулой

где ${{G}_{k}}$ – кинетический момента ДМ с номером k. Значения углов α и β выберем так, чтобы реализуемые системой ДМ максимальные угловые скорости спутника относительно каждой из осей $O{{x}_{i}}$ $(i = 1,2,3)$ были примерно одинаковы [14]. Тогда

При заданном значении ${\mathbf{H}}$ соотношения (5) нельзя единственным образом разрешить относительно ${\mathbf{G}}.$ Для достижения единственности потребуем, чтобы вектор ${\mathbf{G}}$ имел минимальную евклидову норму. Тогда

Здесь ${{D}^{ + }} = {{D}^{T}}{{(D{{D}^{T}})}^{{ - 1}}}$ – псевдообратная матрица матрицы $D$ [15]. Выписанные соотношения (6) задают собственный кинетический момент каждого ДМ. В данной работе предполагается, что предельно допустимые значения собственного кинетического момента каждого ДМ лежат в пределах $\left| {{{G}_{k}}} \right| \leqslant 20$ Нмс $(k = 1,2,3,4).$

5. РЕЖИМ ТРЕХОСНОЙ СОЛНЕЧНОЙ ОРИЕНТАЦИИ СПУТНИКА

Сперва рассмотрим режим трехосной солнечной ориентации спутника, в котором $O{{x}_{2}}$ неизменно направлена на Солнце, ось $O{{x}_{1}}$ лежит в плоскости орбиты, абсолютная угловая скорость спутника мала. Пусть ${\mathbf{s}}$ – орт направления “Земля–Солнце”. Орт ${\mathbf{n}} = {{({\mathbf{s}} \times {{{\mathbf{E}}}_{2}})} \mathord{\left/ {\vphantom {{({\mathbf{s}} \times {{{\mathbf{E}}}_{2}})} {\left| {{\mathbf{s}} \times {{{\mathbf{E}}}_{2}}} \right|}}} \right. \kern-0em} {\left| {{\mathbf{s}} \times {{{\mathbf{E}}}_{2}}} \right|}}$ лежит в плоскости орбиты спутника и ортогонален орту ${\mathbf{s}}.$ Изменение ортов ${\mathbf{s}}$ и ${\mathbf{n}}$ в системе $O{{x}_{1}}{{x}_{2}}{{x}_{3}}$ описывается уравнениями

Закон управления гиросистемой ${{{\mathbf{M}}}_{c}} = {{{\mathbf{M}}}_{c}}({\mathbf{\omega }},{\mathbf{s}},{\mathbf{n}})$ обеспечивающий затухание возмущенного движения спутника в окрестности положения ${\mathbf{\omega }} = 0,$ ${{{\mathbf{e}}}_{1}} = {\mathbf{n}},$ ${{{\mathbf{e}}}_{2}} = {\mathbf{s}}$ с требуемой скоростью имеет вид [16]: где $\hat {W} = {\text{diag}}(1,1,\sqrt 2 ),$ ξ – положительный параметр.

Ниже принято, что ξ = 0.01 с^–1. Полагаем, что закон управления (7) формируется в соответствии с показаниями установленных на спутнике датчиков ориентации и угловой скорости.

Моделирование режима трехосной солнечной ориентации спутника сводилось к численному интегрированию системы (4), (7) с использованием соотношений (5), (6). Компоненты орта ${\mathbf{s}}$ в гринвичской системе координат рассчитывались по приближенным формулам [17]. Начальные условия системы (4), (7) в момент времени $t$ = 0 задавались следующим образом. Спутник в гринвичской системе координат занимает положение ${{{\mathbf{e}}}_{1}} = {\mathbf{n}},$ ${{{\mathbf{e}}}_{2}} = {\mathbf{s}}.$ Кинетический момент системы ДМ ${\mathbf{H}}$(0) = 0, соответственно ${\mathbf{G}}$(0) = 0. Начальные значения компонент угловой скорости ${{{{\omega }}}_{i}},$ задавались с учетом ошибок в их реализации: ${{{{\omega }}}_{1}}(0) = {{{{\omega }}}_{2}}(0) = {{{{\omega }}}_{3}}(0) = 0.01\,\,{{{\text{град}}} \mathord{\left/ {\vphantom {{{\text{град}}} {\text{c}}}} \right. \kern-0em} {\text{c}}}.$ Результаты расчетов движения спутника, полученные в рамках приятой модели для рассматриваемых орбит приведены на рис. 3, 4 (орбита I) и 5, 6 (орбита II). На рисунках приведены графики зависимости от времени углов ψ, θ, φ и величин собственного кинетического момента каждого из ДМ ${{G}_{k}}.$ Все графики на рис. 3–6 построены на интервале времени 13 сут. Переходной процесс (процесс гашения возмущенного движения), обусловленный ошибками в задании начальной угловой скорости спутника, длится менее 10 мин и из-за масштаба на рисунках не виден.

Результаты моделирования показывают, что для двух рассматриваемых вариантов орбит использование закона управления (7) обеспечивает трехосную солнечную ориентацию спутника и затухание его возмущенного движения в окрестности положения ${\mathbf{\omega }} = 0,$ ${{{\mathbf{e}}}_{1}} = {\mathbf{n}},$ ${{{\mathbf{e}}}_{2}} = {\mathbf{s}}.$ Ошибки солнечной ориентации определяются значением угла ${{\sigma }} = \arccos ({\mathbf{s}} \cdot {{{\mathbf{e}}}_{2}})$ (θ = $\vartheta $ при σ = 0), в данном случае σ < 0.4°. Амплитуды установившихся колебаний компонент угловой скорости ограничены значениями:

|ω₁| < 1.0 ⋅ 10^–4 град/c, |ω₂| < 2.0 ⋅ 10^–3 град/c,

|ω₃| < 7.0 ⋅ 10^–4 град/c.

Микроускорения на борту спутника в рассматриваемой точке $P$ по модулю не превышают значения 3.5 ⋅ 10^–6 м/с². Указанный уровень микроускорений является приемлемым для проведения космических экспериментов в области микрогравитации. В работе [18] показано, что наибольший вклад в общий уровень микроускорений ${\mathbf{b}}$ вносит составляющая ${{{\mathbf{b}}}_{g}}$ от действия гравитационного поля Земли, а наименьший – составляющая ${{{\mathbf{b}}}_{r}}$ от вращательного движения спутника относительно центра масс.

На рис. 4 и 6 представлены графики величин ${{G}_{k}}$ $(k = 1,2,3,4),$ на которых видно, что в процессе поддержания солнечной ориентации спутника для всех рассматриваемых орбит кинетический момент каждого из ДМ существенно превысил свое предельно допустимое значение. Разгрузка системы ДМ с помощью магнитных катушек или реактивных двигателей в данной работе не предполагается, поскольку необходимость их проведения сокращает отрезки невозмущенного полета спутника. По этой причине целесообразно рассмотреть возможность использования закона управления, позволяющего не только обеспечить затухание возмущенного движения спутника в окрестности положения покоя с требуемой скоростью, но и дополнительно увеличить продолжительность отрезков полета спутника без проведения разгрузок гиросистемы. Чтобы ограничить рост кинетического момента ДМ не нарушая режим солнечной ориентации спутника в работах [16, 19] предложен соответствующий закон управления, результаты численного моделирования которого будут представлены ниже.

6. РЕЖИМ СОЛНЕЧНОЙ ОРИЕНТАЦИИ СПУТНИКА С ОГРАНИЧЕНИЕМ РОСТА ГИРОСТАТИЧЕСКОГО МОМЕНТА

Рассмотрим возможность использования закона управления, позволяющего не только обеспечить затухание возмущенного движения спутника в окрестности положения ${{{\mathbf{e}}}_{2}} = {\mathbf{s}}$ с требуемой скоростью, но и дополнительно ограничить рост накапливаемого кинетического момента системы ДМ. Чтобы не нарушать режим солнечной ориентации спутника накопление кинетического момента системы ДМ будем контролировать только за счет управления углом φ поворота спутника вокруг орта ${{{\mathbf{e}}}_{2}} = {\mathbf{s}}.$ Закон управления гиросистемой обеспечивающий поддержание режима солнечной ориентации с заданными выше условиями имеет вид [16]:

где χ, ${{{{\kappa }}}_{i}}$ $(i = 1,2,3)$ – положительные постоянные, ${{{{\gamma }}}_{i}}$ – компоненты вектора ${\mathbf{r}},$ ${{K}_{i}}$ – компоненты вектора кинетического момента спутника в его движении относительно центра масс.

Покажем, что выбранный закон изменения управляющего момента (8) действительно обеспечивает солнечную ориентацию спутника и при этом ограничивает накопление кинетического момента системы ДМ. С этой целью вычислим решение системы (4), (8) с использованием соотношений (5), (6), приняв χ = 2ξ, ξ = 0.01 с^–1, ${{{{\kappa }}}_{1}}$ = = ${{{{\kappa }}}_{2}}$ = 1 (Нмс)^–1, ${{{{\kappa }}}_{3}}$ = 3 (Нмс)^–1 с начальными условиями в момент времени $t$ = 0 заданными как в п. 5. Результаты расчетов движения спутника для рассматриваемых орбит показывают, что использование закона управления (8) обеспечивает солнечную ориентацию спутника и затухание его возмущенного движения в окрестности положения ${{{\mathbf{e}}}_{2}} = {\mathbf{s}}.$ Графики зависимости от времени углов ψ, θ совпадают с соответствующими графиками представленными на рис. 3, 5. Точность солнечной ориентации спутника σ < 0.4°. Переходной процесс обусловленный ошибками в задании начальной угловой скорости спутника, в данном режиме длится менее 20 мин. После окончания переходного процесса графики угла ${{\varphi }}(t)$ представляют собой установившиеся колебания с доминирующей частотой $4{{f}_{0}} \approx 0.00072$ Гц, где ${{f}_{0}} = {{{{{{\omega }}}_{0}}} \mathord{\left/ {\vphantom {{{{{{\omega }}}_{0}}} {2{{\pi }}}}} \right. \kern-0em} {2{{\pi }}}}$ – орбитальная частота (${{{{\omega }}}_{0}}$ – среднее движение центра масс спутника, для рассматриваемых в данной работе орбит ${{{{\omega }}}_{0}}$ ≈ 1.09 ⋅ 10^–3 с^–1). Амплитуда установившихся колебаний угла φ ограничена значением $\left| {{\varphi }} \right|$ < 6°. В случае движения спутника по орбите II при значении $\vartheta $ > 70° происходит увеличение амплитуды колебаний угла φ до значений 0° < φ < 15°, максимум которых достигается при максимальном значении угла $\vartheta $ ≈ 88° [16]. В дальнейшем, по мере уменьшения угла $\vartheta ,$ амплитуды колебаний угла φ снижаются до значений $\left| {{\varphi }} \right|$ < 6°. Амплитуды установившихся колебаний компонент угловой скорости ограничены значениями:

|ω₁| < 1.2 ∙ 10^–4 град/c, |ω₂| < 2.2 ∙ 10^–3 град/c, |ω₃| < 7.0 ∙ 10^–4 град/c.

Микроускорения на борту спутника в рассматриваемой точке $P$ по модулю не превышают значения 3.8 ∙ 10^–6 м/с². Указанный уровень микроускорений является приемлемым для проведения космических экспериментов в области микрогравитации.

На рис. 7, 8 на интервале времени 13 сут приведены графики зависимости от величин собственного кинетического момента каждого из ДМ ${{G}_{k}}$ $(k = 1,2,3,4).$ В случае движения спутника по орбите I величины ${{G}_{k}}$ остаются ограниченными значением $\left| {{{G}_{k}}} \right|$ < 10 Нмс $(k = 1,2,3,4)$ (см. рис. 7) на всем интервале движения, что в два раза меньше максимально допустимого значения $\left| {{{G}_{k}}} \right|.$ В случае движения спутника по орбите II, когда Солнце относительно плоскости орбиты поднимается достаточно высоко ($\vartheta $ > 70°$)$, значения $\left| {{{G}_{k}}} \right|$ начинают возрастать (см. рис. 8), достигая своего максимума $\left| {{{G}_{k}}} \right|$ ≈ 48 Нмс при $\vartheta $ ≈ 88°. В дальнейшем, по мере уменьшения $\vartheta ,$ происходит разгрузка системы ДМ до значений $\left| {{{G}_{k}}} \right|$ ≤ 10 Нмс.

7. СПЕКТРАЛЬНЫЙ АНАЛИЗ МИКРОУСКОРЕНИЙ

С целью анализа уровня микроускорений, возникающих на борту спутника в процессе поддержания солнечной ориентации, в работе определены характерные частоты колебаний компонент вектора микроускорений ${\mathbf{b}}$. Частоты были найдены с помощью спектрального анализа, выполнявшегося по следующей схеме [20]. Пусть ${{x}_{n}}$ $(n = 1,\;2,\; \ldots ,\;N)$ – значения какой-либо переменной $x(t)$ исследуемого решения в узлах равномерной временной сетки $\{ {{t}_{n}}\} {\text{:}}$ ${{x}_{n}} = x({{t}_{n}}).$ Во всех рассмотренных ниже примерах шаг сетки $h = {{t}_{{n + 1}}} - {{t}_{n}} = 10$ с. Предположим, что исследуемая функция имеет вид

где $f_{m}^{ \circ } \in (0,{{{{h}^{{ - 1}}}} \mathord{\left/ {\vphantom {{{{h}^{{ - 1}}}} 2}} \right. \kern-0em} 2})$ и ${{{{\alpha }}}_{0}},$ ${{{{\alpha }}}_{m}},$ ${{{{\beta }}}_{m}}$ $(m = 1,\;2,\;\, \ldots \;\,,\;M)$ – постоянные параметры, причем среди частот $f_{m}^{ \circ }$ нет одинаковых. Оценить частоты и амплитуды отдельных гармоник, входящих в $x(t),$ можно, исследуя максимумы периодограммы Шустера на интервале $0 < f < {{{{h}^{{ - 1}}}} \mathord{\left/ {\vphantom {{{{h}^{{ - 1}}}} 2}} \right. \kern-0em} 2}.$ Функция $I(f)$ имеет много максимумов, из которых выбирают несколько наиболее выделяющихся. Если функция $I(f)$ имеет такой максимум в точке ${{f}_{ * }},$ то предполагается, что ${{f}_{ * }}$ близка одной из частот $f_{m}^{ \circ },$ а величина ${{2\sqrt {I({{f}_{ * }})} } \mathord{\left/ {\vphantom {{2\sqrt {I({{f}_{ * }})} } N}} \right. \kern-0em} N}$ является оценкой амплитуды $\sqrt {{{\alpha }}_{m}^{2} + {{\beta }}_{m}^{2}} $ соответствующей гармоники. Периодограмму Шустера удобно преобразовать к виду, называемому амплитудным спектром $A(f) = 2{{N}^{{ - 1}}}\sqrt {I(f)} .$ Выделяющиеся максимумы функции $A(f)$ напрямую оценивают амплитуды отдельных гармоник, но ее максимумы выражены менее наглядно максимумов периодограммы.

На рис. 9 показаны амплитудные спектры ${{A}_{{bi}}}(f)$ величин суммарных микроускорений ${{b}_{i}}(t)$ $(i = 1,2,3),$ в левой части рисунка в случае использования закона управления (7), в правой – при использовании закона (8). Спектры приведены в диапазоне частот от 0 до 0.0015 Гц, при значениях $h = 10$ с и $N = 8640.$ Как и в случае использования закона управления (7), наибольший вклад в общий уровень микроускорений ${\mathbf{b}}$ на борту спутника при использования закона (8) вносит составляющая ${{{\mathbf{b}}}_{g}}$ с доминирующей частотой $2{{f}_{0}} \approx 0.00036$ Гц. Такие колебания возбуждаются градиентом гравитационного поля. При использовании закона (8) на всех рассматриваемых орбитах движения спутника происходит увеличение амплитуд колебаний величин ${{b}_{1}}(t)$ и ${{b}_{3}}(t)$ на частоте $4{{f}_{0}} \approx 0.00072$ Гц, что соответствует доминирующей частоте колебаний величин угла φ и ${{{{\omega }}}_{2}}.$ Однако такое увеличение амплитуд практически не сказывается на общем уровне микроускорений ${\mathbf{b}}$, возникающих на борту спутника.

8. РЕЖИМ КОМБИНАЦИИ ГРАВИТАЦИОННОЙ И СОЛНЕЧНОЙ ОРИЕНТАЦИИ СПУТНИКА

При движении спутника по орбите II в режиме его солнечной ориентации происходит рост величин $\left| {{{G}_{k}}} \right|$$(k = 1,2,3,4)$ до значений, превышающих максимально допустимое значение $\left| {{{G}_{k}}} \right|$ примерно в два раза (см. рис. 8). Это происходит вследствие уменьшения абсолютной величины действующего на спутник гравитационного момента при достаточно большом значении угла θ. Чтобы ограничить рост величин $\left| {{{G}_{k}}} \right|,$ можно ограничить допустимую величину угла ориентации θ используя следующий сценарий вращательного движения спутника с помощью закона управления (8). До наступления условия θ > η спутник стабилизируется относительного положения ${{{\mathbf{e}}}_{2}} = {\mathbf{s}},$ при выполнении условия θ > η спутник стабилизируется относительно положения ${{{\mathbf{e}}}_{2}} = {\mathbf{\tilde {s}}},$ где η – некоторое наперед заданное значение угла между ортом ${\mathbf{\tilde {s}}} = ({{\tilde {s}}_{1}},\;{{\tilde {s}}_{2}},\;{{\tilde {s}}_{3}})$ и плоскостью $C{{Z}_{1}}{{Z}_{3}}.$ Компоненты орта ${\mathbf{\tilde {s}}}$ в системе координат $C{{Z}_{1}}{{Z}_{2}}{{Z}_{3}}$ определяются по следующим формулам:

Здесь ${{s}_{i}}$ – компоненты вектора ${\mathbf{s}}$ в системе координат $C{{Z}_{1}}{{Z}_{2}}{{Z}_{3}}.$ Вычислим решение системы (4), (8) с использованием соотношений (5), (6), по описанному выше сценарию. Примем η = 70°, значения остальных параметров закона управления (8), а также начальных условий в момент времени $t$ = 0 примем как в п. 7. Результаты расчетов движения спутника, полученные в рамках принятой модели для орбиты II приведены на рис. 10, 11, и представлены графиками зависимости от времени углов ψ, θ, φ, а также величин собственного кинетического момента каждого из ДМ ${{G}_{k}}$ $(k = 1,2,3,4).$ Графики построены на интервале времени 13 сут. Переходной процесс, как и в п. 7, длится менее 20 мин. Результаты моделирования показали, что поддержание режима гравитационно-солнечной ориентации спутника позволяет более чем в два раза уменьшить значения $\left| {{{G}_{k}}} \right|$ < 20 Нмс $(k = 1,2,3,4),$ рис. 11. При этом уровень микроускорений на борту спутника в указанной выше точке $P,$ по прежнему, по модулю не превышает значения 3.8 ∙ 10^–6 м/с². Происходит снижение освещенности панелей солнечных батарей спутника примерно на 5% [18], но поскольку при значении $\vartheta $ > 70° спутник не попадает в тень Земли, то суммарный энергосъем будет близок к максимальному.

ЗАКЛЮЧЕНИЕ

В работе показано, что c помощью системы ДМ можно реализовать режим солнечной ориентации спутника с уровнем квазистатических микроускорений, приемлемым для проведения космических экспериментов в области микрогравитации. Стабилизацию режима солнечной ориентации можно обеспечить без накопления кинетического момента ДМ за счет управления углом поворота спутника вокруг направления на Солнце.