Содержание к диссертации
Введение
1 Уровни неоптимальности для итерационных алгоритмов в методе наименьших модулей 6
1.1 Введение к первой главе 6
1.2 Метод наименьших модулей
1.2.1 Постановка задачи 7
1.2.2 Метод Вейсфельда 8
1.3 Анализ вариационных задач 9
1.3.1 Формулировки основных экстремальных задач 9
1.3.2 Двойственные задачи 10
1.3.3 Решение задач гладкой оптимизации 13
1.4 Оценки сверху для уровня неоптимальности 17
1.4.1 Подход, основанный на задаче взвешенного метода наименьших квадратов 17
1.4.2 Подход, основанный на свойствах задачи метода наименьших модулей 18
1.5 Метод наименьших модулей в спутниковой навигации 21
1.5.1 Доплеровские измерения спутниковой навигационной системы 21
1.5.2 Численное решение задачи метода наименьших модулей 23
1.6 Заключение к первой главе 26
2 Задача оценивания смещений в погрешностях БИНС 27
2.1 Введение ко второй главе 27
2.2 Математическая модель БИНС
2.2.1 Уравнения ошибок БИНС 28
2.2.2 Особенности уравнений ошибок при стендовых испытаниях 31
2.2.3 Модели погрешностей чувствительных элементов 33
2.2.4 Особенности начальной выставки 35
2.2.5 Уравнения формирующего фильтра 38
2.2.6 Динамическая система в непрерывном времени 40
2.3 Структура задачи оценивания 40
2.3.1 Дискретизация динамической системы 40
2.3.2 Структура измерений в задаче оценивания 42
2.3.3 Эксперимент с угловой и азимутальной информацией. Декомпозиция динамической системы 49
2.3.4 Моделирование массива измерений 51
2.3.5 Анализ точности приближенных моделей 55
2.3.6 Эксперимент с полной угловой информацией. Второй вариант декомпозиции системы 59
2.3.7 Масштабирование и переход к безразмерным переменным 62
2.4 Вариационные проблемы 66
2.4.1 Оценивание состояний динамической системы 66
2.4.2 Сведение к проблеме линейного программирования 69
2.4.3 Переход к проблеме безусловной минимизации 75
2.5 Заключение ко второй главе 80
3 Решение и численный анализ задач оценивания 81
3.1 Введение к третьей главе 81
3.2 Структурные элементы программной реализации 81
3.3 Оценивание скачков в инерциальных датчиках при наличии скоростной и азимутальной информации 83
3.3.1 Входные параметры вычислительной процедуры 83
3.3.2 Оценивание скачков в показаниях акселерометров 86
3.3.3 Варьирование весовых коэффициентов 92
3.3.4 Сравнение с методом наименьших квадратов 95
3.3.5 Оценивание скачков в показаниях гироскопов 100
3.4 Оценивание скачков в инерциальных датчиках при наличии полной угловой информации 105
3.4.1 Особенности вариационных проблем и их решений 105
3.4.2 Границы применимости метода іі-аппроксимации 108
3.5 Алгоритмические особенности Zi-аппроксимации 111
3.5.1 Методы, основанные на линейном программировании 111
3.5.2 Уровни неоптимальности в проблемах /і-аппроксимации для динамических систем 114
3.5.3 Оконное Zi-сглаживание для больших массивов измерений 121
3.6 Заключение к третьей главе 128
Основные результаты диссертации 129
Список литературы
- Метод Вейсфельда
- Особенности уравнений ошибок при стендовых испытаниях
- Масштабирование и переход к безразмерным переменным
- Методы, основанные на линейном программировании
Введение к работе
Актуальность темы. В задачах прикладной математики и механики часто возникает необходимость оценить значения неизвестных параметров по произведенным измерениям. В данной диссертации представлен один из подходов к задачам оценивания — метод 1\-аппроксимации, также называемый методом наименьших модулей (МНМ).
Известно, что по сравнению с многими другими методами оценивания /і-аппроксимация обладает большей устойчивостью по отношению к аномально большим ошибкам в измерениях. Кроме того, как показано в диссертации, 1\-аппроксимация может быть с успехом применена в динамических задачах оценивания в случае, когда некоторые параметры рассматриваемых систем меняются скачкообразно.
Идея использования метода наименьших модулей в прикладных задачах не нова. У истоков данного направления стояли П.С. Лаплас и Р. Боскович. Среди современных исследований в этой области необходимо упомянуть работы И.Б. Челпанова, В.И. Мудрова и В.Л. Кушко, П. Блумфилда и У. Стейгера. Важное место /і-аппроксимация занимает и в работах С. Бойда и Б.Т. Поляка, посвященных теории выпуклых задач оптимизации. Однако метод наименьших модулей непрост с вычислительной точки зрения в случае большого количества измерений и неизвестных параметров. Поэтому в динамических проблемах оценивания, нередко возникающих в навигации, /і-аппроксимация почти не использовалась. В связи с этим, значительный интерес вызывают алгоритмы численного решения, позволяющие выполнять обработку измерений за ограниченное время и с приемлемой точностью. Этому вопросу посвящена отдельная часть диссертации.
Значительное внимание в данной работе уделено применению /і-аппроксимации для обработки результатов стендовых испытаний бесплатформенных инерциальных навигационных систем (БИНС) — устройств, предназначенных для определения местоположения, скорости и ориентации движущегося объекта при помощи инерциальных датчиков (акселерометров и гироскопов). Решается задача, состоящая в определении моментов и величин скачков в показаниях инерциальных датчиков БИНС. Интерес к данному направлению исследований возник в связи с сотрудничеством лаборатории управления и навигации МГУ с Московским институтом электромеханики и автоматики.
Основы современной теории навигационных задач оценивания заложены в работах А.Ю. Ишлинского, В.Д. Андреева, Г.О. Фридлендера, Е.А. Девянина и Н.А. Парусникова. Большой практический опыт в этой области накоплен в лаборатории управления и навигации МГУ им. М.В. Ломоносова. Так, эффективные алгоритмы численного решения навигационных задач оценивания разработаны Ю.В. Болотиным, Н.Б. Вавиловой, А.А. Голованом и В.В. Тихомировым.
Часть современных исследований по навигации посвящена и обработке информации с аномально большими ошибками. Большое распространение получили робастные статистические методы, представленные, например, в монографии Дж. П. Хьюбера. Также необ-
ходимо отметить диссертацию А.Ю. Невидомского, в которой применительно к навигации и гравиметрии рассмотрены робастные методы оценивания, в том числе и метод наименьших модулей. Еще один способ оценивания в случае скачков в обрабатываемых сигналах основан на так называемом банке фильтров Калмана. Одним из первых методику использования данного подхода описал Д.Г. Лайниотис. Существенное развитие его результатов предложено в работах СП. Дмитриева, О.А. Степанова и Д.А. Кошаева, посвященных методу многоальтернативной фильтрации в навигационных задачах оценивания.
Как правило, подходы к обработке навигационной информации основаны на статистических методах оценивания. Рассмотренный в данной работе метод /і-аппроксимации принципиально отличается от них: он является детерминированным и не требует привлечения гипотез о вероятностных свойствах изучаемых явлений.
Цель работы. Цель диссертации состоит в том, чтобы на основе метода /і-аппроксима-ции разработать математическую формализацию и алгоритмы решения для навигационных проблем оценивания (в том числе, динамических), обладающие высокой эффективностью для систем с аномально большими измерительными ошибками и скачками в оцениваемых параметрах.
Достоверность и обоснованность. Разработанная в диссертации методика применения /і-аппроксиматщи основана на теории выпуклых вариационных задач, а также на принятых в математической теории навигации моделях и методах. Предложенные алгоритмы и все численные эксперименты реализованы в системе Matlab.
Научная новизна. В данной работе представлена новая методика использования 1\-аппроксимации в динамических задачах оценивания, возникающих в навигации. При помощи этого метода решена задача по определению моментов и величин скачков в показаниях чувствительных элементов бесплатформенных инерциальных навигационных систем при стендовых испытаниях.
Кроме того, для одного из алгоритмов численного решения задачи /і-аппроксимации — алгоритма Вейсфельда — предложены оценки уровней неоптимальности, позволяющие контролировать качество приближенных решений.
Теоретическая и практическая ценность. В работе на основе 1\-аппроксимации построены методы идентификации скачков в показаниях чувствительных элементов инерциальных навигационных систем. Также представлены алгоритмы обработки сигналов спутниковых навигационных систем, особенно эффективные в случае аномально больших измерительных ошибок. Построены оценки уровней неоптимальности текущих итераций для алгоритмов решения задач оценивания по методу наименьших модулей, что повышает эффективность соответствующих численных процедур.
Обсуждение работы и публикации. Результаты диссертации докладывались на конференции молодых ученых „Навигация и управление движением" (С.-Петербург, 2008 г.),
на международной научной конференции „Физика и управление" (4th International Scientific Conference on Physics and Control, Катания, Италия, 2009 г.), на 3-й мультиконференция по проблемам управления (С.-Петербург, 2010 г.), на научном семинаре им. А.Ю. Ишлинского (МГУ, 2009, 2010, 2011 гг.), на научном семинаре по теории автоматического управления и оптимизации под руководством Б.Т. Поляка (ИПУ РАН, 2011 г.).
Основные результаты диссертации опубликованы в пяти работах, еще одна работа (в трудах 18-го конгресса ИФАК) принята к публикации. Список публикаций приведен в конце автореферата. Работа над диссертацией выполнялась при поддержке РФФИ (проект № 08-08-00904-а).
Структура и объем диссертации. Работа состоит из введения, трех глав, заключения и списка литературы (57 наименований). В работе приведено 42 рисунка. Общий объем диссертации составляет 133 страницы.
Метод Вейсфельда
Пусть по произведенным измерениям необходимо оценить значения неизвестных параметров, которые связаны с измерениями соотношением: z = HTq + г, где z Є RjV — вектор измерений, г Є R — вектор ошибок измерений, q є Rn — неизвестный векторный параметр, Н — (Hi,..., HN) — заданная матрица размерности п х N {N п). Метод наименьших модулей состоит в решении вариационной задачи: /( ?) = f - tf l - jsL (1-і) г=1 т. е. в минимизации Zi-нормы вектора невязки.
Данный метод позволяет снизить влияние измерений со сбоями при получении искомой оценки параметров. В этом состоит отличие МНМ от МНК, обладающего хорошими свойствами осреднения, но не обладающего помехоустойчивостью по отношению к аномально большим ошибкам измерений.
Отметим также, что на практике различные измерения Z{ могут иметь разный "вес". Иными словами иногда необходимо, чтобы различные измерения по-разному учитывались в целевом функционале: N ад = р;к-# 1- mRfn, (1.2) г=1 Pi О - заданные весовые коэффициенты. Обычно они выбираются тисходя из представлений о масштабах измеряемых величин. Легко заметить, что проблемы (1.1) и (1.2) очень похожи: достаточно домножить zi,Hi из второй задачи на pi, чтобы считать все весовые коэффициенты равными единице. Однако с вычислительной точки зрения так не всегда удобно делать, что будет пояснено в следующем подразделе.
Существует несколько способов решения вариационной проблемы (1.1). Среди них можно указать подход, основанный на сведении данной проблемы к задаче линейного программирования [12, 47], которая может быть решена как при помощи широко известного симплекс-метода [12], так и при помощи более современного метода внутренней точки [47]. Однако решение задач большой размерности с помощью этих алгоритмов не совсем тривиально. Например, использование программ из пакета Matlab при количестве измерений N 104 проблематично, так как требует чрезмерно большого выделения памяти. Для решения задачи (1.1) используется и метод вариационно-взвешенных квадратических приближений (называемый также алгоритмом Вейсфельда), который будет описан далее. Все эти три алгоритма являются итерационными, поэтому возникает необходимость получения критериев их остановки, особенно при обработке большого количества измерений. 1.2.2 Метод Вейсфельда
Одним из алгоритхмов решения задачи (1.1) является метод вариационно-взвешенных квадратических приближений [31], идея которого была предложена Вейсфельдом в 1937 г. [55, 56]. Эта идея заключается в следующем. Вместо минимизации негладкой функции I(q) производится последовательность итераций, на каждой из которых ищется вектор, минимизирующий специальную квадратичную по q форму где к — номер итерации, W = (jzi — H7q k \) — весовой коэффициент, соответствующий г-й компоненте вектора невязки на предыдущем шаге, r/fc_1) — вектор, полученный на предыдущей итерации. В качестве начального вектора q (при к = 0), берется некоторое априорное (возможно, грубое) значение неизвестного параметра, которое требуется уточнить. Иными словами, строится последовательность векторов {q }, такая что g(fc) = argmin Q(q,q{k l)) Отметим, что задача нахождения минимума функции Q{q,q k 1 ) представляет собой проблему взвешенного МНК, для которой существуют очень простые и хорошо известные алгоритмы решения. Этим объясняется популярность алгоритма Вейсфельда в инженерной практике. Если на каком-либо шаге выполнится условие: Q{q,q ) Q(q(k\ q ) для всех q, то можно доказать, что q — точное решение проблемы (1.1) [31]. Если на каждом шаге приведенное выше условие не выполняется, то последовательность будет бесконечной.
Указанный способ получения квадратичной формы Q{q,q ) наталкивается на вычислительные сложности при малых значениях компонент вектора невязки z — HTq . А если одна из его компонент становится равной нулю, то вычислить весовой коэфи-циент W{ = (zi — Hfq№) 1 на данном шаге невозможно. Для решения этой проблемы используется прием регуляризации. Пусть а 0 — достаточно малое число, тогда г-е слагаемое в Q(q,q ) задается так [31]: 2а ГТ„\2 + Уг 0 г Ч) , если \Zi - H7qM\ а, I\ -HT%Y если г " HTq(")l а Или в терминах весовых коэффициентов: если \zi — Hfq(k )\ а, если \zi — Hfq№\ а. Тогда функция Q(q,q ) является квадратичной по q формой, весовые коэффициенты которой могут зависеть как от q k\ так и от а. При этом всегда выполняется неравенство W$k) 1/а, что обеспечивает необходимую ограниченность весовых коэффициентов.
Вопрос о сходимости алгоритма вариационно-взвешенных-квадратических приближений остается открытым, особенно в случае регуляризации. Поэтому представляется важным оценить, насколько неоптимальным является полученный на каком-либо шаге вектор q(k\ т. е. насколько I(q ) больше неизвестного оптимального значения IQ — I(q), где q — решение задачи (1.1).
Уровнем неоптимальности к—й итерации метода вариационно-взвешенных квадрати-ческих приближений назовём величину [54] где N «« = argmin J(q), J(q) = ]Г W (z{ - Hjqf. г=і Оценки неоптимальности такого типа первоначально появились в гарантирующем оценивании [9, 25, 26, 54]. Точное значение уровня неоптимальности неизвестно, так как неизвестно оптимальное значение целевого функционала IQ. В данной главе предлагаются два способа построения оценки А0 уровня неоптимальности А, один из которых может быть использован и для других итерационных алгоритмов решения задачи (1.1).
Эти оценки могут оказаться полезными при решении прикладных задач. Если на данном шаге А Д0 и А0 достаточно близко к единице, то I(q ) и IQ тоже мало отличаются. Тогда q можно взять в качестве удачного приближенного решения задачи (1.1). Таким образом, получение оценки сверху уровня неоптимальности для текущей итерации дает надежный критерий остановки вычислительного процесса.
Особенности уравнений ошибок при стендовых испытаниях
Опишем примеры использования метода Вейсфельда и построения оценок уровней неоптимальности для решения задачи оценивания параметров по большому количеству измерений (N 104). Рассмотрим задачу оценки скорости объекта по реальным доплеровским измерениям спутниковой навигационной системы GPS [И]3- Не вдаваясь в технические подробности, кратко опишем структуру доплеровских измерений GPS.
Пусть известен вектор-столбец R = (Ri,R2,Rs)T декартовых гринвичских координат приемника спутникового сигнала. Кроме того заданы гринвичские координаты М видимых в данный момент спутников (М 10): Rs&tl = (Rf1 1, Rs2at , Rlatl)T, а также их скорости l/sat = (VT\ V2sa\ VU)T, і = 1,..., M. Измеряемой величиной является рг — проекция радиальной скорости "объект — спут ник с номером г": г, = pt + Арт + Арг, г = 1,...,М, где zx — доплеровский сигнал, соответствующий спутнику с номером г, Арт — ошибка часов приемоиндикатора, установленного на подвижном объекте, Арг — погрешность измерения.
Радиальная скорость связана со скоростями спутников и искомой скоростью объекта V следующей формулой [11]: рг = hJ(VS!ltt — V), где векторы-столбцы ht задают направление "объект — спутник": (Rs - R) . Ы = и V, , = г, г = 1,.... М. Д8а «-Л2 3Автор признателен А.А. Головану и Н.Б. Вавиловой за ценные консультации по теории спутниковых навигационных систем и указание на содержательную работу [11] Таким образом, имеют место соотношения: Zi = hJ(Vs&u -V) + Арт + АРі і = 1,..., М. Исключение ошибки часов приемоиндикатора производится за счет вычитания одного равенства, например первого, из всех остальных: Zi-Z! = h\ (Vs u -V)-h{ (Vs&ti -V) + A - Apu i = 2,...,M. Следовательно, в фиксированный момент времени t скорость объекта V связана с измерениями соотношением z (t) = (H )TW + r (t), где z (t) = (hfVs - hJV )(t) - (St - 50(0, H (t) = (hi-h1)(t), r (t) = (Api-Apl)(t), i = 2,...,M(t), M(t) — количество доступных в данный момент спутников. Чтобы продемонстрировать работу описанной в данной главе методики для случая большого количества измерений, строилась следующая модельная задача. Обрабатывался массив реальных измерений системы GPS, соответствующий неподвижному объекту (скорость объекта, приемника сигнала, V = 0). При этом учитывались сигналы, поступившие на протяжении длительного интервала времени (і = 1,..., L, L 1000) в предположении, что скорость объекта оставалась постоянной. Кроме того, исходные измерения в некоторые моменты времени содержали аномально большие ошибки, что, как будет показано ниже, сильно искажает оценку, полученную при помощи МНК. Целью примера было убедиться, что предлагаемый алгоритм правильно определяет скорость, которая заранее известна, и при этом протестировать предложенный алгоритм на большом массиве измерений. , L, в один Группируя все измерения, соответствующие моментам времени t = 1, вектор-столбец, получим: ( 2 (1) (2) / г (1) \ г (2) Z = (H (l),...,H (L))TV + {L) ) { г Щ J следовательно, вектор измерений z имеет размерность iV = 2t=i(M(t) — 1), вектор неизвестных параметров (скорости объекта) — трехмерный, п = 3, а Н = (Я (1),..., H (L)) — матрица размерности 3 х N. Таким образом, задача оценивания формулируется в стандартном виде (1.1), где q = V Є R3. 1262 1261 1260 1259 1258 1257
В рассматриваемом здесь примере суммарное количество измеренрпї составило N — 10006, поэтому решение при помощи сведения вариационной проблемы к задаче линейного программирования оказалось затруднительным4. Зато метод вариационно-взвешенных квадратических приближений вполне реализуем на персональном компьютере и, будучи усовершенствован согласно описанной в данной главе методике, работает достаточно эффективно.
Для инициализации итерационного процесса в алгоритме вариационно-взвешенных квадратических приближений необходимо задать два параметра: малую величину а, используемую при регуляризации, и пороговое значение оценки уровня неоптималыю-сти Aend) по достижении которой итерационный процесс останавливается. Обе величины выбираются исходя из конкретной задачи. В данном случае они были выбраны следующими: а = Ю-10, Aend = 1 + Ю-7.
В ходе итерационного процесса считались оба варианта оценки уровня неоптималь-иости — А0 , AQ И критерий остановки имел вид: min{Ag , AQ } АопЛ. Результат работы алгоритма показан на следующих рисунках. На рис. 1 продемонстрировано изменение значения целевого функционала I(q) в ходе итерационного процесса. При этом явно наблюдается свойство алгоритма Вейсфельда — на каждом шаге значение целевого функционала монотонно уменьшается. Динамика оценок уровней иеоптимальности представлена на рис. 2, причем для наглядности показаны величины A J — 1, AQ — 1.
На последней итерации (с номером 110) оценки уровней неоптимальности были сле 4В принципе, методы линейного программирования пригодны для реализации МНМ и при больших N. Речь идет о том, что прямое их использование не всегда возможно. Иногда это требует некоторых дополнительных усилий и определенной математической культуры исследователя, например для выбора удобной формы представления задачи.
Как видно, на протяжении почти всего процесса первая оценка была точнее, ближе к единице, однако на некоторой итерации вторая оценка, строящаяся на основе знаков невязок, скачкообразно приблизилась к единице, гарантируя высокое качество решения на данной, 110-й итерации. Подчеркнем еще раз, что обе оценки носят гарантирующий характер: неоптимальность текущего решения не может быть хуже величины тіп{До , До } При этом вторая оценка на "поздних" итерациях оказывается намного точнее. Ее использование позволило раньше остановить вычислительный процесс, тем самым избежать заведомо лишних итераций. Однако на практике и первый вариант оценки оказывается вполне полезным: Д0 также достигает необходимой точности, но за большее (в полтора-два раза) количество шагов, при этом с вычислительной точки зрения данный подход намного проще.
Для сравнения приведем результаты применения метода наименьших модулей и классического метода наименьших квадратов. Так, оценка скорости, полученная с помощью алгоритма Вейсфельда для МНМ, имеет вид, хорошо согласующийся с исходной информацией о неподвижности объекта: VWeis = (-0.002, 0.0009, -0.0002)г [м/с], а МНК-оценка скорости Vls = (—1.43, 16.47, —4.22)т [м/с]. Таким образом, МНМ, в отличие от МНК, оказался мало восприимчив к аномально большим ошибкам в измерениях, а соответствующая вариационная проблема эффективно и с высокой точностью решается с помощью метода вариационно-взвешенных квадратических приближений. При этом в задачах с таким числом измерений использование линейного программирования без дополнительной модификации алгоритмов затруднительно.
Масштабирование и переход к безразмерным переменным
Итак, функции погрешностей Xj(tk) будем считать колебаниями около некоторых „средних" значений, причем иногда это „среднее" значение скачкообразно меняется. На рис. 2.1 представлен характерный вид моделируемой погрешности x\{tk) [м/с2].
Теперь рассмотрим вопрос выбора начального значения фазового вектора системы (2.30). Для этого кратко напомним об особенностях инициализации процесса вычислений при стендовых испытаниях БИНС. При инициализации задаются начальные значения координат х хН и скоростей Vx4, г = 1,2, модельной точки М , а также приближенные значения углов ориентации приборного трехгранника. Так, значение угла курса берется из показаний датчика угла стенда:
Начальные значения компонент вектора скорости в силу неподвижности основания считаются равными нулю, то есть скорости модельной и идеальной точек в начальный момент времени совпадают: Vx, = VZx = 0. В качестве координат модельной точки берутся координаты места проведения эксперимента, полученные из некоторых дополнительных источников, например, с помощью спутниковой навигационной системы [11]. Иными словами, снова используется сторонняя информация: х хЧ(0) = 2 (0), г = 1,2. Из формул (2.34)-(2.38) при таком способе инициализации следует, что Для(0) = х хЧ(0) - xx i(0) = хе(0) - ххЧ(0) = гхН, то) = Vx;i(0)-VZxl(0)=0, г =1,2. Следовательно, задав конкретные значения для соответствующих ошибок, получим и начальные значения АХІ : ГХ \ —6 М, 7V2 = 4 м, Да;1(0) = -6м, Дх2(0)=4м, 6Vi(0) = 0 м/с, 6V2(0) = 0 м/с, (2.59) где значения гхч выбраны исходя из того, что точность сторонней информации о координатах составляет несколько метров. Величины /5,-(0), j = 1,2,3, выбираются в соответствии с зависимостями (2.43): А(0) = rJ _e1_5JM а д д ft(0) _ _!Vi + l + MM, а д д &(0) = - tgv-бф; причем значения ГХЦ,ЄІ заданы ранее (см. (2.59)), широта места проведения эксперимента также считается известной, возьмем, например, ір = 56. Погрешность датчика угла стенда 5ф имеет порядок нескольких угловых минут, примем ее равной 10_3, а значения Sfi(0) также известны в рамках данного примера: 5/j(0) = а{ cos(u /0) = а{. Следовательно значения /3j(0), j = 1,2,3, вычисляются по формулам (2.43) на основе заданных (в модельной ситуации) параметров погрешностей приборов: А(0) = -8.7 1(Г5, /?2(0) = 1.8 Ю-4, /33(0) = -1.0 10"3. Таким образом, становится возможным смоделировать набор значений (Дхь Ах2, 5VU 5V2, A, fo, fo)T(tk), к = 0,..., К, соответствующий первым семи уравнениям системы (2.30), заданным начальным условиям и неоднородностям в правой части данной системы.
Продолжительность одной серии измерений при стендовых испытаниях БИНС обычно не превышает 20 — 30 минут. Интервал между двумя последовательными измерениями примем равным одной секунде, получим, что количество дискретных моментов времени в данной задаче составляет порядка тысячи, например, К = 1200. Теперь сформируем массивы измерений zxi(tk),zVi(tk), г = 1,2, z,p{tk), z0(tk),zK(tk) для модельной ситуации, для этого воспользуемся полученными ранее формулами (2.40), (2.41) и (2.48): Zxiitk) = Axi(tk)+rxi(tk), zvi(tk) = SVi(tk) + rVi(tk), г = 1,2, Z {tk) = A}(ifc)+7 (fc), zo(tk) = Pi(tk)+re(tk), zK{tk) = fo(tk)+rK(tk), k = 0,...,K.
Погрешности в координатах rxi(tk) складываются из случайной составляющей — погрешности округления r x4{tk) и постоянной ошибки сторонней информации гхц, возникающей из-за неточного знания координат объекта. Величины гхч в этом примере задаются соотношением (2.59). Положим, что г .,1(0) = г х,2(0) = 0, а гхН(Ьк), г = 1,2, к = 1,... ,К —независимые нормально распределенные случайные величины с нулевым математическим ожиданием и среднеквадратическим отклонением Ю-1 м: Щг хн(ік)} = 0, M[r xll(tk) r 2( i)] = О, (Щг х,Мг х,ШУл = Ю-1 &н (м), і = 1, 2, М = 1,
Ошибки rVi(ifc) в измерениях zvi{tk) вызваны округлением, то есть ограниченностью разрядной сетки при выводе данных на экран или в соответствующий файл. Моделировать их также будем как независимые нормально распределенные случайные величины с нулевым математическим ожиданием и среднеквадратическим отклонением 10 3 м/с2: M[rvi(tk)} = 0, M[rVi{tk)rv2(ti)} = О, (М[гу ( к)ту4( ,)])з = 1(П3 бы (м/с), г = 1,2, к,1 = 1,...,К. Погрешность в измерениях углов Д описывается соотношением (2.49) тФ(ік) = 5ф-5ф {ік)-Г- -Щ \ Lb r0(tk) = 69-59 (tk) + -, rK(tk) = SK-SK )-7 -, k = l,...,K. Lb Причем моделирование шума rxi(tk) описано ранее, а постоянные ошибки датчиков углов стенда 5ф,5в,5к примем равными 10 3. Погрешности вывода информации 5ip (tk), 59 (tfc), 6n (tk) будем считать независимыми нормально распределенными случайными величинами с нулевым математическим ожиданием и среднеквадратическим отклонением Ю-4: Щ5ф {ік)] = M[50 (tk)} = М[5к (ік)} = О, Щ5ф (ік) 5ф\и)\ = М[6в (ік) 59 ( ,)] = M[6K (tk) 6K (tt)] = 0, кфі {(М[(5ф (ік))2]У = (М[( ( ))2]) = (М[(5к (ік))2})" = 10-4, k,l = l,...,K. Таким образом, описаны характерные значения величин для моделей погрешностей, соответствующие наблюдаемым на практике уровням шумов в чувствительных элементах БИНС, а также указаны способы моделирования измерений для задачи оценивания.
После того, как рассмотрены основные аспекты, связанные с определением характерных значений параметров системы уравнений ошибок, возникает необходимость оценить, насколько корректны использованные при получении системы (2.50) допущения. Напомним, что речь шла, во-первых о переходе к дискретной системе по методу Эйлера, а во-вторых — о подстановке вместо некоторых фазовых переменных соответствующих измерений. Анализ точности при переходе к дискретной системе может быть произведен путем сравнения результатов моделирования состояний динамической системы (2.25) при различных значениях шага дискретизации и при различном количестве слагаемых в представлении экспоненты степенным рядом. Значения начальных условий и неодно-родностей в правых частях можно задать так, как это описано в предыдущем разделе.
Так, например, сравним результаты моделирования для системы (2.30) при At = 1 с и системы, полученной с помощью дискретизации уравнения (2.25) при At = 0.25 с, где матричная экспонента для формулы (2.29) вычислялась с высокой точностью при помощи пакета Matlab. При этом дополнительные шумы, описанные в предыдущем разделе не добавлялись к выходному сигналу, так как цель этого примера — оценить точность дискретизации.
Оказалось, что разница между двумя моделями очень небольшая: порядка Ю-1 м для Ахи Ю- 1 м/с для 8VX и около Ю-8 для углов Рі,р2,Рз, что заведомо приемлемо, поэтому можно считать, что дискретная модель (2.30) с высокой точностью аппроксимирует исходную динамическую систему (2.25). Если же шаг дискретизации в системе (2.30) сделать больше, например, At = 10 с, то расхождение с более точным вариантом модели (при At = 0.25 с) возрастет. Соответствующие значения Да;г будут отличаться на 1-2 метра, значения 8Ц, — приблизительно на 10_3м/с, а углы /?і,/32,Д5 — на величины порядка Ю-7. Как будет показано в третьей главе, даже использование в задаче оценивания модели (2.30) при At = 10 с позволяет с вполне приемлемой точностью выявлять скачки в показаниях чувствительных элементов БИНС.
Второй этап упрощения системы — это ее декомпозиция, описанная в разделе 2.3.3. Главная ее особенность — подстановка в правую часть системы (2.30) измерений соответствующих фазовых переменных. При этом погрешности измерений могут существенно влиять на точность, с которой получившаяся укороченная разностная схема (2.50) аппроксимирует исходную систему дифференциальных уравнений (2.25). Для того, чтобы оценить влияние погрешностей измерений, проведем сравнение наборов данных, полученных двумя разными способами: с помощью процедуры, описанной шагами а)-г) в прошлом разделе (причем с учетом погрешностей измерений), а также с помощью упрощенной системы (2.50).
Динамика изменения фазовых переменных во обоих случаях приведена на следующих графиках 2.2-2.5, сплошной линией показан результат моделирования, соответствующий методике из предыдущего раздела, а штриховой — динамика упрощенной модели (2.50).
Как видно, компоненты SVzftk) и fii(tk) очень неточно описывается упрощенной моделью. Причина этого заключается в том, что на них оказывает влияние неучтенное слагаемое «2? ( fc)- Дело в том, что ошибка r (tk) имеет достаточно высокий порядок (около Ю-3), а согласно упрощенной модели (2.50), измерение r ,(ifc), содержащее эту ошибку, подставляется вместо значения /?з(/с)- Если при формировании массива измерений не добавлять слагаемое 5ф, ошибку датчика угла стенда, то сравнительный анализ
Методы, основанные на линейном программировании
В данной главе рассматриваются вычислительные аспекты метода /і-аппроксимации в навигационных проблемах оценивания. Основное внимание уделено решению поставленной в предыдущей главе задачи определения моментов и величин скачков в показаниях инерциальных датчиков (акселерометров и гироскопов) БИНС. Эта задача решается в нескольких модификациях с привлечением различных наборов измерений.
Значительное внимание в третьей главе уделено исследованию численных особенностей получаемых решений. В частности, рассмотрен важный с прикладной точки зрения вопрос о чувствительности искомых оценок к вариациям весовых коэффициентов. Также проводится сравнение /і-аппроксимации с более распространенным методом оценивания — методом наименьших квадратов (/г-аппроксимацией).
Наконец, применительно к динамическим задачам оценивания вычисляются полученные в первой главе оценки уровней неоптимальности. Показано, что эти оценки позволяют характеризовать качество приближенных решений вариационных проблем, и, тем самым, контролировать точность вычислений.
Таким образом, главная цель третьей главы — продемонстрировать наиболее значимые свойства метода /і-аппроксимации в динамических задачах оценивания.
Рассмотрим основные структурные элементы компьютерных программ, при помощи которых производится моделирование массивов измерений и решаются проблемы оценивания.
В качестве среды разработки для указанные программ использовался пакет Matlab (версия 7.0). В принципе, соответствующие процедуры могут быть реализованы в раз личных языках программирования, например, С, C++. Главными преимуществами системы Matlab являются, во-первых, богатый выбор возможностей при операциях с матрицами, и во-вторых, большой набор встроенных в систему численных методов. В частности, следует упомянуть пакет приложений CVX, предназначенный для решения выпуклых вариационных проблем. Он также реализован в языке программирования Matlab. Пакет CVX разработан в Стэнфордском университете на базе идей, изложенных в книге [47]; он выложен в сети Internet для свободного некоммерческого использования1.
Разработанный в рамках данного исследования набор программ может быть представлен в виде нескольких модулей (процедур), изображенных на следующей схеме. Моделирование \ измерений ) Динамические модели Масштабир ование Результаты оценивания Численное решение Вариационные задачи Здесь использованы следующие обозначения. Моделирование измерений — процедура, обеспечивающая построение всех модельных данных для последующей обработки. Логика этой процедуры основана на подходе, изложенном в разделе
Динамические модели — этап, на котором вычисляются матрицы, определяющие правые части динамических моделей, и производится декомпозиция уравнений ошибок. Этому шагу соответствуют формулы (2.50)-(2.57) или (2.63)-(2.67) из предыдущей главы.
Масштабирование — простой, но принципиально важный для численной реализации шаг, на котором производится обезразмериваиие и масштабирование параметров задачи так, как это изложено в разделе 2.3.7.
Вариационные задачи — процедура, производящая построение структуры вариационных задач согласно утверждениям 2.1 и 2.2. В случае перехода к проблеме ЛП формируются матрицы и векторы для функционала (2.93) и ограничений (2.94)-(2.95). В случае перехода к задаче МНМ (2.101) формируются матрица Н и вектор Z. Если возникает необходимость сравнить результаты использования -аппроксимации (2.84),(2.85) и
На книгу [47] и пакет программ CVX автору любезно указал Б.Т. Поляк. /о-аппроксимации (2.86),(2.85), то на этом же шаге задаются матрицы и векторы задачи /о-аппроксимации (2.86),(2.85).
Численное решение — набор процедур, предназначенных для численного решения соответствующих вариационных задач. Для решения проблем ЛП используется метод внутренней точки, реализованный в пакете CVX. Заметим, что традиционным способом решения задач ЛП является симплекс-метод, однако, как показывают численные эксперименты, более современный метод внутренней точки позволяет найти решение намного быстрее, поэтому далее использоваться будет именно он. Для решения проблемы безусловной оптимизации (2.101) используется либо подход, основанный на линейном программировании, либо алгоритм Вейсфельда, реализованный автором диссертации в системе Matlab. Для алгоритма Вейсфельда при помощи отдельной процедуры вычисляются оценки уровней неоптимальности текущей итерации, согласно полученным в первой главе результатам. Поиск решения проблемы /г-аппроксимации (2.86),(2.85), как и в случае с ЛП, производится средствами встроенных в Matlab функций.
Результаты оценивания — процедура построения оценок неизвестных параметров — значений вектора состояния исходной динамической системы (2.50)-(2.57) (или (2.63)-(2.67)), на основе решений вариационных задач, полученных на шаге Численное решение. Также на данном шаге производится визуализация результатов оценивания — построение соответствующих графиков.
Таким образом, описаны основные шаги программной реализации, соответствующие изложенной в предыдущей главе методике. В следующих разделах будут продемонстрированы результаты различных численных экспериментов и обсуждены принципиально важные свойства полученных решений.
В предыдущей главе в общем виде описана вариационная задача (2.84)-(2.85), при помощи которой ищутся оценки неизвестных погрешностей БИНС. Эта проблема может быть сформулирована для частной модели (2.74)-(2.77). В качестве входных данных в этой задаче используется априорная информация о начальном значении фазового вектора гарг и набор масштабированных измерений z(k), к = 1,..., К. Также, известными являются матрицы динамической системы F,G и матрицы для измерений и априорной информации Н и Яарг соответственно.
