Вычислительные методы алгебры и оценивания. И.В. Семушкин (2011) (1185350), страница 52
Текст из файла (страница 52)
Этап экстраполяции: совпадает с этапом экстраполяции базовогометода, т. е. РККИФ см. формулу (14.19).32814.9 Скаляризованный модифицированный квадратно-корневой И-фильтр14.9Скаляризованный модифицированныйквадратно-корневой информационный фильтрAлгоритм СМККИФI. Этап обработки измерений (фильтрация):••(0)−T /2T(0), xt = x̂−Положить: St = P̃tt .Для k = 0, 1, .
. . , m − 1 вычислять:0(k+1)St(k+1)K̄p,t−T(k+1)St−(k) = Ot,1 (k+1)htT(k+1)−ēt(k+1) (k+1)xtΦTt St(k+1)/σt(k)St=0−T(k)StΦTt(k+1)(k+1)/σt−zt(k) (k)St x t,(14.28)(k)где Ot,1 — ортогональное преобразование, которое приводит либок нижнему треугольному виду первый (блочный) столбец, либо кверхнему треугольному виду второй (блочный) столбец матрицы,стоящей в правой части формулы (14.28).•−T /2Положить: P̂t(m)= St−T /2 +x̂tи P̂t(m) (m)= St x t .II. Этап экстраполяции: совпадает с этапом экстраполяции базовогометода, т.
е. МККИФ, см. формулу (14.22).14.10Скаляризованный комбинированныйквадратно-корневой фильтрAлгоритм СКоККФI. Этап обработки измерений (фильтрация):•(0)Положить: St1/2(0)= P̃t , xt = x̂−t .32914 Ортогонализованные блочные алгоритмы•Для k = 0, 1, . . . , m − 1 вычислять:T(k+1)(k+1)(k+1)0−ētK̄p,t re,t−T −T=(k+1) T(k+1)(k+1)(k+1)Φt0StStStxtT(k+1)(k+1)(k+1)(k+1)(k+1)/σt− ht−zt/σt0σt(k) −T−T= Ot,1 ,(k) T(k) (k+1)(k)(k)(k)St Φ tSt htStStxt(14.29)(k)•где Ot,1 — ортогональное преобразование, приводящее либо к верхнему треугольному виду первых два (блочных) столбца, либо книжнему треугольному виду третий (блочный) столбец матрицы, стоящей в правой части формулы (14.29).(m)(m)1/2Положить: P̂t = St и x̂+t = xt .II.
Этап экстраполяции: совпадает с этапом экстраполяции базовогометода, т. е. КоККФ, см. формулу (14.25).Теоретические обоснования четырех алгоритмов, приведенных выше вподразд. 14.7, 14.8, 14.9 и 14.10, помещены в приложение A.1 (стр. 339) вавторской редакции М. В. Куликовой [51].14.11Задание на лабораторный проект № 10ВведениеВ данном лабораторном проекте мы предлагаем к изучению изложенныевыше ортогонализованные блочные алгоритмы калмановской фильтрациина материале реальной прикладной задачи. В качестве такой задачи возьмем задачу получения оценок высоты и вертикальной скорости летательногоаппарата (ЛА) по показаниям двух приборов: инерциального датчика вертикального ускорения ay (t) и барометрического датчика высоты, которыйобладает собственной инерционностью в своих показаниях h(t), характеризуемой постоянной времени τ .Введем обозначения физических величин, участвующих в математической модели движения объекта (ЛА) по высоте и в модели наблюдения за33014.11 Задание на лабораторный проект № 10этим движением.
В соответствии с принципиальной моделью — вторым законом Ньютона, движение центра масс объекта характеризуем тремя величинами:1) y = y(t) — вертикальная координата (высота над Землей),2) vy = vy (t) — вертикальная скорость,3) ay = ay (t) — вертикальное ускорение.Будем рассматривать случай движения с постоянной силой тяги двигателяЛА.
В этом случае следует считать, что ускорение ЛА тоже постоянно (пренебрегая изменением массы ЛА), однако, оно неизвестно, то есть подлежит оцениванию по показаниям приборов. В роли оценивателя должен бытьиспользован фильтр Калмана, который, как известно, оптимален, если обемодели — модель состояния и модель наблюдения — линейные.Прочитанные показания приборов (измеренные величины) обозначим:1) za = za (t) — измеренное ускорение a(t),2) zh = zh (t) — измеренная барометрическая высота h(t).Считывание этих показаний происходит в дискретные моменты времени tiс аддитивными погрешностями v1 (ti) и v2(ti ), соответственно.
Это означает,чтоza (ti) = a(ti ) + v1(ti ) , zh (ti) = h(ti ) + v2(ti ) .Погрешности датчиков между собой независимы и в отдельные моментывремени ti являются гауссовыми случайными величинами с нулевыми средними значениями и с постоянными дисперсиями σ12, σ22 , соответственно.Истинная высота y(t), на которой находится ЛА, воздействует на барометрический датчик (БД) высоты так, что барометрическая высота h(t) подчиняется уравнениюdτ h(t) + h(t) = y(t) .dtТеперь есть все сведения, чтобы записать уравнения состояния обобщенного динамического объекта, включающего ЛА и БД.
Для этого введемосновные обозначения: 1) вектор состоянияTx = x(t) = x1 ≡ y(t) x2 ≡ vy (t) x3 ≡ ay (t) x4 ≡ h(t)и 2) вектор измеренийTz = z(t) = z1 ≡ za (t) z2 ≡ zh (t)33114 Ортогонализованные блочные алгоритмыTc погрешностью v(t) = v1(t) v2(t) .чаем непрерывную модель состояния0 1 0 0dx(t) = F x(t) , F = 0 0dtα 0и дискретную модель наблюденияz(ti ) = Hx(ti) + v(ti) ,Из принятых предположений полу0 01 0 ,0 0 0 −αH=α , 1/τ0 0 1 00 0 0 1.(14.30)(14.31)1. ЗаданиеA. Построить дискретную модель обобщенного объекта (14.30), (14.31).Для этого найти резольвентную матрицу Φ(s) , (Is − F )−1, где s — комплексная переменная преобразования Лапласа [75]. Доказать, чтоs−1s−2s−30s−1s−200.Φ(s) = −100s0s−1(s + α)−1α s−2 (s + α)−1α s−3 (s + α)−1α (s + α)−1Совершить отсюда обратное преобразование Лапласа [75] и тем самым найтипереходную матрицу состояния на интервале времени t:a1 (t) = 1 − e−αt ,1tt2 /20 0a2 (t) = (αt − 1 + e−αt )/α ,1t0 2,Φ(t) = −αt010 0a3 (t) = (1 − αt + (αt))/α2 ,2 −ea1 (t) a2 (t) a3 (t) a4 (t)a4 (t) = e−αt .Задать τs — постоянный интервал выборки (sampling interval), т.
е. темпсчитывания данных с датчиков; при этом исходить из условия γ , τs /τ < 1,например, γ = 1/10 или меньше. Подстановкой t = τs в Φ(t) определитьпостоянную Φ , Φ(τs) — переходную матрицу состояния дискретной моделиМодель А:xt+1 = Φxt ,zt = Hxt + vt ,(14.32)которая является частным случаем общей модели (14.1), (14.2), где нижнийиндекс в записи {·}t есть индекс (номер) дискретного момента времени для33214.11 Задание на лабораторный проект № 10переменной {·}. Таким образом, исходная дискретная модель общего вида(14.1), (14.2) конкретизирована для исследуемого объекта в виде (14.32).Б. Допуская некоторое усложнение модели (14.32), ввести в нее «шумовую» составляющую Gt wt по формуле()xt+1 = Φxt + Gt wt , zt = Hxt + vtМодель Б:(14.33)TGt = 0 1 0 0= constс независимым дискретным белым шумом wt : wt ∼ N (0, Qt = qτs = const),где q — коэффициент диффузии соответствующего процесса броуновскогодвижения (винеровского процесса), м2 /c5 .
Таким образом может быть учтенанестационарность силы тяги двигателя ЛА или «скомпенсированы» другиенеточности идеализированной математической модели А (14.32), но при этомпотребуется задавать q. При q = 0 возвращаемся к модели А.В. Спроектировать, отладить и продемонстрировать в действии программу решения задачи оптимального оценивания состояния применительнок модели А и также модели Б в соответствии с вашим номером и содержанием проекта по табл. 14.1 (стр. 334) и вариантом исходных данных для негопо табл. 14.2 (стр.
334).Г. Спланировать вычислительный эксперимент с целью сравнить поведение оценок состояния объекта в модели А и в модели Б. При планированиипредусмотреть вариации параметров по табл. 14.2 на ±50%.Д. Получить от преподавателя восьмизначное число N , в котором каждая цифра ni ∈ {1, 2, 3, 4, 5, 6} обозначает номер варианта исходных данных по табл. 14.2. Начальный вектор состояния выбрать самостоятельноили спросить у преподавателя.Е.
Сделать выводы по выполненной работе. Выводы делаются на основании аналитического исследования и по результатам компьютерного моделирования. Они должны отражать результаты решения следующих вопросов:•построение дискретной модели,•оценка точности использованных алгоритмов оценивания,•оценка времени сходимости к установившемуся режиму фильтрации,•аналитическое определение установившегося решения P̄ алгебраического уравнения Риккати (14.9); для этого надо в уравнении (14.9) поло−жить P̄ , Pt+1= Pt− ,•оценка влияния вариаций параметров модели:33314 Ортогонализованные блочные алгоритмыТаблица 14.1. Номер и содержание проекта № 1012РКККФРККИФ3203223Номер проекта4567Разновидность алгоритма по разд.
14МККИФ КоККФ СКККФ СККИФ СМККИФ8СКоККФСтраница, где начинается описание алгоритма324326327328329329Таблица 14.2. Варианты исходных данных для проекта № 10№12345678••••ПараметрДиффузия q,м2 /c5Постоянная времени τ, cДисперсия σ12 ,(м/c2 )2Дисперсия σ22 ,м2Дисперсия P11 (0),м2Дисперсия P22 (0), (м/c)2Дисперсия P33 (0), (м/c2 )2Дисперсия P44 (0),м2Вариант исходных данных1234563000 340 400 300 3500 33500.05 0.65 0.80 0.90 0.10 0.121.00 2.25 4.00 6.25 6.00 5.50403025354550102030405060605040302010152025303540454035253015диффузии q,постоянной времени τ,дисперсий σ12, σ22,а также дисперсий P11 (0), P22 (0), P33(0), P44 (0)на скорость сходимости параметров фильтра к значениям установившегося режима оценивания.
Выводы обосновать протоколами компьютерного моделирования. Для сравнения рекомендуем использоватьmatlab.14.12Варианты задания на лабораторный проект № 10Данный проект — групповой, т. е. рассчитан на выполнение одного проекта двумя (максимум тремя) студентами.
Свое задание и содержание проекта студенты определяют по номеру назначенного проекта в табл. 14.1 иисходные данные для него берут из табл. 14.2.334ЗаключениеДанное учебное пособие разделено на две части, первая из которых содержит материал начального уровня, составляющий основу всех вычислительных методов алгебры, а вторая часть представляет собой курс повышенноготипа и, соответственно, может использоваться как учебный материал дляспециальных разделов вычислительной математики.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
