1. Определение квазистатической составляющей микроускорения

по телеметрической информации

Полет МКС после пристыковки к ней Служебного модуля (СМ) практически неизменно проходил в ориентированном состоянии. В основном поддерживалась неизменная ориентация  станции либо в орбитальной системе координат, либо в абсолютном пространстве. Первые несколько месяцев поддержание ориентации и ее изменение выполнялись с помощью реактивных двигателей СМ. Затем эти двигатели стали применяться, как правило, только для выполнения разворотов и разгрузки используемых для поддержания ориентации американских гиродинов. Поскольку система управления ориентацией работала непрерывно, в БЦВМ имелась полная информация о вращательном движении станции. Несколько раз в сутки часть этой информации передавалась на Землю для разного рода проверок. По специальному заказу можно организовать сбор этой информации во время проведения научных экспериментов.

Используемая телеметрическая информация представляет собой значения в дискретные моменты времени кватерниона ориентации, вектора угловой скорости станции, а также моменты включения и продолжительность работы реактивных двигателей системы управления ориентацией. При пролете станции вблизи одного из наземных приемных пунктов  можно получить указанные значения с шагом около 1 с в течение не более 10 мин. Это так называемый режим непосредственного сброса информации. К сожалению, получаемые данные относятся только к участку пролета МКС над наземным приемным пунктом. Измерительную информацию можно также записывать в течение продолжительного времени (обычно около витка) в специальное запоминающее устройство и только потом в удобной ситуации сбрасывать на Землю. Однако при этом шаг по времени перечисленных выше величин составляет 1.5 – 2 мин. Для расчета квазистатической составляющей микроускорения такие продолжительные отрезки информации с редкими по времени данными наиболее удобны.

Самым ценным видом данных являются значения кватерниона. Они наиболее точны и достаточны для расчета квазистатической составляющей. Методики расчета описаны в [1] и успешно применялись еще в случае станции «Мир». Данные угловой скорости менее точны, особенно, если движение станции происходит с малой угловой скоростью [2]. Эти данные целесообразно использовать для контроля. Применение остальных видов данных требует привлечения специальных и довольно сложных моделей функционирования органов управления станцией. В данной работе применялись простейшие модели движения станции, основанные на кинематических уравнениях [1].

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

.                  (1)

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

Восстановление фактического вращательного движения станции выполнялось по значениям кватерниона, задающего ее ориентацию относительно инерциальной системы координат. Пусть на каком-либо временном интервале моменты со значениями такого кватерниона расположены достаточно часто. Эти значения сглаживаются кватернионной функцией времени, имеющей непрерывную вторую производную. Дифференцированием этой функции находятся  и , в результате получается аппроксимация вращательного движения станции. Движение центра масс последней считается кеплеровым. Элементы этого движения определяются по данным траекторных измерений. Знание движения центра масс станции и ее движения относительно центра масс позволяет рассчитать первые три слагаемых формулы (1). Вычисление четвертого слагаемого требует еще задания способов вычисления  и . Первая из этих величин рассчитывается согласно модели атмосферы [4], вторая считается постоянной и определяется по данным траекторных измерений.

В качестве инерциальной системы координат, относительно которой задается ориентация станции, используется правая система , связанная с земным экватором. Точка  – центр Земли, ось  направлена в точку весеннего равноденствия эпохи 2000.0, ось  направлена в соответствующий северный полюс мира. Под ориентацией станции понимается ориентация жестко связанной с ее корпусом строительной системы координат . Эта система также правая. Точка  центр масс станции, ось  параллельна продольной оси СМ и направлена от его переходного отсека к агрегатному отсеку, ось  перпендикулярна оси вращения солнечных батарей СМ. Положение системы  относительно системы  задается с помощью кватерниона , имеющего единичную норму:  Матрицу перехода от системы  к системе  обозначим , где  – косинус угла между осями  и . Элементы этой матрицы выражаются через компоненты  с помощью известных формул , ,  и т. п. Ниже компоненты векторов и координаты точек указываются в системе .

Телеметрическая информация собирается на временном интервале не более нескольких часов и содержит последовательность моментов времени и кватернионов

.                           (2)

Здесь ,  – значение кватерниона , вычисленное на момент . Как правило,  мин, . Кватернион, задающий ориентацию станции, определен с точностью до знака. Знаки  и момент  в (2) выбираются из условия

.

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

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

В [5] отыскание  в предположении, что эта функция дважды непрерывно дифференцируема, сводится к решению вариационной задачи

,          .                     (3)

Здесь  – заданное положительное число. Решением задачи (3) является кубический сплайн. В [5] приведена программа на алголе-60 вычисления коэффициентов этого сплайна по величинам , ,  . Эта программа, переписанная на турбо-паскаль, применялась в данной работе.

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

,

,

.

Продифференцировав последние уравнения по времени и подставив в полученные выражения первую и вторую производные нормированной кватернионной функции, можно найти . Теперь все готово для расчета микроускорения по формуле (1).

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

,      ,

,      ,          (4)

,      ,      .

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

,      .        (5)

Здесь  – положительный параметр.

Смысл функционала (5) состоит в следующем. Для кватернионной функции , удовлетворяющей уравнениям (4) и аппроксимирующей данные измерений (2), выполнены соотношения  . Положения станции, задаваемые кватернионами  и , близки, поэтому можно считать, что положение  получается из положения  бесконечно малым поворотом c вектором . Этот факт выражается соотношениями

,

,

,

,

с помощью которых  в (5) можно преобразовать к виду

.

Следовательно,  характеризует близость расчетного и фактического движений станции. Последняя формула служит также обоснованием способа аппроксимации фактического движения посредством покомпонентного сглаживания последовательности кватернионов (2).

Второе слагаемое в правой части  представляет собой штраф за слишком большое угловое ускорение станции, которое может возникнуть при точной аппроксимации. Таким образом, функционал (5) является компромиссным. Этот функционал похож на функционал, возникающий при решении изопериметрической вариационной задачи (3) методом неопределенных множителей Лагранжа [5]. В [5] множитель Лагранжа, аналогичный параметру  в (7), выбирается из условия выполнения второго (изопериметрического) условия (3). При численном решении новой вариационной задачи  считался независимым параметром и уменьшался до тех пор, пока  в (5) продолжало монотонно убывать [1]. Решение задачи, отвечающее минимальному значению , считалось аппроксимацией фактического движения станции.

Примеры реконструкции вращательного движения станции по информации (2) и расчета на основании этой реконструкции квазистатической компоненты микроускорения приведены на рис. 1 – 4. На рис. 1, 3 правые графики иллюстрируют сглаживание телеметрических данных (маркеры) решениями уравнений (4) (непрерывные кривые) в рамках второго способа, графики слева представляют зависимость от времени компонент и модуля микроускорения в точке (5 м, 2 м , 0). Детали реконструированных движений представлены на рис. 2, 4. Здесь слева непрерывными кривыми изображены графики углов Крылова [3], задающих ориентацию системы  относительно системы  (рис. 2) и орбитальной системы координат [1] (рис. 4), маркерами указаны значения тех же углов, рассчитанные по данным (2). Способ введения углов в данном случае не важен. Важны только диапазоны изменения углов, которые характеризуют точность поддержания ориентации. В средней части рис. 2 и 4 приведены графики компонент угловой скорости станции. Непрерывные кривые – результат реконструкции, маркеры – данные измерений. Справа на рисунках приведены графики компонент углового ускорения станции, горизонтальные прямые указывают нуль.

Рисунки 1 – 4 иллюстрируют движение станции во время поддержания  гиродинами неизменной ориентации в инерциальной (рис. 1, 2) и орбитальной (рис. 3, 4) системах координат. Это достаточно мягкое движение – в данных примерах не было включений реактивных двигателей системы ориентации. Для сравнения приведем пример более сложной ситуации.

Рисунки 5, 6 иллюстрируют режим поддержания неизменной ориентации станции в абсолютном пространстве с помощью реактивных двигателей системы управления. Эти рисунки аналогичны рис. 1 – 4, но они получены первым способом, и для задания ориентации системы  относительно системы  на них использованы другие углы. Пики на графиках микроускорений (рис. 5) стали несколько резче, но уровень микроускорений практически не изменился. На рис. 7 приведены результаты реконструкции того же движения по той же измерительной информации с помощью динамической модели вращательного движения станции, учитывающей работу реактивных двигателей. Содержание этого рисунка в точности аналогично содержанию левой и средней частей на рис. 6. Сравнение обоих рисунков показывает, несколько худшую точность аппроксимации данных измерений динамической моделью. Не столь уж большую, как кажется на первый взгляд из рисунка, всего лишь в несколько десятых долей градуса. Заметим, что в первом способе суммарная погрешность аппроксимации исходных данных –  в (3) – входной параметр. В данном примере при сглаживании каждой компоненты кватерниона принималось .

Главное отличие рис. 6 и 7 – в графиках угловой скорости. Вертикальные отрезки графиков компонент угловой скорости на рис. 7 – это скачки, вызванные срабатыванием реактивных двигателей системы управления. Формально этим скачкам в угловой скорости отвечают бесконечные значения (дельта-функции) углового ускорения и линейного ускорения. Фактически же в моменты срабатываний двигателей на борту станции возникают конечные, но весьма значительные возмущения. Величина углового ускорения достигает 0.0006 с. При правильном расчете это должно привести к узким пикам на графиках микроускорений, в 500 – 1000 раз превышающим пики на рис. 5. При применении методик [1] к имеющейся измерительной информации соответствующие импульсы длительностью 0.1 с распределяются по отрезку времени более 1 мин (интервал между соседними значениями кватерниона). В результате узкие пики на графиках компонент микроускорения скрадываются.

Расчет микроускорения в последнем примере содержит еще одну неточность, которая тесно связана с только что указанной. Некоторые связки двигателей ориентации СМ не сбалансированы, и при их включении к станции прикладывается не чистый момент (пара сил), а еще и сила. Эта сила имеет импульсный характер и может привести к изменению указанных выше пиков.

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

2. Измерение вибрационной составляющей микроускорения акселерометрами ИМУ-128

Вибрационная составляющая микроускорения на борту МКС измеряется с помощью 10 трехкомпонентных акселерометров ИМУ-128. Пять из них расположены на Функциональном грузовом блоке (ФГБ) и пять – на СМ. Основные технические характеристики этого прибора следующие:

амплитудный диапазон измеряемых ускорений –       0.1 – 0.0001 м/с,

частотный диапазон измеряемых ускорений –            0.01 – 20 Гц,

порог чувствительности –                                              0.00002 м/с,

скорость выборки –                                                         50 отсчетов в сек.

Для увеличения точности измерений в области малых амплитуд измерения по каждой компоненте передаются по двум телеметрическим каналам. По одному каналу передается полный диапазон амплитуд 0.01 м/с, по другому – в 10 раз меньший – 0.001 м/с. Из-за «ухода нулей» более точный канал не работал ни у одного из датчиков. Менее точный канал работал примерно у 4 однокомпонентных составляющих прибора ИМУ-128 на ФГБ. На СМ в штатном режиме полета станции работал только один "грубый" канал. Причиной этого могли быть вибрации от работы расположенного на СМ компрессора системы кондиционирования воздуха. Поэтому был проведен сеанс измерений при выключенном компрессоре. Это позволило получить информацию по девяти "грубым" (0.01 м/с) и трем "точным" (0.001 м/с) каналам.

Полученные данные измерений представляют собой отрезки временных рядов длительностью несколько минут. Анализ этих данных включал решение двух задач: 1) выделение циклических трендов, 2) оценивание спектральной плотности ряда остатков, получающегося из исходного ряда удалением значимых циклических трендов. Ниже описаны методы решения этих задач, использованные в данной работе. Обоснование и более подробное изложение этих методов можно найти в [6 – 8].

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

Циклическим трендам отвечают скачки в спектральной функции ряда. Если такие тренды удалить, то спектральную функцию ряда остатков можно считать абсолютно непрерывной и вместо этой функции рассматривать ее производную – спектральную плотность. Поиск циклических трендов выполнялся на основе анализа периодограммы ряда [6,7].

Пусть  – натуральное четное число,   – данные измерений какой-либо компоненты микроускорения, полученные в моменты времени , . Периодограммой называется функция , рассматриваемая на отрезке , ,и определенная соотношениями

,       .

Здесь  – пробная частота,  – частота Найквиста. Если моменты времени  выражены в секундах, то  и  измеряются в герцах. В рассматриваемом случае  с,  Гц.

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

Обычно периодограмма представляет собой пилообразную функцию с большим числом локальных экстремумов. Расстояния между ее соседними максимумами . Разыскивая циклические тренды, следует выбирать только значимые максимумы периодограммы, которые существенно больше ее соседних максимумов. Для выделения значимых максимумов использовался критерий Шустера [6,7]. Правило его применения состоит в следующем. Над данными измерений выполняется дискретное преобразование Фурье

.

Имеют место соотношения (напомним )

,

позволяющие свести поиск значимых максимумов периодограммы к изучению величин . А именно, находятся локальные максимумы последовательности  , т. е. элементы , удовлетворяющие неравенствам , . Локальный максимум  считается значимым, если выполнено неравенство

,

где  – приемлемая для исследователя вероятность совершить ошибку. В данной работе принималось . Выписанное неравенство представляет собой асимптотику критерия Шустера при , . Оценка частоты циклического тренда, отвечающего значимому локальному максимуму , вычислялась по формуле

.

Эта формула задает абсциссу вершины квадратной параболы, проходящей через точки ,  и . Найденная оценка затем уточнялась при построении наилучшей в смысле метода наименьших квадратов аппроксимации данных измерений выражением . Здесь , ,  и  – параметры. Значения параметров находились методом Гаусса - Ньютона.

После определения значимых частот  , отвечающий им тренд искался в виде

,                                    (3)

где  – постоянные величины. Значения постоянных определялись методом наименьших квадратов из условия наилучшей аппроксимации выражением (3) данных измерений.

После того, как тренд из исследуемого отрезка ряда удален, следует провести спектральный анализ остатков. Ряд остатков можно считать рядом с непрерывным спектром и оценивать его спектральную плотность. Для получения оценки спектральной плотности использовалось дискретное преобразование Фурье. Пусть последовательности чисел  и   представляют собой соответственно ряд остатков и его дискретное преобразование Фурье. Справедливо соотношение (равенство Парсеваля)

.

Периодограмма ряда остатков в точках  равна , а оценками значений спектральной плотности в этих точках могут служить величины  . Однако эти оценки несостоятельны. Чтобы получить состоятельную оценку спектральной плотности, следует разбить отрезок  на несколько равных интервалов и выполнить усреднение указанных значений по этим интервалам. Найденные средние значения будут служить оценкой спектральной плотности на соответствующих интервалах. В результате спектральная плотность оценивается кусочно-постоянной функцией. А именно, пусть натуральное число  – делитель числа : . Отрезок  делим на  интервалов одинаковой длины, и на интервале   спектральную плотность  ряда остатков оцениваем постоянной

.

Заметим, что в силу равенства Парсеваля найденная описанным способом оценка  спектральной плотности удовлетворяет соотношению

.

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

.

В данной работе использовалось окно Ханна [6 – 8]

.

Графики оценки спектральной плотности  не дают наглядного представления о распределении амплитуд микроускорения по частотам. Вместо этих графиков строились графики функции . Функция  – также кусочно-постоянная. Ее значение на каждом отрезке постоянства представляет собой среднюю квадратическую на этом отрезке частот амплитуду колебаний с непрерывным спектром.

Данные измерений микроускорения на ФГБ вдоль оси  и их амплитудные спектры приведены на рис. 8, 9. Эти рисунки устроены следующим образом. Верхние графики представляют собой ломаные, ординаты вершин которых – данные измерений с удаленными средними значениями, а абсциссы образуют равномерную сетку с шагом 0.02 с. Вторые сверху – графики функции . Они также представляют собой ломаные. Их вершины находятся в точках ,  . Они построены при  для данных, чьи начальные отрезки представлены на верхних графиках. Использование функции  вместо периодограммы  удобно по той причине, что значения  в точках максимума являются оценками амплитуд соответствующих циклических трендов. Недостаток функции  состоит в том, что ее значимые максимумы выражены менее наглядно значимых максимумов периодограммы.

Анализ значимых максимумов  и последующий поиск циклических трендов методом наименьших квадратов позволили обнаружить циклические тренды со следующими частотами  (выражены в герцах) и амплитудами  (выражены в 0.001 м/с): , , , , , , ,  в случае данных представленных на рис. 8; , , , , , , ,  в случае данных представленных на рис. 9. Оба тренда с частотами  могли быть вызваны упругими колебаниями корпуса станции или действиями экипажа. Остальные тренды были вызваны функционированием бортовых систем.

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

Как видно из рисунков, периодограмма имеет значимые максимумы при  Гц. Однако акселерометр ИМУ-128 имеет резонанс на частоте 25 Гц и при  Гц его амплитудная характеристика отличается от 1. Циклический тренд с указанной частотой был удален из первоначальных данных, но он не имеет отношения к анализу микрогравитационной обстановки на станции.

Данные измерений микроускорений на СМ вдоль осей , ,  и их амплитудные спектры приведены на рис. 10 – 12. Эти рисунки аналогичны рис. 8, 9 и построены точно таким же способом. Приведенные данные были получены, когда не работали установленный на СМ компрессор системы кондиционирования воздуха и некоторые другие агрегаты. По этой причине уровень измеренных микроускорений оказался достаточно мал. В данных измерений были найдены следующие частоты и амплитуды циклических трендов: , , , , ,  в случае данных, приведенных на рис.. 10; , , , , ,  в случае данных на рис. 11; , , , , , в случае данных на рис. 12. Здесь использованы те же единицы измерения величин  и , что и выше. Найденные тренды оказались несущественными. Их исключение из данных практически не повлияло на амплитудные спектры (ср. графики функций  и  на рисунках). Это означает, что микроускорения на СМ в то время имели непрерывный спектр. Микроускорения с таким спектром могли быть вызваны слабыми упругими колебаниями корпуса станции и деятельностью экипажа.

Как показала обработка полученной информации, квазистатическая составляющая на расстоянии не более 5 м от центра масс станции не превышает обычно 0.00001 м/с, вибрационная составляющая на ФГБ не превышала 0.01 м/с в диапазоне частот от 0 до 20 Гц. Во время спокойной обстановки на СМ вибрационная составляющая микроускорения также не превосходила указанного предела. При включенном компрессоре системы кондиционирования воздуха уровень микроускорений на СМ повышался более чем на порядок.

Литература

1.     Сазонов В.В., Беляев М.Ю., Ефимов Н.И., Стажков В.М., Бабкин Е.В. Определение квазистатической компоненты микроускорения на станции «Мир». Космические исследования, 2001, т. 39, № 2, с. 136-147.

2.     Бабкин Е.В., Беляев М.Ю., Ефимов Н.И., Сазонов В.В., Стажков В.М. Неуправляемое вращательное движение станции «Мир». Космические исследования, 2001, т. 39, № 1, с. 27-42.

3.     Сазонов В.В., Комаров М.М., Полежаев В.И., Никитин С.А., Ермаков М.К., Стажков В.М., Зыков С.Г., Рябуха С.Б., Асеведо Х., Либерман Е. Микроускорения на орбитальной станции "Мир" и оперативный анализ гравитационной чувствительности конвективных процессов тепло-массопереноса. Космические исследования, 1999, т. 37, № 1, с. 86-101.

4.     Модель верхней атмосферы для баллистических расчетов. ГОСТ 22721-77. М., Изд-во стандартов, 1978.

5.     Reinsch C.H. Smoothing by spline functions. Numerische mathematik, 1975, B. 24, No. 5, S. 383-393.

6.     Хеннан Э. Анализ временных рядов. М., Наука, 1964.

7.     Теребиж В.Ю. Анализ временных рядов в астрофизике. М., Наука, 1992.

8.     Отнес Р., Эноксон Л. Прикладной анализ временных рядов. М., Мир, 1982.