Метод приближенного корректирующего оператора в задачах вариационного оценивания параметров движения космических аппаратов Текст научной статьи по специальности «Математика»

УДК 629.191

В. И. Миронов, Ю. В. Миронов, Р. М. Юсупов

МЕТОД ПРИБЛИЖЕННОГО КОРРЕКТИРУЮЩЕГО ОПЕРАТОРА В ЗАДАЧАХ ВАРИАЦИОННОГО ОЦЕНИВАНИЯ ПАРАМЕТРОВ ДВИЖЕНИЯ

КОСМИЧЕСКИХ АППАРАТОВ

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

Ключевые слова: навигационное оценивание, краевая задача, метод приближенного корректирующего оператора.

Введение. В настоящее время основными при определении орбит космических аппаратов (КА) являются методы, базирующиеся на совместной обработке результатов наблюдений по полной выборке. Они позволяют успешно решать широкий круг важных и сложных прикладных задач [1—8 и др.]. В динамических задачах оценивания условий непосредственно применяются метод максимального правдоподобия (ММП) и метод наименьших квадратов. По смыслу они представляют необходимые условия оптимальности, характерные для прямых методов оптимизации.

Вместе с тем метод максимального правдоподобия может быть реализован на основе использования условий оптимальности оценок вариационного типа. Вопросы обоснования и разработки соответствующей вариационной технологии рассматривались в работах авторов [9, 10] применительно к оцениванию состояния нелинейных динамических систем.

Настоящая статья посвящена вариационному подходу к решению задач навигационного оценивания параметров движения КА по данным измерений установленной на его борту навигационной аппаратуры потребителя (НАП), работающей по сигналам спутниковой радионавигационной системы. При этом в качестве данных измерений рассматриваются оценки векторов параметров текущего состояния КА, получаемые НАП на заданном мерном интервале времени, которые подлежат совместной обработке. С целью сокращения объема вычислений при решении краевых задач вариационного оценивания применяется конструктивный метод приближенного корректирующего оператора, предложенный в работе [11].

В статье определяются и конкретизируются необходимые условия оптимальности оценок вариационного типа применительно к модели движения КА в нормальном гравитационном поле, учитывающем полярное сжатие Земли. Приводятся результаты численных расчетов.

Постановка задачи. Рассмотрим задачу оценивания параметров движения динамического объекта, которая заключается в наилучшем определении п-мерного вектора его исходного состояния Хо на начальный момент времени * = ¿о по результатам измерений, проведенных в N точках Ц , заданных на интервале измерений т = Т - ^ .

Задача. Пусть динамика объекта описывается векторным дифференциальным уравнением

Х = ф(x, ¿X Х(^о) = Хо, * Т],

измеряется да-мерный вектор

) = у[ х(*)].

Измеренное значение вектора у в момент обозначим как у (г.) = у . и представим модель измерений в виде

у (г.) = у[х(гг-)] + 8., I = 1(1) N; г. Т ].

Здесь 8. — да-мерный вектор случайных ошибок измерений, стохастическое изменение которого зададим некоторым многомерным непрерывным дифференцируемым распределением / (8., а.) с параметрами а., отличающимся в общем случае от нормального распределения.

Требуется найти такую оценку вектора х0 , которая обеспечивает минимальное значение функционала:

N

I = Е Р. {у(г.), У[х(г.)], а.},

.=1

где

Р. = 1п А {у (г.) - у[х(4), а. ]}; I = 1(1)N.

Функции ф(х, г) и у[х(г.)] будем считать однозначными, ограниченными, непрерывными и дифференцируемыми по всем аргументам во всей области их определения. Нетрудно видеть, что функционал I есть не что иное, как логарифмическая функция правдоподобия.

Для решения поставленной задачи в работе [9] были получены условия оптимальности оценок вариационного типа, которые заключаются в следующем: оптимальная оценка вектора хо и соответствующая ей оптимальная траектория являются решением краевой задачи для следующей системы дифференциальных уравнений:

х = ф(х, г); к =

дх

при граничных условиях

к (г+) = к (г-) + [у., у(х.), г1 ]; дх

к (го) = 0; к(Т) = 0; г = 1(1) N.

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

Приведенные выше условия оптимального оценивания можно конкретизировать применительно к заданному виду распределения вектора случайных ошибок измерений. Так, если для вектора 8. принимается нормальным распределение N(0, К§.) с нулевым вектором математического ожидания и корреляционной матрицей К§., что, как правило, имеет место на

практике, решение задачи оптимального оценивания сводится к решению краевой задачи

д—Т

х = —(х, г); к = ——к;

дх

к (%) = 0; к(Т) = 0;

Т

к ) = к (г-)^ {у. - у[х(гг)]}; . = 1(1) N.

Согласно этим условиям, для получения оптимальной оценки вектора Хо необходимо решить краевое уравнение

Х(Хо, Т) = 0,

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

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

Ч к+1 = М

-Е * (я?, т)

г=о

к = 1,2,3,

где М [•] — приближенный корректирующий оператор (ПКО), представляет собой алгоритм приближенного решения приближенного краевого уравнения

* (Ч, т) = о.

Конструктивной особенностью метода ПКО является возможность синтеза семейства алгоритмов численного решения краевой задачи, различающихся сложностью и скоростью сходимости в зависимости от выбора приближенной модели *(ч, Т) и принятого способа приближенного решения приведенного уравнения. Условия сходимости данного метода определены в [8] на основе принципа сжимающих отображений.

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

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

Значение оператора М [•] может изменяться или уточняться в ходе вычислительного

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

В целом метод ПКО выражает достаточно общую идеологию конструирования алгоритмов решения краевых задач и нелинейных уравнений.

Определение параметров орбиты космического аппарата по данным спутниковой навигации. Рассмотрим особенности применения вариационного метода максимального правдоподобия на примере решения задачи статистического оценивания параметров движения КА по результатам текущих навигационных оценок, полученных его бортовой навигационной аппаратурой, работающей по сигналам спутниковой навигационной системы. Движение КА будем рассматривать в нецентральном гравитационном поле Земли. Соответствующие уравнения движения в абсолютной геоцентрической системе координат имеют вид [12]

* = ^; у = vv ; г = уг ;

(

ух =

(

-—+р I *; уу =

V г )

|

— + р I у; л>г = V г3 2

Г

3 I

Р = 2"Г ^20

2 г 3

Г г2 V 1 - 5 т

V г )

Яе Г + 15 I г

~ I + тзи 40

V

-"Г + Р + Др Iг; V г3 )

2

I

Др = 3^3 Л

Г яе

20

-м +—^ л

15 |

40

Г г2 г4 V Я

-1 +14— - 21 —

2 4

V г Г )У

Г 28 г2 ^ - - - 4„

Г яе

г = \[х2 + у2 + г2:

Л20 =-1 082 627 • 10

-9 .

У40 = 23 71• 10

-9

т т

где г = [ х, у, г ] и V = [V*, Уу, у2 ] — векторы координат и скорости движения КА соответст-

3 2

венно; Яе = 6378,136 км — экваториальный радиус Земли; ц=398 600,44 км /с — постоянная притяжения Земли.

Проведем прямые полные дискретные измерения элементов векторов г и V, так что модель измерений принимает вид

у 1 = ) + 5г.; у 2 = v(ti) + 5vг.; I = 1(1) N, где 5гг-, 5vг■ — ошибки оценок НАП, полученных в моменты времени .

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

Аг =%ку-4(гггК]; АV = -аг;

г г

А = [А,г, Ау]Т = [Х*Лу,Xг,Хух,хуу,\г ]Т.

Теперь в соответствии с приведенными выше вариационными условиями оптимально-

т

сти для оценивания вектора начального состояния КА q = [г, V] необходимо решить двухточечную краевую задачу для приведенных систем уравнений движения и сопряженных уравнений с учетом следующих граничных условий:

X (¿0) = 0; А(Т) = 0;

А г (¿+) = А г ) + К -1 [у 1 - г )]; А у (¿+) = X у ) + К -V [у 2 - v(íf)]; I = 1(1) N, где К8г и К 8у — корреляционные матрицы ошибок измерений векторов координат и скорости движения КА соответственно.

Таким образом, решение задачи сводится к поиску корней краевого уравнения

А Т) = 0.

Применение метода Ньютона приводит к следующему итерационному алгоритму:

дА (^ Т)"

q к+1 = q к

дq

А ^ к, Т); к = 0,1,2,...

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

r = V;

v = -ш r; X r = ш X v;

Xv = —Xr ;

X (¿7) = X (¿7) + К -1 [X у - х(^.)]; ] = 1(1) N; X (¿о) = X (Т) = 0,

где ш — угловая скорость орбитального движения спутника в начальный момент ¿о :

Ш=Ш(Хо) = ^ого3 ; го =\1 хо + Уо + • Решение этой краевой задачи сводится к определению корня Хо нелинейного уравнения

Ni

Z V (T, tj )K—1 [Х j — U(tj, to)xo

j=i

= 0.

где

U(tj, to) =

V(T, tj) =

cos ш (tj — t0)E ш 1 sin ш (tj — t0)E —ш sin ш (tj — t0)E cos ш (tj — t0)E

cos ш(Т — tj )E ш sin ш(Т — tj )E —ш—1sin ш(Т — tj )E cos ш(Т — tj )E

E — единичная (3х3)-матрица. Нелинейность этого уравнения определяется зависимостью фундаментальных матриц U и V от параметра ш, который, в свою очередь, зависит от неизвестного вектора Х0.

Проанализировав структуру приведенных выражений, приходим к заключению, что для решения исходной краевой задачи вариационного оценивания параметров движения КА можно применить итерационную вычислительную процедуру метода ПКО с корректирующим оператором

Mk =

Ni

j=1

Z V(T, tj / X0k )K—1 U(tj, t0 / X0k)

-1

Приведем некоторые результаты вычислений. Расчеты проводились для спутника, находящегося на орбите с высотой h = 900 км и эксцентриситетом 0,003. С помощью датчика случайных величин по нормальному закону распределения на мерном интервале T = 100 с моделировалась с шагом Át = 1 с статистическая выборка прямых измерений вектора текущего состояния КА. При этом предельные ошибки измерений задавались значениями 100 м по элементам вектора координат и 1 м/с — по элементам вектора скорости. Некоторые результаты расчетов приведены в табл. 1 и 2.

Таблица 1

Результаты оптимального оценивания

Величина x0, км Уo, км z0, км vx , км/с x0 vy , км/с y0 v , км/с z0 '

Точное значение 0 -7249,926 0 0,904 95 0,005 75 7,370 23

Начальное приближение 50 -7199,926 50 0,954 95 0,055 75 7,820 23

Оптимальная оценка -0,006 -7249,928 -0,007 0,904 97 0,005 80 7,370 30

Ошибка оценивания -0,006 -0,002 -0,006 -0,000 01 0,000 05 0,000 06

Таблица 2

Сходимость вычислительного процесса

Номер итерации Ошибки оценивания

5 x0, м 5 Уo, м 5 z0, м 5 Vx , м/с 5 Vy , м/с у0 5 К , м/с

0 50 000 50 000 50 000 50 50 50

1 -7,6 -2,8 -7,9 0,05 0,04 -0,01

2 -6,3 -1,9 -6,5 0,02 0,05 0,06

3 -6,3 -1,9 -6,5 0,02 0,05 0,06

В табл. 1 даны точные значения параметров начального фазового состояния КА в абсолютной геоцентрической системе координат, принятое начальное приближение элементов уточняемого вектора, полученные в результате вариационной обработки измерений оптимальные оценки, а также характеристики точности оценивания. В табл. 2 представлены значения ошибок оценивания по итерациям. Приведенные в таблицах данные свидетельствуют о высокой точности и скорости сходимости вычислительного процесса. Расчеты показали, что применение метода ПКО в принятых условиях позволяет сократить объем вычислений примерно в четыре-пять раз.

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

Работа выполнена при поддержке РФФИ (проект № 09-08-00259).

СПИСОК ЛИТЕРАТУРЫ

1. Аким Э. Л., Энеев Т. М. Определение параметров движения космических аппаратов по данным траекторных измерений // Космические исследования. 1963. Т. 1, № 1. С. 5—50.

2. Брандин Н. К., Разоренов Г. Н. Определение траекторий КА. М.: Машиностроение, 1978. 216 с.

3. Космические траекторные измерения / Под ред. П. А. Агаджанова, В. Е. Дулевича, А. А. Коростелева. М.: Сов. радио, 1969. 504 с.

4. Линник Ю. В. Метод наименьших квадратов и основы теории обработки наблюдений. М.: Физматгиз, 1958. 350 с.

5. Основы теории полета космических аппаратов / Под ред. Г. С. Нариманова и М. К. Тихонравова. М.: Машиностроение, 1972. 608 с.

6. Шебшаевич В. С., Дмитриев П. П., Иванцевич Н. В. и др. Сетевые спутниковые радионавигационные системы. М.: Радио и связь, 1993. 408 с.

7. Статистические методы обработки результатов наблюдений / Под ред. Р. М. Юсупова. М.: МО СССР, 1984. 563 с.

8. ЭльясбергП. Е. Определение движения по результатам измерений. М.: Наука, 1976. 416 с.

9. Миронов В. И., Миронов Ю. В., Юсупов Р. М. Вариационное оценивание состояния нелинейной динамической системы по критерию максимального правдоподобия // Мехатроника, автоматизация, управление. 2009. № 11. С. 2—6.

10. Миронов В. И., Миронов Ю. В., Юсупов Р. М. Метод наименьших квадратов в задачах вариационного оценивания состояния нелинейных динамических систем // Информационно-управляющие системы. 2009. № 6. С. 2—6.

11. Миронов В. И., Миронов Ю. В., Юсупов Р. М. Синтез итерационных алгоритмов решения краевых задач и нелинейных уравнений // Изв. вузов. Приборостроение. 2010. Т. 53, № 1. С. 9—14.

12. Мамон П. А., Половников В. И., Слезкинский С. К. Баллистическое обеспечение космических полетов. Л.: ВИКИ им. А. Ф. Можайского, 1990. 622 с.

Сведения об авторах

Вячеслав Иванович Миронов — д-р техн. наук, профессор; Санкт-Петербургский институт информа-

тики и автоматизации РАН; E-mail: [email protected] Юрий Вячеславович Миронов — д-р техн. наук, доцент; Санкт-Петербургский институт информатики и

автоматизации РАН; E-mail: [email protected] Рафаэль Мидхатович Юсупов — д-р техн. наук, профессор; Санкт-Петербургский институт информатики и автоматизации РАН; директор; E-mail: [email protected]

Рекомендована институтом Поступила в редакцию

07.07.10 г.

УДК [517.938 + 519.713/718]: 621.398

А. В. Ушаков, Е. С. Яицкая

РЕКУРРЕНТНОЕ СИСТЕМАТИЧЕСКОЕ ПОМЕХОЗАЩИТНОЕ ПРЕОБРАЗОВАНИЕ КОДОВ: ВОЗМОЖНОСТИ АППАРАТА ЛИНЕЙНЫХ ДВОИЧНЫХ ДИНАМИЧЕСКИХ СИСТЕМ

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

Ключевые слова: двоичная динамическая система, систематический помехо-защищенный код, проверочная и образующая матрицы, образующий модулярный многочлен, помехозащитное преобразование кодов.

В работах [1, 2], посвященных преобразованию кодов (ПК), представленных в виде кодовых последовательностей, получен конструктивный инструментарий модельного представления процедур ПК в виде линейных последовательностных машин (ЛПМ). Так как термин ЛПМ выпадает из общей теории систем, то со временем он был заменен понятием „линейная двоичная динамическая система " (ЛДДС). Таким образом, ЛДДС, использующая в основном векторно-матричные модельные представления (ВММП) аппарата пространства состояний (АПС), со временем стала инструментом ПК, в том числе и помехозащитного преобразования систематических кодов в задачах кодирования и декодирования [3—5]. Из теории АПС [6] известно, что одной из проблем формирования векторно-матричного модельного представления динамических процессов над бесконечными и конечными полями является поиск базиса представления, в котором матричные компоненты ВММП обладают желаемыми исследователю свойствами. Причем возможность выбора базиса представления при построении ВММП процессов общего вида практически не ограничена. К сожалению, этого нельзя сказать о ВММП ЛДДС, использованных в задачах помехозащитного преобразования кодов.

Известно [7, 8], что систематическое помехозащитное преобразование кодов (ППК) в задачах кодирования и декодирования может быть осуществлено несколькими способами, основанными на различном описании процесса ППК. Тем не менее эти способы имеют единое