Увеличение грубости к параметрической неопределенности при решении задачи фильтрации в робастных системах Текст научной статьи по специальности «Математика»

УДК 62-506

Н. В. Воробьев1, О. А. Ремизова2,

В. В. Сыроквашин3, А. Л. Фокин4

Введение

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

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

Построение динамических наблюдателей, к которым относится и фильтр Калмана, широко используется при разработке информационных систем [7], так как решение задачи оценивания не измеряемых координат объекта расширяет возможности при проектировании систем автоматизации. Например, в [8, 9] так решена проблема оценивания температурного поля протяженного термодинамического объекта с целью снижения процента брака и уменьшения количества аварийных ситуаций.

Задача оптимальной фильтрации на основе построения фильтра Калмана также оказывается востребованной при решении проблем диагностики [10]. Поэтому актуальна проблема повышения грубости систем обработки информации, построенных на основе фильтра Калмана. Для решения этой проблемы при диагностике используются минимаксные методы [11], логико-динамический подход [12, 13], интервальные [14] и другие методы.

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

УВЕЛИЧЕНИЕ ГРУБОСТИ К ПАРАМЕТРИЧЕСКОЙ НЕОПРЕДЕЛЕННОСТИ ПРИ РЕШЕНИИ ЗАДАЧИ ФИЛЬТРАЦИИ В РОБАСТНЫХ СИСТЕМАХ

Санкт-Петербургский государственный технологический институт (технический университет) 190013, Санкт-Петербург, Московский пр., д. 26

Рассмотрена методика уменьшения чувствительности процесса оптимальной фильтрации к влиянию параметрической неопределенности в математической модели объекта

Ключевые слова: неопределенность, оптимальное оценивание, состояние, фильтрация, фильтр Калмана, расширенная модель

Фильтр Калмана для расширенной модели объекта

Рассмотрим линейную модель объекта с параметрической неопределенностью вида

х0 = (А0 + АЛ)х + (В0 + АБ)ы + у0, х0(о)= х°, (1)

Уо ()= С0х0 ((), (2)

где А0, В0, Со известные номинальные значения матриц, ДА, ДВ, ДС параметрическая неопределенность в уравнении состояния, у0, -^о - некоррелированные возмущения типа белого шума в уравнении объекта и в уравнении наблюдения, х0, у0, еЛп, и еЯ", у0, w0 еЛ1.

Для начальных условий и возмущений известны вероятностные характеристики

М[хо0] = х°0 , М[(хЦ - хо0)(хо0 - Хо0)т] = Ро,

м К ()] = 0, м |у0 (Оу/ (г')] = бо (г Жг -г'X

мК (0] = 0, М[^0(0^0Т (г')] = ко-*') ,(3) где М - операция математического ожидания, 20,Р0>0,

Ло>0, я0 > 0.

Оценка Калмана Х0 для вектора состояния х0 является

несмещенной и оптимальной в смысле минимума среднего квадрата ошибки

/0 = М[(х0 ()- х0 ()У (х0 ()- х0 ())], (4)

где х0() - оптимальная оценка состояния.

Для номинальной модели она вычисляется по известным формулам [15]

хо = АоХо + В0и + К0 (0 - Сох0), *0(°)= , (5)

К о = ОСТ0 Щ1, (6)

О = А0О + ОАТ0 - ОСТЯ01С0О + бо, Ф) = Р0, (7) где о = М[е0()еТ ()]- дисперсионная матрица ошибки е0() = х0()- Х0(?) размерности п х п, Ко -, матрица наблюдателя размерности п х I.

1 Воробьев Николай Вячеславович, аспирант каф. автоматизации процессов химической промышленности, enicolas@yandex.ru

2 . Ремизова Ольга Александровна, канд. техн. наук, доцент каф. автоматизации процессов химической промышленности, remizova-oa@yandex.ru

3 Сыроквашин Владислав Викторович, канд. техн. наук, ст. преподаватель каф. автоматизации процессов химической промышленности

4 . Фокин Александр Леонидович, д-р техн. наук, профессор каф. автоматизации процессов химической, fokin_sa@mail.ru

Дата поступления - 27 июня 2011 года

Оптимальная оценка (5), вычисленная на основе Х0 М = дх () + х1 ()■ (16)

номинальной модели объекта (а0, В0, Со) имеет большую / \ /\

чувствительность к возможным вариациям этих матриц Вект°ры оценки а^() и ^ () получаются как

(ДА, ДВ, ДС). Чтобы построить наблюдатель, позволяю- результат декомпозиции вектора х() из (13). Если систе-

щий получить грубую по отношению к неопределенности оценку вектора состояния, расширим номинальную (ДА =0, ДВ = 0, ДС = 0) модель (1), (2) состояния путем разложения вектора состояния х0(г) на две составляющие при помощи системы не минимально фазовых фильтров вида

x1 (р)= Wo (p)x0 (Р) =

1 - ТФ1 р

(8)

1 + ТФ2 Р

где XI - опорная траектория, ТФ1, ТФ2 - постоянные времени фильтра.

В результате получим разложение вектора состояния х0(г) на два вектора х1(г) и вектор

Ах1 ( ) = х0 ()- х1 (), (9)

где Дх1(г) - разностный векторный сигнал.

ма (1), (2) рассматривается с параметрической неопределенностью (ДА, ДВ, ДС), то оценка (16) становится смещенной, так как нарушается взаимная компенсация векторов Дс1 (?) и X ()■

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

Л (17)

~о ( ) = Dx(t )= [ D2 ]•

Ax () xi().

где ~0 скорректированный вектор оценки для вектора

В сумме векторы х^г) и Дх^г) всегда составляют состояния хо(0 х0 ((), Б = [А °2] - матрица настраиваемых

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

(10) (11)

х = Ах + Ви + v' x(o) = (x00 )T 0

Уо = Cx + w0'

где

A = О A П • _ T-1 I T Ф21 0 A n • , B = ■ вВо ' , x = Ax1

№2 (і - — T-1T A T Ф2Т Ф1A0 _ _- ТФ2ТФ1В0 _ _ x1 _

расширенный 2n вектор состояния, р = і + Тф2Тф

C = [Co C0 ]'

в ■ Тфг^оі1.

Ковариационные матрицы для возмущения V и начальных условий х(0) расширенной модели состояния (10) будут иметь следующие значения

д

Q =

(12)

Qo [ - Т-^/]' о о

Для расширенной номинальной модели (10)' (11) можно также получить оптимальную оценку Калмана вида

г~ і \ ~~\Т

x = Ax + Bu + K(0 - Cx)' jc(o) = (x00 )T 0 '(13)

K = LC T Rq1 ' (14)

L = AL + LAT - LCTR01CL + Q' L(o) = P' (15)

где x() = [ () xf ()] - оценка расширенного 2n век-

тора состояния x(t)' Axj (t)- оценка n вектора Дхі(г) и jij ()- оценка n вектора xi(t)' L = M[e()eT ()] - дисперсионная матрица ошибки e() = x()-x() размерности 2nx2n' K - матрица наблюдателя размерности 2nx1.

Основной результат

В работах [3' 16] было показано' что основным механизмом минимизации квадратичного функционала при использовании расширенной модели объекта является взаимная покоординатная компенсация векторных составляющих x1(t) и Дх^) расширенного вектора состояния x(t).

Фильтр Калмана также является результатом минимизации квадратичного функционала (4) и при использовании формул (13) - (15) оптимальная оценка исходного вектора состояния х0 () системы (1) при номинальных

значениях матриц А0' В0' С0 получается после частичной взаимной компенсации составляющих оценки расширен-

(19)

параметров размерности п х 2п.

Физически введение матрицы Б означает преобразование Дс1 () с матрицей Б1 и преобразование ^ () с

матрицей Б2 с тем, чтобы улучшить полученную ранее оценку х0 (г) при наличии неопределенности. Для оптимальной настройки рассмотрим задачу минимизации по Б квадратичного функционала

*1 = МИуо()- соОх{ї)]■ [у()- соОх{ї)]Т}. (18) В каждый момент времени г происходит измерение вектора выхода у0(г) и вычисляется вектор оценки х() для расширенной модели объекта. Далее происходит сравнение истинного измеренного значения / х 1 вектора у0(г) и его аналога, полученного на основании оценки расширенного вектора состояния (13) по формуле у 0 () = С0 ~0 () = С0 . На основании (18) получим

М [ у о (г) ут0 (і )]-

3 = ҐГ•! -2С0Бх()М[уі (і)^ +

+С0 Бх (ґ )хт (ґ )БТСІ

Так как на матрицу Б ограничения не наложены, то продифференцируем (19) по Би приравняем нулю, получим

3 = -2 х )М [ УоТ )] С0 + . (20)

+2 х (ґ )хт (ґ )БТСТС0 = 0

Отсюда, используя операцию псевдообращения матрицы, получим

Б = (сТ Со)+ СТм()]гт ()[^(^;ЕГ ()] . (21)

После этого новая скорректированная оценка вычисляется по формуле (17). В (21) входит математическое ожидание М[у0(г)]. Оно вычисляется в реальном времени при помощи вспомогательного фильтра. В качестве такого фильтра можно использовать, например, рекуррентный алгоритм вычисления оценки математического ожидания вида

У 0

(к ):

k -1 Jo (» -1) + ІУ0 (к)' k = У-" (22>

к к' где к - дискретное время, у0 (к)- оценка математического

ожидания вектора у0, которое входит в формулу (21).

В формуле (21) присутствует матрица (с^ С0),

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

ного вектора с°ст°яния х(), что будет показано в даль- некоторые строки матрицы Б состоят из нулевых элемен-

нейшем на примере. При этом по аналогии с (9) справед- тов. Тогда дополнительная коррекция оценки вектора

лива формула состояния по формуле (17) возможна не для всех коорди-

v

0

нат вектора состояния, а только для тех, которые соответствуют не нулевым строкам матрицы О.

Это является недостатком предлагаемого метода при решении задачи оценивания. Хотя в некоторых случаях такое решение все-таки возможно и при наличии нулевых строк в матрице О. Например, в работах [8, 9] рассматривается задача оценивания температурного поля протяженного объекта при наличии измерения температуры.

Но это может быть достоинством при решении задач диагностики. Тогда на основании сравнения двух фильтров Калмана (5) - (7) и (13) - (15), (17), (21), (22) можно сделать вывод о наличии неисправности, если она вызвана изменением коэффициентов модели (1), (2). Действительно, в первом случае оценка наблюдаемых координат будет смещенной, а в другом несмещенной. Если эти оценки различаются, то это свидетельствует о наличии неисправности.

Если неисправность вызывает изменение сигнала управления, и факт такого изменения нужно установить по изменению координат вектора состояния, то здесь в условиях действия неопределенности нужно использовать только несмещенные оценки, следовательно, можно использовать не все координаты вектора состояния, а только те, для которых существует скорректированная оценка (17).

Пример

Рассмотрим линейный инерционный объект с одним входом и одним выходом, который описывается передаточной функцией

25 . (23)

(1), (г)

W° (P (30 p + l)(2G p + l)2(l0 p +1) Соответствующая номинальная модель имеет вид

¿G — AgXG + BGU + BG VG '

Уо = C0 x0 + w0 ,

(г4)

(2s)

где vo, w0 єRl - соответственно помехи объекта и наблюдения,

Ao =

- 8.33 -10-6 - 0.0007 - 0.0192 - 0.2333 0 0 0

0.2 -10-

С0 = [1 0 0 0].

Скалярное управление и в данном случае может быть любым. Здесь использовался пилообразный сигнал с амплитудой, равной 0.33, и периодом в 50с. Помехи у0, '№0 представляют собой белые шумы с интенсивностями, равными 0.5 и 0.4 соответственно. Для этой номинальной модели был построен фильтр Калмана (5) - (7) с вектором наблюдателя

К0Т = [0.0236 0.0003 -1.8-10 -5 - 4.64 • 10-7 ].

На рисунке 1 для примера штриховой линией показано изменение первой координаты х01(0 вектора состояния х0(0, а сплошной линией первая координата х01 () вектора оценки х0() вида (5) - (7).

Далее для сравнения на рисунке 2 представлены аналогичные временные зависимости для неопределенной реальной модели, у которой коэффициент передачи в передаточной функции (23) увеличен в 8 раз и теперь равен 8 • 2.5 = 20 . Хорошо видно, что оценка, сделанная по номинальной модели с коэффициентом передачи, равным 2.5, стала смещенной.

Чтобы устранить этот недостаток далее рассматривается расширенная номинальная модель (10), (11), на основе которой построен фильтр Калмана (13) - (15). Вектор наблюдателя будет

Кт = [0.148 -0.0088 -0.0003 3.4М0-5 -0.125 0.009 0.0002 -3.46-10~5]

Рисунок 1. Сравнение истинного движения и оценки Калмана для номинальной модели

Рисунок 2. Сравнение истинного движения и оценки Калма-на для номинальной модели

При этом происходит покоординатная взаимная компенсация векторных составляющих дх1 (t ) и (t ), как

это видно из рисунке З. Здесь представлены первая x (t ) (сплошная линия) и пятая x5(t) (штриховая линия) вектора оценки x(t ), которые являются первыми координатами векторов AXj (t) и x t).

Рисунок 3 Частичная взаимная компенсациялервых координат векторов Ах1 () и Х1 ().

Оценка исходного вектора состояния при этом получается по формуле (16) и она совпадает с полученной ранее оценкой, показанной на рисунке 1. Так же точно, как и раньше, в модели с неопределенностью происходит смещение оценки, как показано на рисунке 2.

Чтобы сделать оценку несмещенной используется дополнительная коррекция по формуле (17). Матрица В

G

l

G

G

вычисляется для каждого момента времени по формуле (21). Вместо математического ожидания М[уо(0] используется его оценка, которая вычисляется по формуле (22).

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

Рисунок 4. Сравнение истинного движения и оценки Калмана для реальной модели с параметрической неопределенностью задания коэффициента передачи после дополнительной коррекции.

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

Для сравнения была также рассмотрена 4x4 матрица наблюдения полного ранга Со = I. В этом случае все координаты наблюдаемы и для них оценки, полученные по формуле (17), существуют и являются несмещенными для модели с неопределенностью.

Заключение

В работе рассмотрена проблема обеспечения грубости оптимальной оценки вектора состояния для линейного динамического объекта, полученной при помощи фильтра Калмана. Эта проблема актуальна при построении информационных систем и систем диагностики.

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

Для компенсации влияния неопределенности дополнительно решается специально поставленная задача оптимальной коррекции смещенной оценки по наблюдаемым переменным модели объекта. Это позволяет получить несмещенную оценку вектора состояния.

Список литературы

1. Янушевский Р.Т. О грубости решения задачи аналитического конструирования регуляторов // Автоматика и телемеханика. 1960. № 3. C. 18-25.

2. Бахилина И.М., Степанов С.А. Синтез грубых линейных квадратичных гауссовских регуляторов // Автоматика и телемеханика. 1998. №7. С. 96-106.

3. Бороздин П.А., Сыроквашин В.В., Фокин А.Л. Робастное управление линейным инерционным объектом // Изв. РАН Теория и системы управления. 2008. № 4. С. 41-49.

4. Ремизова О.А., Рудакова И.В., Фокин А.Л. Синтез робастной системы стабилизации на основе квадратичной теории // Изв. СПбГТИ(ТУ) 2009. № 6 C. 71-75.

5. Климов А.П., Ремизова О.А., Рудакова И.В., Фокин А.Л. Уменьшение чувствительности Иг оптимальной системы к влиянию неопределенности модели объекта // Изв. РАН Теория и системы управления. 2010. № 3. С. 2732.

6. Климов А.П., Ремизова О.А., Рудакова И.В., Фокин А.Л. Достижение робастности системы стабилизации, синтезированной на основе квадратичной теории // Изв. вузов. Приборостроение. 2010. № 7. С. 18-26.

7. Игнащенко Е.Ю, Панков А.Р, Семенихин К.В. Минимаксно-статистический подход к повышению надежности обработки измерительной информации // Автоматика и телемеханика. 2010. № 2. С. 76-91.

8. Блинников A.A., Бойков В.И. Выбор параметров конечно-разностной аппроксимации модели термодинамического объекта // Изв. вузов. Приборостроение. 2007. № 11. С. 43-49.

9. Блинников А.А. Анализ свойств и синтез наблюдающего устройства для восстановления состояния протяженного термодинамического объекта. Дис. ... канд. техн. наук. СПб. 2009. 113с.

10. Кошаев Д.А. Многоальтернативный метод обнаружения и оценки нарушений на основе расширенного фильтра Калмана // Автоматика и телемеханика. 2010. № 5. С. 70-83.

11. Панков А.Р., Платонов Е.Н., Семенихин К.В. Робастная фильтрация процесса в стационарной разностной стохастической системе // Автоматика и телемеханика. 2011. № 2. С. 167-182.

12. Жирабок А.Н., Кучер Д.Н., Филаретов В.Ф. Обеспечение робастности при диагностировании нелинейных систем // Автоматика и телемеханика. 2010. № 1.

С. 159-173.

13. Жирабок А.Н. Нелинейные соотношения паритета: логико-динамический подход // Автоматика и телемеханика. 2008. № 6. С. 160-174.

14. Guanrong Chen, Jianrong Wang, Leang S. Shieh. Interval Kalman Filtering // IEEE Trans. on Aero. Elect. Sys. 1997. Vol. 33, Р. 250-259.

15. Ким Д.П. Теория автоматического управления. Т. 2. Многомерные, нелинейные, оптимальные и адаптивные системы 2- е изд. испр. и доп. М.: ФИЗМАТЛИТ, 2007. 440 с.

16. Бороздин П.А., Сыроквашин В.В., Фокин А.Л. Синтез робастной системы управления методами прямого поиска экстремума // Изв. вузов. Приборостроение. 2007. № 5. С. 25-34 .