Метод снижения вычислительной сложности ансцентного фильтра Калмана Текст научной статьи по специальности «Математика»

УДК 519.6

Р.В. Белов1, Д.А. Кляпнев1, К.О. Огородников1, 2

МЕТОД СНИЖЕНИЯ ВЫЧИСЛИТЕЛЬНОЙ СЛОЖНОСТИ АНСЦЕНТНОГО ФИЛЬТРА КАЛМАНА

ПАО «АНПП «Темп-Авиа», г. Арзамас1,

Нижегородский государственный технический университет им. Р.Е. Алексеева2

Рассматривается возможность уменьшения количества математических операций, осуществляемых при проведении ансцентного преобразования, с помощью перестановки элементов вектора состояния фильтра Кал-мана. При использовании ансцентного преобразования вычисление множества сигма-точек производится с применением разложения Холецкого P = SST. Поскольку матрица S треугольная, то в наборе сигма-точек могут присутствовать одинаковые элементы. Это обстоятельство может быть использовано для уменьшения вычислительной сложности ансцентного преобразования.

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

Ключевые слова: ансцентное преобразование, вычислительная сложность, фильтр Калмана.

Введение

В настоящее время для решения навигационных задач часто используется фильтр Калмана (ФК) в его традиционной интерпретации, позволяющий проводить оценку вектора состояния динамической системы. Однако во многих случаях уравнения состояния системы и модели наблюдения оказываются нелинейными, и поэтому при решении задач оценивания практикуется применение ансцентного ФК (UKF - Unscented Kalman filter), позволяющего учитывать нелинейный характер задачи [1-8]. Кроме того, UKF не требует вычисления Якобиана, который в некоторых случаях оказывается затрудненным.

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

Вычислительная сложность ансцентного преобразования в значительной степени определяется количеством сигма-точек. Известно, что минимальный набор состоит из n+1 сигма-точек [1]. Вместе с тем, дополнительное уменьшение количества математических операций, осуществляемых при проведении ансцентного преобразования, оказывается возможным за счет перестановки элементов вектора состояния ФК. В работе показано, что объем выполняемых вычислительных операций зависит от того, в каком порядке составлен вектор состояния и, соответственно, от записи уравнений преобразования.

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

© Белов Р.В., Кляпнев Д.А., Огородников К.О., 2019.

Ансцентное преобразование

Предположим, что задано нелинейное преобразование:

У = / (х), (1)

где х - случайный вектор размерности п с математическим ожиданием ~ и матрицей ковари-аций Рх; / (•) - известная вектор-функция размерности I.

Требуется рассчитать оценки ковариации Ру и математического ожидания ~ вектора у . В соответствии с процедурой ансцентного преобразования эти моменты вычисляются как

~ = 1ЪКХ,), Ру = п £ (/(х,) - у)/(х, ) - ~У, (2)

т 1=1 т 1=\

где X, - сигма-точки, с известными свойствами [1-8]; т - их общее количество.

Вычисление множества сигма-точек производится с использованием факторизации Рх = SST, где S является треугольной. Следует учесть, что возможно несколько вариаций построения набора сигма-точек:

- т = 2п. Тогда они вычисляются следующим образом:

[х21 - = ~ + ^

< ' ~ /, 1 = 1,..., п, (3)

I х21 = ~ - ^

где х21_1з х2, - две симметричные сигма-точки, соответствующие 1-му столбцу матрицы S; 8,

- 1-ый столбец матрицы 8;

- т = 2п+1. Тогда к предыдущему набору добавляется центральная точка:

х2„+1 = ~; (4)

- т = п+1. Для нахождения сигма-точек строится симплекс размерностью п. Процедура получения сигма-точек описана в [1, 4].

Идея снижения вычислительной сложности

Поскольку матрица 8 треугольная, то в наборе сигма-точек могут присутствовать одинаковые элементы. К примеру, если вектор х двумерный и центрированный, то при т = 2п = 4 в результирующем наборе будут присутствовать две точки, равные нулю. Это обстоятельство может быть использовано для уменьшения вычислительной сложности ансцентного преобразования. Покажем это на простейшем примере для случая т = 2п. Рассмотрим нелинейное преобразование:

/ х2 ) = х2 ) + g2(x2), (5) где х1, х2 - случайные величины, для простоты центрированные; /(•), (•), з2 (•) - скалярные функции.

Если применить разложение Рх = ££т для вектора [х1з х2 ]т , то при вычислении моментов (2) значения функций з 1, з 2 требуется рассчитать во всех сигма-точках. Если же произвести перестановку элементов х1 и х2, то для полученного вектора [х2, хх ]т значение функции З 2 нужно дважды вычислить в точке 0. Таким образом, в первом случае нужно 4 раза рассчитать з 2, а во втором достаточно трех вычислений.

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

Q = ^, (6)

где ех = (1 1 1 — 1) - вектор-строка размерности I; W - весовая матрица размерности I • п , выражающая зависимость /-го уравнения преобразования от j-го параметра и определяющая сложность повторного вычисления /-го уравнения при изменении j-го элемента вектора состояния; £ - вектор-строка размерности п, ^ый элемент которого определяет количество различных ^ых элементов в наборе сигма-точек.

Вектор-строка £ принимает значения в зависимости от выбранной вариации построения набора сигма-точек:

%2п =[1 11 ■■■ 11 0] + 2[1 2 3 — п-2 п-1 п] для т = 2п, (7) £2и+1 = 1 + 2[2 3 4 ■■■ п - 2 п -1 п] для т = 2п+1, (8)

£и+1 =[1 11 — 11 о]+[1 2 3 ■■■ п - 2 п -1 п] для т = п+1. (9)

Уравнение (6) позволяет достаточно просто построить процедуру синтеза оптимальной перестановки элементов вектора х и, как следствие, строк матрицы W, обеспечивающей минимизацию (6). К примеру, считая каждую из функций g 1 и g 2 в (5) одинаково сложной, составим весовую матрицу W . От параметра х1 зависит только функция g1, следовательно, при изменении х1 требуется повторно вычислить только g 1, поэтому назначается вес, равный 1. От параметра х 2 зависят как g 1, так и g 2, поэтому назначается вес, равный 2. В итоге получаем W = [1 2] и, следовательно, Q = [1][1 2][3 4]т = 11. Если же элементы вектора поменять местами и произвести соответствующую перестановку для W = [2 1], то в таком случае Q = 10.

Следует отметить, что описанный подход применим не всегда. Если в (5) функция g 2 будет зависеть не только от х2, но и от х1, то перестановка элементов вектора [ х1, х2 ] не позволит уменьшить вычислительную сложность ансцентного преобразования ввиду того, что значение функции g2(х1, х2) будет отличаться при переходе к каждой последующей сигма-точке. Записав W = [2 2] для этого примера, можно видеть, что перестановка элементов не приведет к изменению оценки сложности. Кроме этого, необходимо иметь в виду, что в рассмотренном примере функции g 1 и g 2 считались одинаково сложными. В общем случае при формировании матрицы W необходимо учитывать трудоемкость вычисления нелинейной функции / .

Минимизацию вычислительной сложности Q ансцентного преобразования обеспечивает перестановка, полученная путем сортировки вектора elW по убыванию, т. к. вектор £ представляет собой неубывающую последовательность.

Пример синтеза перестановки уравнений

Рассмотрим алгоритм беспоисковой корреляционно-экстремальной навигационной системы (КЭНС) по рельефу местности, способной проводить коррекцию решения инерциаль-ной навигационной системы (ИНС) посредством сопоставления наблюдаемого поля высот с заранее подготовленным эталоном. В данном случае беспоисковая КЭНС основана на анс-центной калмановской фильтрации, позволяющей оценивать координаты и вектор скорости путем статистической линеаризации информативного поля.

Уравнения объекта описывают линейную зависимость между параметрами:

д' 1 т т! 0 0 0 0 " д

ду 0 1 2 т 0 0 0 0 ду

да 0 0 1 0 0 0 0 да

дух = 0 0 0 1 0 0 0 дух

дУ2 0 0 0 0 1 0 0 дуг

дх 0 0 0 т 0 1 0 дх

д к 0 0 0 0 т 0 1 д2

(10)

J к-1

где I - период дискретизации; дк,ду,да - ошибки по высоте, вертикальной скорости и ускорению соответственно; 5vx,5vz - ошибки по горизонтальным проекциям вектора скорости; дх, д - ошибки по горизонтальным координатам.

Модель измерений выглядит следующим образом:

= Кнс -кпп[хиНС +дx,zиHC + &]+&, (П)

где кинс - высота, полученная от ИНС; хинс, гинс - координаты, полученные от ИНС; крв - высота радиовысотомера; Ъпп [хинс + дх, гинс + &] - функция выборки высоты подстилающей поверхности из эталонного массива высот по координатам хинс + дх, гинс + д.

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

Определив зависимость одних параметров от других, запишем W, учитывая сложность проводимых вычислений:

W = [1 00002 2] (12)

Заметим, что наибольший вес имеют ошибки по координатам, что связано с высокой сложностью выборки высоты подстилающей поверхности.

Запишем вектор е1№ , который в данном случае идентичен W:

е1Ж = [1 0 0 0 0 2 2] (13)

Осуществив сортировку вектора e1W по убыванию, получим Wp с учетом перестановки, вектор е1Жр :

=[2 2 1 0 0 0 0],

(14)

и, соответственно, вектор перестановки р:

р = [6 7 1 2 3 4 5] (15)

Для подтверждения найденного решения был реализован ансцентный ФК, а ансцентное преобразование проводилось по 2п сигма-точкам. Результатом проведения перестановки стало перемещение элементов вектора состояния, от которых зависит вычисление указанных выражений, в начало вектора состояния.

ду да дух дvz дх дг]т [дх д д ду да дух дуг] (16)

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

Аналогичные действия были проведены для набора, состоящего из п+1 сигма-точек. Здесь вычисления необходимо проводить для 4 из 8 и 5 из 8 сигма-точек соответственно.

Результаты, полученные путем математического моделирования, подтверждают корректность разработанного метода перестановок элементов вектора состояния для увеличения эффективности ЦКР (рис. 1).

100 90 80 70 60 50 40 30 20 10 0

1

2п

2 3 4

2п с перестановкой п+1 п+1 с перестановкой

Рис. 1. Время выполнения ансцентного преобразования

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

Проведем вычисление оценки сложности преобразования сигма-точек без учета и с учетом полученной перестановки:

0 = е,^ = 57; (17)

б = е1Шре2п = 23. (18)

Соответственно, время выполнения ансцентного преобразования должно уменьшиться на 59,6%:

О - О 57 - 23

^^ • 100% = ' б 57

Сравнение оценок времени при вычислении ансцентного преобразования до и после проведения перестановки элементов вектора состояния показало, что уменьшение времени выполнения ансцентного преобразования составило 55,7%:

(58мс -16.9мс) - (35. Ыс -16.9мс) _ шо% ^ 55 ?% _ (20)

(58мс -16.9мс)

Оценка сложности преобразования также была вычислена и для набора, состоящего из п+1 сигма-точек. Соответственно, уменьшение времени выполнения ансцентного преобразования должно составить 53,3%:

е = еЩТп+1 = 30; (21)

(22)

• 100% = ——— • 100% « 59.6%.

(19)

<2 = е^^ = 14;

О - О 30 -14

р • 100% = —--100% « 53.3%.

б 30

В этом случае сравнение оценок времени показало, что время выполнения ансцентного преобразования уменьшилось на 51,3%:

(33.6мс -10.6мс) - (21.8мс -10.6мс) 100^ 51 3% (24)

(33.6мс -10.6мс)

Оценки сложности ансцентного преобразования, полученные для наборов, состоящих из 2п и п+1 сигма-точек, соответствуют реальным полученным результатам и, таким образом, позволяют получить предварительное заключение о возможности увеличения эффективности ансцентного ФК.

Заключение

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

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

Библиографический список

1. Julier, S.J. Reduced sigma point filters for the propagation of means and covariations through nonlinear transformations / S.J. Julier, J.K. Uhlmann // IEEE American Control Conference. - 2002. - Vol. 2. -P. 887-892.

2. Julier, S.J. Unscented Filtering and Nonlinear Estimation / S.J. Julier, J.K. Uhlmann // Proc. IEEE. -2004. - Vol. 92(3). - P. 401-422.

3. Kandepu, R. Applying the unscented Kalman filter for nonlinear state estimation / R. Kandepu, B. Foss, L. Imsland // Journal of Process control. - 2008. - Vol. 18, № 7-8. - P. 753-768.

4. Moireau, P. Reduced-order unscented Kalman filtering with application to parameter identification in large-dimensional systems / P. Moireau, D. Chapelle // COCV 17. - 2011. - P. 380-405.

5. Simon, D. Optimal State Estimation. Kalman, Нда, and Nonlinear approaches / D. Simon. - Hoboken, 2006. - 550 p.

6. Wan, E.A. The square-root unscented Kalman filter for state and parameter estimation / E.A. Wan, R. Van der Merwe // Acoustics, Speech, and Signal Processing. Proceedings of the IEEE International Conference. - 2001. - Vol. 6. - P. 3461-3464.

7. Wan, E.A. The unscented Kalman filter for nonlinear estimation / E.A. Wan, R. Van der Merwe // Symposium on Adaptive system for signal processing, communication and control (AS-SPCC). - 2000. - Vol. 6. - P. 153-158.

8. Степанов, О.А. Линейные оптимальные алгоритмы в задачах оценивания с нелинейными измерениями. Связь с алгоритмами калмановского типа / О.А. Степанов, А.Б. Торопов // Известия ТулГУ. Технический вестник. - 2012. - № 7. - С. 172-189.

9. Markley, F.L. Fundamentals of Spacecraft Attitude Determination and Control / F.L. Markley, J.L. Cras-sidis. - New York: Springer, 2014. - 486 p.

Дата поступления в редакцию:11.01.2019

R.V. Belov1, D.A. Klyapnev1, K.O. Ogorodnikov1, 2

COMPUTATIONAL COMPLEXITY REDUCING METHOD IN APPLICATION

OF UNSCENTED KALMAN FILTER

PJSC «Arzamas research and production enterprise «Temp-Avia»1, Nizhny Novgorod state technical university n.a. R.E. Alekseev2

Purpose: This paper considers the possibility of reducing the number of mathematical operations performed in unscented transform that is done by permuting the elements of the Kalman filter state vector.

Approach: Using unscented transform we calculate a set of sigma-points with Cholesky factorization P = SST. Equal elements in set of sigma-points could be found since the matrix S is triangular. This circumstance can be used to reduce the computational complexity of the unscented transform.

Results: The considered method allows us to calculate the most complex expressions fewer times by carrying out the state vector elements permutation. It is applicable to all of three known variants of sigma-points construction mentioned in this article. However, the possibility of computational complexity reducing depends on the sequence of elements in state vector, and, correspondingly, on the equations themselves.

Key words: unscented transform, computational complexity, Kalman filter.