УДК 629.78 : 681.51
аналитический синтез программного движения космических аппаратов наблюдения
© 2004 Е.И. Сомов1, С.А. Бутырин2
1 Центр исследований устойчивости и нелинейной динамики института машиноведения РАН, г. Москва 2 НИИ проблем надежности механических систем, г. Самара
Многоаспектная проблема гиросилового управления пространственным движением маневрирующих космических аппаратов (КА) наблюдения более 20 лет находится под пристальным вниманием авторов [1-7]. Представляется решение двух актуальных задач аналитического синтеза программного углового движения КА - пространственного поворотного маневра с краевыми условиями общего вида и маршрутного движения. Это позволяет вычислять потребные изменения всех координат силовых гироскопов по аналитическим соотношениям и обеспечивать слабое возбуждение колебаний упругих элементов конструкции КА.
Введение
В гиросиловых системах управления движением (СУД) КА наблюдения традиционно используется координатно-временной метод углового наведения. Программа маршрутного движения КА, где выполняется оптико-электронная съемка, получается на основе векторного сложения всех элементарных движений в геодезическом базисе при задании орбитального движения КА, движения Земли, начальных координат и азимута сканирования с учетом зенитного угла Солнца и текущей перспективы наблюдения исходя из главного требования: изображение заданного участка Земли должно перемещаться требуемым образом на фотоприемной структуре в фокальной плоскости телескопа.
Для перехода между маршрутными движениями КА совершает поворотные маневры с краевыми условиями общего вида, рис. 1. Слабое возбуждение упругих колебаний конструкции КА в процессе таких маневров достигается за счет аналитически рассчитанного [4] гиросилового управления при заданном программном движении корпуса КА.
Рис. 1. Схема поворотных маневров КА наблюдения
Задача аналитического синтеза программы поворотного маневра
Задача аналитического синтеза программы поворотного маневра (ПМ) КА на заданном интервале времени ( е Тр = [(Р, ],
= + Тр, состоит в определении явных функций времени: кватерниона ориентации А(?) связанного с корпусом КА базиса в относительно известного инерциального ба-
зиса I, векторов угловой скорости ю (г), ускорения е (г) и его производной б (г) = е* (г) + ю (t) х Е(г ). Кватернион Л = (Л,0,X),X = (Л,{}, векторы ю = (ш } , б = (вг-} = со и вектор производной ускорения
Е = (вi} = б* + ю х б должны удовлетворять следующим краевым условиям на левом ( t = гр) и правом ( t = гр ) концах траектории ПМ:
Л(0 = Л о; ю (О = ю 0- е0ш Юо;
Е (го^ = Ео = ео 8о ;
Л(гр ) = Лу; ю (гр) = ю у = е00 Юу ;
(1)
Е (гг) =Е г =е /8 г;
Е (гр) = Еу = еВ* 8* + юу х еу. (2)
Последнее условие в (2) представляет требования к гладкости сопряжения ПМ с последующим участком маршрутного движения (МД) КА.
Подход к решению этой задачи основывается на необходимом и достаточном условии разрешимости классической задачи Дарбу - определения Л(г) в аналитическом виде (в том числе в квадратурах) из уравнения
Л(г) = 1 Л(г)©ю (г) при известных Ло и
ю (г). Введем базис Ео, фиксированный в
инерциальном базисе I кватернионом Л о
(1), и подвижные базисы Е к (к = 1,...п), где
базис Е п совпадает со связанным базисом в. Необходимое и достаточное условия разрешимости задачи Дарбу состоят в возможности представления вектора угловой скорости ю (г) в виде
ю (г) = ю п^) + ю пЛ(г) + •••+ю (), где вектор ю к(г) имеет неизменное направление в базисе Е к1 и является вектором угловой скорости базиса Е к относительно базиса Е к1, т.е. в виде
ю (г) = ю п4(0+Л п(0©(ю П-2(0+Л п-1(?)©(ю п:2(?)+Л пМ)©>
©(• + ю \(г)) • ©Л п-2(г) )©Л п-1(0 )©Л а(г) '
где вектор-столбец ю к-1(г), к = 1,...п составлен из проекций вектора ю к(г) фиксированного направления в базисе Ек-1, а Л к(г) является кватернионом ориентации базиса Е к
относительно базиса Е к1.
Решение поставленной задачи представляется как результат сложения в общем случае шести одновременно происходящих элементарных поворотов "вложенных" базисов
Ек вокруг ортов ек осей Эйлера, положение которых определяется из краевых условий (1), (2) исходной пространственной задачи. Краевые условия всех 6 элементарных поворотов относительно ортов е к приведены в табл. 1, где для первых пяти элементарных
Таблица 1. Краевые условия
к Тип элементарного движения Краевые условия на левом конце траектории (г = гр) Краевые условия на правом конце траектории (г = гр)
ф ш в в ф ш в в
1 Гашение ускорения во о о во вО ф 1 о о о
2 Гашение скорости ш о о шо о в ° ф 2 о о о
3 Позиционный переход о о о • о в 3 f ф 3 о о о
4 Разгон до скорости ш у о о о в 4 ф 4 ш у о о
5 Разгон до ускорения в у о о о в ? ф 5 о в у о
6 Разгон до производной ускорения в у о о о в 6 ф 6 о о * 8/
движений требуется обеспечить равенство нулю локальной (собственной) производной ускорения на правом конце траектории.
Кватернион А(() ориентации КА в базисе I определяется произведением А(0 = А„©Л^)©А2(()©Аз(/)©А4(/)©А5(/)©А«(/). (3) Здесь индексы 1-6 кватернионов А к(() соответствуют их номерам в табл. 1, причем кватернион
Ак(() = (^(Фк(()/2), eк • sin( Фк(()/2)), где ф к(() и e к - текущий угол и орт оси Эйлера к -го поворота. В силу неподвижности орта e к в базисе E к-1 имеем
ю к(() = ф к(() •e к; ё к(()=ф к(() •e к и
ё к(() =фк(() •e к. Вектор угловой скорости ю (() , векторы углового ускорения ё(() и его производной ё(() при начальных обозначениях векторов
ю (1)(() = ю 1((), ё(1)(() = ё 1((), ё(1)(() = ё1(() определяются аналитически по рекуррентному алгоритму:
для верхних индексов к = 2:6 последовательно вычисляются
ю к(() =: А к(( )©ю (к> )©А к(();
ю (к)(() = ю к(() + ю к((); ё ^ (() =: Ак(() ©ё (к-1)(() ©А к((); ё(к)(() = ёк(() + ё ^ (() + ю к(() х ю к((); ё ^ (() =: Ак(() ©ё (к-1)(() ©А к(();
ё(к)(()=ё к(()+ё5 (()+(2ё5 (()+ю5 ^) х
х ю к(()) х ю к(() + ю 5 (() х ёк(() ,
(4)
двух движений выберем функции ф к((). к=1,2 с краевыми условиями
Рк(0) =0; Фк(0) =шю0|;
фк(^) = ¿0 =| ё Фк(<рг) = 0; Фк((р) = 0; Фк(<р) = 0
0
(5)
искомые векторы получаются как ю(() = ю(6)((), ё(() = ё(6)((), ё(() = ё(6)(() .
Функции ф к((), представляющие в аналитическом виде углы элементарных поворотов, выбираются в классе полиномов (сплайнов) соответствующей степени.
Гашение начальных значений угловой скорости и углового ускорения. Для первых
в виде сплайнов фк(т) 5-ой степени нормированного времени т = ((- (0р)/ Тр е [0,1]:
фк(() = ¿к(() = -6 (а0- ат + 2а2т 2)/Тр;
фк(/) = ек(/) = е0 -т (6а0 - 3а1т + 4а2т2);
фк(() = ш к(() = ®с + Трт [ес -т (3а0 -а1т + а2т')];
фк(/) = Трт{&0 + Трт [е0 - т(20а0 - 5а1т + 4а2т2)/10]/2},
(6)
где коэффициенты а.. определяются так:
а0 =2 ш 0/Тр + е0; а =8 ш 0Т + 3 е0 и а2 = 3 ю0 /Тр + е0.
Здесь формально гашение начального углового ускорения получается по соотношениям (6) при ш 0 = 0, а гашение начальной угловой скорости - при в 0 = 0. На левом конце траектории производные ускорений
вк((0) = вк = фк((0) = -6а0/Тр, а углы фк на
правом конце траектории принимают значения фк((р) - фк = Тр[8(юг/Тр) + ег ]/20 .
Разгон до заданных значений угловой скорости и ускорения. Для четвертого и пятого движений функции ф к(( ),к = 4,5 должны удовлетворять краевым условиям
Фк((Р) = 0; Фк ((0р) = 0;
Фк((р) = 0;
Мр) = ю/|; Мр) = 8 ё/ \ Фк((р) = 0,
(7)
поэтому фк(() также выбираются в виде сплайнов 5-ой степени
Ф'к(г) =в к(г) = 6(ао + ат+2а 2т2)/Тр;
ф к(г) = в к(г) = т(6ао + 3ахт + 4а2т2);
«ФДО = шк(г) = Тр*2(3ао + ат+а2^2); (8) рк(г) = Т2ртъ(2оа0 + 5ат + 4а 2т2)/2о,
где нормированное время
т = (г-^)/Тр с [о,1]
и коэффициенты а г. вычисляются по соотношениям
ао = 2юу/Т - 8у; а1 =-8юу/Т + 5 8у;
а2 = 3юу/Т - 2 8у.
На левом конце траектории собственные производные ускорений
вк(гор ) = вк = фк(гор ) = 6 ао/Тр , а углы ф к принимают конечные значения
фк(гр ) = ф к = Тр2[12(ю у/Тр ) - 3 8у ]/2о,
к = 4,5.
Разгон до заданного значения локальной производной углового ускорения. Для
шестого элементарного движения функция
ф 6(г) должна удовлетворять краевым усло-
виям
ср6(гр) = о; ф6(гр) = о;
ф №) = о;
Ф 6(гр) = о; ф 6(гр) = о; ф,б(гр) = 8 **
(9)
она представляется сплайном 5-ой степени нормированного времени т с [о,1] в виде
8* (г) = Ь6[1 - 6т + 6т2]/Тр; в 6(г) = Ь6Т[1 - 3т + 2т2]; ш6(г) = Ь6Трт2[1 - 2т + т2]/2; (р6(г) = Ь6Тр2т3[1о - 15т + 6т 2]/6о,
86(гор ) = в° =ф 6(гор ) = 8*, угол ф 6 поворота относительно орта е6 принимает конечное
значение ф 6(гр ) = ф6 = Тр Ь6 / 6о = Тр 8* / 6о.
Позиционный переход. Функция позиционного перехода ф 3(г) по углу поворота в третьем замыкающем элементарном движении должна удовлетворять краевым условиям
ф 3^) = о; ф 3(гр) = о; ф ) = о; (11)
ф 3(гр) = ф 3; ф 3(гр) = о; ф ) = о; ф ¿у) = о. (12) Здесь угол ф 3 = ф* = 2агссо$(Л*о) определяется краевыми условиями (1) и (2), Х*о -скалярная часть кватерниона
Л *= (Хо, X *) = Л2(г; )©л,(г; )© ©Л о ©Л у ©~6(гр )©Л5(г; )©Л4(г;) с ортом оси Эйлера е3 = X* / sin(ф*/2), кватернионы Л к(гу ), к=1,2,4,5,6 однозначно определяются углами ф:к, представленными выше в явном виде, см. табл. 1, и ортами
_ у _ а .
е1 = е о ; е 2 = е о ;
е 4 = Л 6 (гр) ©Л 5 (гр) ©е> Л5 (гр) © Л6 (гр); е 5 = Л 6(гр )©е£ © ~6(гр); (13)
У* * / I * I
е 6 = еу = Еу /| Е у |,
где в соответствии с последним условием в
* в* *
(2) вектор еу = е у 8у = бу - ю у х еу.
На модуль скорости движения в позиционном переходе может накладываться ограничение с заданной константой ю* вида
тах | ф 3(т) |< ш*
те[о,1]
(1о)
где Ь 6= 8*Тр. На левом конце траектории
(14)
(14)
Рассмотрим сначала случай позиционного перехода без ограничений (14). Весь интервал такого перехода в нормированном времени те [о,1] разделим на 2 участка: первый участок длительностью т < 1 и второй участок длительностью 1-т. На каждом участке
введем свое нормированное время
т1 = (( - (0)/Т^ Т1 -ДТр
и Т2 = О - (0 - Т1 )/Т2 , Т2 - (1 - • Тр .
Угловое движение на первом участке
представим полиномом 4-ой степени ф31) (т1)
с краевыми условиями (11) на левом конце, а на втором участке - полиномом 5-ой степени ф 32)(т 2) с краевыми условиями (12) на правом конце. В момент абсолютного времени ( = (0 + дТр для гладкого сопряжения участков с достижением максимальной скорости шт = 10ф*/[Тр (4 + д)] элементарного движения такого позиционного перехода потребуем выполнение равенств
ф31)(1) = ф32)(0), Ф 31)(1) = Ф32)(0) = Ф 31)(1) = Ф32)(0)
где т1 = (( - (0)/Т2 = (( - (0 - Т1 )/Т2 ,
Вт = 6ш т/Т12; Вт = 6ш т/Т1 ; ф31 = ШтТ1/2 и
коэффициенты а0 = 2шт/Т2 ; а1 = 8шт/Т2 ;
а2 = 3шт/Т2. На левом конце траектории собственная производная ускорения
в3((0р ) - ¿0 = ф3((0р ) = 6 а0/Тр, на правом конце в3((р) - в3 = ф3((р) = 0 в соответствии с
(12), а угол ф 3 в позиционном переходе принимает конечное значение ф^) = ф3 = ф* = ®т /[(Т /2) + (2Т2/5)] = 10Шт/[Тр(4 + д)] Рассмотрим теперь вариант наличия условия (14). В этом варианте решение заключается в том, что между первым и вторым участком "вставляется" участок движения
длительностью Тс с постоянной скоростью
(15) ф3(() = ш* =const ("полка"). При этом
Параметр де (0,1) является свободнымТ его можно использовать для перераспределе-с ния интенсивности ПМ на участках: уменьшение т увеличивает интенсивность ПМ на первом участке и уменьшает ее на втором. При определении параметра т из условия равенства производной ускорения
в3(т) = ф 3(т) при сопряжении участков (т.е. из условия ф31)(1) =ф32)(0) ) получается
Д
= 42 -1 и функции позиционного перехо-
да ф 3)((), ф 3)(() принимают вид: участок 1:
в 31)(() = Ф3(1) (() = в т (1 - 2Т1); в 31)(() = ф3\() = в тТ(1 -Т1); ш31)(() = Ф3(1)(() = ®тт2(3 - 2т1); ф31)(() = ®тТТ13(2 -Т1)/2;
(16)
участок 2:
в32)(() = ф32)(() = 6(-а0 + а1т2 - 2а2т2)/Т2; в32)(() = ф32)(/) = -т2(6а0 - 3а1т2 + 4а2т2); ш32)(() = ф32)(() = шт -Т2т\(3ай -а1Т2 + а2т2);
ф32)(/) =ф31 -Т2т2(шт -Т2т22(20а0 -5а1т2 + 4а2т2)/20),
Т = дТр /(1 + 5Д); Т2 = (1 - д)Тр/[1 + 5(1 - д)]; = 5[1 - д(1 - д)(2 - д)]Тр/[1 + 5(1 + 5д(1 - д))],
(18)
* »-»
а параметр 5 = 5 , определяющий длительность Тс "полки", находится из квадратного уравнения 52 + Ь5 + с = 0, где при обозначениях z = д(1 -д)(10а - 1) и а = 10ш*Трф* коэффициенты Ь = [(10а -1) -11ад(1 - д)]^ и с = [(4 + д)а -1]^ . Если значение 5 < 0, то при шт < ш* позиционный переход будет без "полки", а при шт >ш* требуемое движение неосуществимо. Если же с[ > 0, то позиционный переход обязательно имеет "полку" длительностью Тс. Явное описание позиционного перехода на участке 1 (от начала движения до времени появления "полки") и участке 2 (с момента времени схода с "полки" до завершения движения) по-прежнему дается (16) и (17), но с подстановкой в них
значений ш* вместо шт, нормированного времени на 2 участке
т2 = (г - (го + Т + Тс ))/Т2 с [о,1]
и значения ф31 = ш* (Тс + Т1 /2). На участке
"полки" параметры движения КА определяются формулами
в3 (г) = ф3 (г) = о; в3с)(г) = ф3с)(г) = о;
ш3с)(г) = ф3с) (г) = ш*;
ф3с) (г) = Тс ш* (1/2 + те),
где нормированное время
Тс = (г - го - Т1)/Тс с [о,1].
При этом угол ф 3 имеет значение
(19)
) = = у =
а
(Т1/2) + Тс + (2Т2/5)
Юа*[1 + q(1 + q^(1 - ^))]
Задача аппроксимации программного маршрутного движения
(2о)
kp > 2. Получаемые при этом кватернион М(г) = (|до(г),д(г)),д(г) = (г)} и вектор
р (г) в явной зависимости от времени г е Тп, используемые в БЦВМ в моменты времени г, е Тп , должны быть близки к кватерниону
Л(г) и вектора ю (г) У г е Тп, соответственно. Кроме того, требуемая аппроксимация должна обеспечивать получение в явном виде векторных функций q (г) = р(г) и q (г), соответствующих функциям е (г) и е* (г). Значения кватерниона М(г" ) = М о = Л о и векторов р (г"), q (г"), ч (г") необходимы для
гладкого сопряжения МД с предыдущим участком ПМ КА, а значения кватерниона
М(г") = М„ и векторов р(г") = р„,
Ч (г" ) = чп требуются для задания начальных
условий последующего участка ПМ.
Подход к решению данной задачи осно-
Задача аппроксимации программного вывается на экстраполяции дискретно задан-
маршрутного движения (МД) КА на заданном интервале времени г е Тп = [г", г" ], где
г" = го" + Тп, состоит в расчете функций времени, приближенно определяющих кватернион Л(г) = (X о(г), Х(г)), Х(г) = (X. (г)}, векторы угловой скорости ю (г), углового ускорения Е(г) = (в. (г)} = ю (г) и производной углового ускорения Е(г) = б*(г) + ю (г) х Е(г) по значениям векторов ю , = ю (г&), заданным в диск-
ной траектории вектора программной угловой скорости КА ю к вектором р (г) У г е Тп при условиях р (гк ) = рк = ю к, k = о : (" -1), Р (гп) = р (г") = ю" с помощью " векторных
сплайнов 3 порядка рк (т) , k = о : (" -1) в нормированном времени
Т = (г-гк)/Та е[о,1], и далее аналитическом получении высокоточной аппроксимации как вектора программного углового ускорения 4 (г) = р (г) с его ло-
ретные моменты времени гб, е тп с периодом кальной производной ч (г), так и кватерниона программной ориентации М(г). При обо-
Тч = ^ - г,, ^ = о,1,2.." = о: , = Т" / Тч
и начальному значению кватерниона значениях рк (о) = рк и р'к (о) = р'к , где
Л(г") = Л о. Решение задачи для бортовой реализации должно получаться экстраполяцией значений векторов ю к = ю (гк ), определенных в моменты времени гк е Тп с шагом
Та = гк+1 - гк, к = о : " , " = Т" / Та , при кратности периодов кач = Та / Тч = 2кр степени
р' к (т) = ёрк (т)/ёт в нормированном времени т = (г - гк )/ Та е [о,1], векторный сплайн рк (т) на сегменте т = к +1, где индексы
к = о : (" -1), представляется в матричном виде
pк (т) = F(т) • Gк;
F(т) = [ВДЛ(тШт), F4(т)];
Gk = {Pк,Pк+^р'к, р'к+1^ (21) где нормированные к длине сегмента Та кубические весовые функции Эрмита р (т) = т2 (2т - 3) +1; F2 (т) = -т2 (2т - 3); Fз(т) = Тат(т-1)2; F4(т) = Тат2(т-1) (22) и использованы обозначения строки [•] и столбца {•} . При введении обозначений
Т н(т) = [1, т, т2, т3];
A н =
1 0
0 0
2 - 2 Т
0 0
Та 0
2Та - Та
Т Т
Fl(t) = 6т(т-1)/Та; F2(t) = -6т(т - 1)/ Та; F з(t) = т(3т - 4) +1;
F4(t) = т(3т- 2); (27)
2
, (23) , (23)
матрица-строка весовых функций имеет вид F(т) = Т н(т) • Aн . В нормированном времени т векторный сплайн pк (т) на т -ом сегменте аппроксимации и его производные по т представляются в явном виде
P к (т) = Ц(т) p к + Р2(т) pк+1 +
+ Fз(т) p' к + F4(т) p' к+1
Чк (Т) = P'к (Т) = Р'1 (Т) Pк + Р'2 (Т) P к+1 + + Б'3(т)р'к +Б '4 (т)р'
к+1
Ч'к (т) = F"l (т) Рк + (т) Рк+1 + + F"з (т) р'к + Р"4(т) р'к+1 , (23)
Дифференцирование весовых функций Эрмита (22) в нормированном времени т дает
F '1 (т) = 6т(т -1); F '2 (т) = -6т(т -1);
F'з(т) = Та [т(3т- 4) +1];
F '4 (т) = Тат(3т- 2), (26) На т -ом сегменте аппроксимации в абсолютном времени t е Тк - [Хк, tk+1],
tk+1 = tk + Та, к = 0 : (п -1) , учитывая что
dт / dt = 1/ Та, имеем производные этих функций
F l(t) = 6(2т -1)/ Та; F ) = -6(2т -1) / Та Fз(t) = 2(3т- 2)/Та; F4(t) = 2(3т-1)/ Та, (28)
где формально т = ^ - tk ) / Та е [0,1]. В абсолютном времени t е Тк производные векторного сплайна рк ^) имеют вид
Ч к 0) = рк 0) = ) рк + Р 2(0 рк+1 + Р з(0 р'к+Р 4(() р'к+;
(29)
Ч к 0) = рк 0) = Р1(*) рк + Р 2(0 рк+1 + Р 3(0 р'к +Р 4(*) р'к+Р
(30)
где использованы представления производных весовых функций (27) и (28). Очевидно, что векторы
р 0 = р 0(0) = ю 0 = ю(^п ) и
р п (0) = р п = р„-1(1) = ю„ = ю (^),
а векторы р'0 = р'0 (0) и р'п = р'п-1 (1) ,
представляемые в абсолютном времени t с
учетом (29) как p(ton) = р'0 и р(^) = р'п , соответствуют векторам ё 0 = ё(^) и
ёп = ё (^ ).
Аналитический расчет векторов р(^ ) и
р(^) на границах заданного интервала Тп выполняется с помощью аппроксимации вектора ю ^) на основе классической интерполяционной формулы Лагранжа порядка р е (3,4,5) на интервалах времени
t е К1, ^ + рТ5 ] и t е [^ - рТч, ^ ] по значениям векторов ю 5 = ю (1В) при ^ = 0: р и
^ = п5 - р : п5 соответственно. Такая аппроксимация в окрестности левой границы интервала Тп имеет вид
ю ^) « I ^) = 10 + 11т 5 + 12 т5 +•••+ I р тр ,
то = ('" )/Тд е [0,р],
и в окрестности правой границы этого интервала представляется как
ю (0 - г (0 = гс + Ггхч + г2^ +•••+ гртр,
то = (' + рТч - )/Т е [0, р], р = 3,4,5.
В этих соотношениях всегда векторы 10 = ю 0
и г0 = ю Пч-р, а остальные векторы вычисляются по двум однотипным блочно-матрич-ным соотношениям
(11,...1 р} = Вр • {ю 0*,...юр};
(г1,...гр} = Вр • {ю "0-р,...ю 0о}, р = 3,4,5, где постоянные прямоугольные матрицы
"-11 18 - 9 2 " 6 -15 12 - 3 -1 3 - 3 1
нозначно определяются из матричного уравнения
в а
1 6
в а
1 24
в а
1
120
- 50 96 - 72 32 - 6
35 -104 114 - 56 11
-10 36 - 48 28 - 6
1 - 4 6 - 4 1
- 274 600 - 600 400 -150 24 225 - 770 1070 - 780 305 - 50
- 85 355 - 590 490 - 205 35 15 - 70 130 -120 55 -10 -1 5 -10 10 - 5 1
"1 0
1 4 1 0
0 1 4 1 0
0 1 4 1 0
0 1 4 1
0 1
р'0 р'0
р'1 3(р2 - р0)/ Та
р'2 3(р3 - р1)/Та
р'к = 3(рк+1 - р к-1)/ Та
р'п-2 3(рп-1 - рп-3)/ Та
р'п-1 3(рп- -рп-2)/ Та
р'п _ р'п
Порядок р интерполяции выбирается из условия малости "краевых эффектов" для производных векторной функции ю (') . Учитывая dт 0 / dt = 1/ Т0, производные векторной функции р (') = ю (') на границах интервала
Тп в абсолютном времени ' рассчитываются так:
Р('") = 11/ То; р('" ) = (г + 2 р • г2 + 3р2 Г3 +•••)/Тч .(31) При векторах р к = ю к , к = 0: п и Р'о = Р('" ) , Р'" = Р('") , входящие в состав составных векторов Gk (21) векторы р'к од-
(32)
так как (п +1) х (п + 1) постоянная ленточная трехдиагональная матрица заведомо не вырождена и ее обращение выполняется специальным методом исключения Гаусса.
Значения векторов р ('" ), q ('" ) = р ('" )
и Ч ('" ) = Р ('" ), соответствующих векторам
ю 0= ю('"), 8 0 = в('") = ю('") и 80 = б*('"), которые необходимы для гладкого сопряжения МД с предыдущим участком ПМ КА, вычисляются по формулам
р ('С) = р 0 = ю 0; ч ('С) = р('С) = р\>;
Ч ) = р('С) = -6(р 0- р 1)/ Т2 - 2(2р'0 + р'1)/ Т,
(33)
а значения векторов р ('" ) = р п, ч ('" ) = ч п,
которые нужны для задания начальных условий последующего участка ПМ КА, и вектора и Ч ('" ) = Чп - по формулам
Ч ) = р('") = р'п;
р) = рп = ю ";
Ч ('") = р('") = 6(р"-1-р")/Т1 + 2(р'"-1+ 2р'")/Та.
(34)
Компактный вид векторного сплайна рк (т) на т -ом сегменте ( т = 1: п ) аппроксимации
рк (т) = п0 + т< +т2 п2 +т3 п^
к = 0: (п -1),
(35)
пк = р к;пк = Т, р'к; п2 = -3(рк - р к+1) - Т, (2р'к +р'к+1);пк = = 2(рк - рк+1) + Та (р'к +р'к+1)
(36)
следует из векторно-матричного соотношения
рк (т) = ВД • Gк = Тн • Ан • {рк , рк+1 , р'к , р'к+1
Г1 _ _ 2 _ 3 п (__ к __ к __ к __ к) = [1,т,т ,т ] • n1, П 2 , П 3 }
Прямая проверка на основе формулы (35) и следствий из нее
Чк (t) = рк (t) = (пк + 2тп2 + 3т2 п3)/Та;
Чк (t) = рк (t) = 2(п2 + 3т п3)/Та2 (37) с учетом (36), показывает справедливость соотношений (33) и (34).
При п0 ^ 0, пк ^ 0, п0 х пк ^ 0 и обок I к к\1/2 к к значениях п0 = (п0,п0) ; тк = п0 хп1 ;
тк = к, тк}12 фиксированные в базисе в
орты элементарных поворотов е= 1,2,3 и
скалярные коэффициенты а к, Ьк, ск вычисляются по явным формулам
ек = п0/п0; е3 = Шк/тк; е2 = е3 х ек ;
I к к к\ тк к к\
2 = п 2 , е2/; Ь3 = (п к , е2 ;
С2к = (п2,е3 ; сзк = п3,е3
Угловые скорости элементарных поворотов в абсолютном времени t е Тк с учетом
т = ^ - tk )/ Та е [0,1] вычисляются по соотношениям
ф к ^) = пк + ак т + ак т2 + ак т3;
ф 2(0 = т2(Ьк + Ьк т); ф 3(0 = т2(с2 + сзк т) ,(38) при этом векторный сплайн
рк (0 = п0 + тпк +т2 п2 +т3 п3 , t е Тк,
т = (t - tk)/Та, к = 0:(п -1) (39) представляется в виде
рк (t) = Г1к^) + г2к(t) + Гзк(t), t е Тк,
к = 0:(п -1), где г^ (t) = е^ фV(0, V = 1,2,3
- векторы угловых скоростей ф V(t) (38) элементарных поворотов относительно фиксированных ортов еV в базисе в, т.е. фиксированных в ССК. Аналитический расчет ква-
} терниона М^) на интервале времени t е Тп
с начальным значением М(^) = М0 = А0
может основывается на его рекуррентном представлении
М(tk) = Мк; Lk (0 = Ьк(0©Ьк(0©Ьк(0; Мк (0 = Мк ©Ьк (t) для времени t е Тк = [Хк, tk+1 ], к = 0 : (п -1) . При этом задача сводится к вычислению кватерниона Ьк (0 , соответствующего представлению сплайна (39) в виде, удовлетворяющем условиям разрешимости задачи Дарбу
р к ^) = Гзк ^) + ) © (г2к (t) +
+ ) © Г к ^) © Lk2(t)) © Ьк; (),
к = 0: (п -1),
где Ьк(0 = (13к(0, IVС*)>, IV = {1кк,12к,13к} ,
V = 1,2,3 - кватернионы элементарных поворотов. Однако с точки зрения вычислительных затрат кватернион М к ^) на каждом из п сегментов времени t е Тк = [tk, tk + Та ],
к = 0 : (п -1) наиболее просто определяется численным интегрированием по методу трапеций по явным рекуррентным соотношениям в дискретные моменты времени ts с периодом Т5, и последующей интерполяцией
значений кватерниона М к ^) в произвольные моменты времени с использованием полиномов низкого порядка.
Численные результаты
Представленные выше методы аналитического синтеза программного движения КА успешно прошли компьютерную апробацию на типовых задачах управления ориентацией КА наблюдения [5].
На рис. 2 приведены все кинематические параметры поворотного маневра корпуса КА на интервале времени Тр = [0, 85] сек с краевыми условиями
О 20 40 60 80 0 20 40 60 80
Рис. 2. Кинематические параметры поворотного маневра
Л 0 = {0.92667,-0.019725,0.37420,-0.030397} ; Л / = {0.92095,-0.092125,-0.37859,-0.0052309};
ю 0 = {-0.9,0.04,0.7}°/с;
ю / = {-0.9,-0.01,-0.7}°/с;
Б0 = {-0.01, 0,0.005}°/с2;
Б г = {-0.0119549,-0.00106716,-0.0089966}° /с,2
причем рис. 2 а соответствует варианту движения в позиционном переходе без ограничения, а рис. 2 Ь - варианту при наличии
ограничения (14) с параметром ш* = 1.5 0/с.
На рис. 3 для интервала времени Тп = [85,133]сек последующего маршрутного движения КА (длительностью 48 секунд) представлены компоненты вектора угловой
скорости ю (') и погрешности сплайновой аппроксимации при значениях периодов
Т0 = 0.25 с и Та = 2 с как угловой скорости
(компоненты вектора 5ю (') = р(') - ю (') ), так и углового положения - компоненты векторной части кватерниона
8Л(') = Л(')©М(') . Из рисунка следует, что
согласно разработанному методу в данном 177расчетном варианте погрешности аппроксимации МД гарантированно не превышают
Рис. 3. Компоненты вектора угловой скорости и погрешности сплайновой аппроксимации
1.510-7 0/сек
по угловой скорости и
210-9рад = 0''.0004 по угловому положению, причем наибольшие значения погрешностей обусловлены "краевыми эффектами" интерполяции по методу Лагранжа.
Заключение
Представлены методы аналитического синтеза основных программ углового движения КА при решении задач землеобзора, а
также результаты их компьютерной апробации. На основе таких программ по аналитическим соотношениям непосредственно на борту КА рассчитываются потребные изменения всех координат силовых гироскопов
[4].
Работа поддержана РФФИ (проект 0401-96501), Президиумом РАН (Программа фундаментальных исследований "Управление механическими системами") и Отделением энергетики, механики, машинострое-
ния и процессов управления РАН (Программа 19).
СПИСОК ЛИТЕРАТУРЫ
1. Сомов Е.И. Оптимизация экстенсивного управления при переориентации летательного аппарата с управляющими силовыми гироскопами // Оптимизация процессов в авиационной технике. Казань: КАИ. 1981.
2. Сомов Е.И., Бутырин С.А., Герасин И.А. Динамика пространственных поворотных маневров КА при цифровом управлении избыточной системой гиродинов с электромагнитным подвесом ротора // Динамика управляемых космических объектов. Иркутск: ИрВЦ СО АН СССР. 1991,
3. Somov Ye.I. Nonlinear spacecraft gyromoment attitude control // Proceedings of 1st International conference on Nonlinear Problems in Aviation and Aerospace. Daytona Beach: ERAU. 1997.
4. Сомов Е.И., Герасин И.А. Оценка реализуемости поворотного маневра космичес-
кого аппарата, управляемого избыточной системой гиродинов // Управление движением и навигация летательных аппаратов. Самара: ПФ Российской академии космонавтики им. К.Э. Циолковского. 1998.
5. Somov Ye.I., Butyrin S.A., Matrosov V.M., Anshakov G.P. et al. Ultra-precision attitude control of a large low-orbital space telescope // Control Engineering Practice. Vol. 7. №. 7. 1999
6. СомовЕ.И., Бутырин С.А., Сорокин А.В., Платонов В.Н. Управление силовыми ги-рокомплексами космических аппаратов // Труды X Санкт-Петербургской Международной конференции по интегрированным навигационным системам. С.- Петербург, ЦНИИ "Электроприбор". 2003.
7. Аншаков Г.П., Сомов Е.И., Бутырин С.А. Динамика прецизионных гиросиловых систем управления космическими аппаратами землеобзора // Труды XI Санкт-Петербургской Международной конференции по интегрированным навигационным системам. С.- Петербург, ЦНИИ "Электроприбор". 2004.
analytical synthesis of a programmed motion for remote sensing spacecraft
© 2004 Ye. I. Somov1, S.A. Butyrin2
1 Stability and Nonlinear Dynamics Research Center, IMASH RAS, Moscow 2 Research Institute for Mechanical Systems Reliability, Samara
A solution for two actual problems of analytical synthesis by spacecraft attitude programmed motion -a spatial rotation maneuver with general boundary conditions and a course motion, is presented. That solution allows to compute the required variations for all coordinates of the moment gyros by analytical relations and ensures a weak excite of the spacecraft construction element's oscillations.