УДК 517.95
ПАРАМЕТРИЧЕСКАЯ ИДЕНТИФИКАЦИЯ ДИНАМИЧЕСКИХ МОДЕЛЕЙ МОРСКИХ СУДОВ Н.А. Жабко
Рассматривается задача идентификации линейных моделей морских судов с использованием методов, базирующихся на минимизации ошибки предсказания и применении эталонной модели, в условиях действия возмущений
Ключевые слова: параметрическая идентификация, адаптивная идентификация, эталонная модель, внешнее возмущение
Введение
Эффективность решения задач управления динамическими объектами существенно зависит от того, насколько адекватное математическое описание используется при формировании законов управления. В качестве математических моделей, которые описывают динамику морских подвижных объектов в разных режимах функционирования, обычно принимаются системы дифференциальных уравнений, зависящие от условий плавания и характеристик конкретного судна. При этом аналитические зависимости, определяющие динамические модели, то есть вид и структуру соответствующих дифференциальных уравнений, часто можно считать известными, в то время как параметры объекта, входящие в уравнения, неизвестны или рассчитаны недостаточно точно. В таком случае заданные параметры оцениваются или уточняются в результате решения задачи параметрической идентификации на основе экспериментов с объектом в нормальном или активном (специально организованном тестовом) режиме функционирования. Наличие внешних возмущений в виде волн, ветра, течения, а также шумов в измерениях при проведении экспериментов осложняет получение точных оценок параметров рассматриваемых моделей, поэтому изучение возможности идентификации параметров с учетом возмущающих воздействий представляет особый интерес. Предполагается, что эксперименты с объектом проводятся на тихой воде, однако на объект действуют медленно изменяющиеся с течением времени или постоянные (низкочастотные) возмущения, вклад высокочастотных возмущений считается приемлемым для получения адекватного результата идентификации. Как известно, необходимым этапом при выполнении идентификации является предварительная обработка экспериментальных данных с целью возможного исключения влияния возмущающих факторов.
Жабко Наталия Алексеевна - СПбГУ, канд. физ.-мат. наук, доцент, e-mail: ztasha@mail.ru
Однако такая обработка не всегда может привести к улучшению данных, кроме того, идентификация может проводиться в режиме реального времени с применением адаптивных схем на основе постоянно обновляемого набора данных, что исключает их предварительную обработку. В этом случае для улучшения результатов идентификации в модель, используемую для идентификации, часто вводится модель шума того или иного вида и, соответственно, дополнительный набор параметров, подлежащих оценке совместно с параметрами модели. В статье обсуждается применимость в такой ситуации известных схем идентификации для определения коэффициентов линейных моделей движения, используемых при решении многих задач управления морских судов [2,5,6]. В частности, рассматривается возможность адаптивной идентификации параметров в линейных моделях динамики морских судов при наличии внешних возмущений.
Минимизация ошибки предсказания
Алгоритмы параметрической идентификации, базирующиеся на минимизации ошибки предсказания, предполагают обычно [1, 2, 4] что объект функционирует в дискретном времени в соответствии с уравнениями
у^) = фт ()0 + е(Г). (1)
Здесь ^ - текущий момент времени, у(^) е Ек -вектор выхода объекта в момент ^, е(^) - вектор шума, 0 е Ер - вектор оцениваемых параметров, (рт (^) - матрица, определяемая в каждый момент времени заданной вход-выходной последовательностью {и(^), у(^)} , ^ е [1, N ]
значений входного и(^) е Ет и выходного сигналов в рассматриваемые дискретные моменты времени.
Тогда модель, представляющая собой прогноз значения выходного сигнала в момент I, задаётся уравнениями
уО, ö) = j (t )ö,
(2)
а задача идентификации сводится к минимизации функционала
1 N
v (°) = N S1 (у(< ) - У(' ■ °)).
t=1
где l ( х) - скалярная неотрицательно-
определенная функция:
V(0) ® min .
BeD
Здесь D - допустимое множество оценок параметров. В частности, для метода наименьших квадратов имеем
1 N1 ~ т
V (0) = й 2 ö (у (t) - y(t ■ ö))t (у (t) - y(t. ö)).
N t=i 2
Если оценка искомых параметров, входящих в уравнения (1), определена путем применения одной из схем, используемых при идентификации на основе минимизации ошибки предсказания, и не является удовлетворительной, она может быть улучшена путем введения дополнительных слагаемых в (2) и их оценивания одновременно с параметрами модели объекта.
Будем считать, что объект функционирует в соответствии с уравнениями
У (t) = jT (t)ö + f (t,7) + e(tX где f (t, g) - низкочастотные возмущения,
ge Eq - неизвестные параметры.
В простейшем варианте можно допустить, что возмущающее воздействие как функция времени на промежутке идентификации представимо в виде
f (t, g) = at + at +... + a,
где a}. e Ek, j = 1, k - постоянные векторы неизвестных коэффициентов, объединенные в вектор дополнительных параметров g.
Дополнительные параметры входят в уравнения объекта линейно, поэтому ошибка предсказания в такой ситуации может быть представима в виде (2) с новым вектором параметров
( 0 ^
0* = 0 .
\g)
Предположим теперь, что известна модель возмущений f(t, g) , которые действуют на объект, то есть аналитическая зависимость от параметров, которая является нелинейной. В этом случае можно выделить линейную часть функции f(t, g) в окрестности некоторого начального вектора g = g0, считая t параметром, и перейти к соотношениям, аналогичным (2). При-
менение соответствующей схемы идентификации для такой модифицированной модели может привести к получению улучшенной оценки параметров модели 0 . Осуществляя подобную процедуру последовательно несколько раз для сформированного вновь на предыдущем этапе начального вектора у = у0 , можно получить дополнительное улучшение оценки параметров 0 модели объекта.
Адаптивная идентификация параметров
Адаптивная идентификация параметров объекта может быть осуществлена с использованием рекуррентных схем идентификации, основанных на методе минимизации ошибки предсказания [1, 2, 4]. Такие схемы предполагают, что оценка параметров уточняется на каждой итерации при переходе от текущего момента времени к следующему, то есть происходит адаптация оценок параметров под изменяющиеся условия функционирования объекта. В частности, для подстройки под текущий уровень возмущений в модель могут быть введены дополнительные слагаемые, как указано в предыдущем пункте, параметры в которых следует также уточнять совместно с параметрами модели на каждой итерации.
Альтернативный подход связан с адаптивной подстройкой к уровню возмущений на основе адаптивного метода с использованием эталонной модели, опирающегося на применение функций Ляпунова. В соответствии с [2,3] предположим, что линейная математическая модель объекта управления представляется системой дифференциальных уравнений
х = Ах + Ви, (3)
где х е Еп - вектор состояния объекта, и е Ет - входной сигнал, подаваемый на объект, А и В - неизвестные матрицы с постоянными компонентами соответствующих размеров, матрица А гурвицева.
Введем также эталонную модель, которая описывается линейной системой
Хт = А тХт + Вти , (4)
где х т е Е п - вектор состояния эталонной модели, и е Ет - входной сигнал эталонной модели, который совпадает с подаваемым на вход объекта (3). Матрицы А т , В т имеют те же
размерности, что и соответствующие матрицы в системе (3), и могут выбираться произвольно, что позволяет решать задачу идентификации, основанную на выборе матриц Ат , Вт , обеспечивающих наименьшее расхождение между
динамикой объекта (3) и его модели (4), мерой которого является разность
e(t) = Хm (t) - X(t^ (5)
векторов состояний модели (4) и объекта (3). При этом движение объекта и модели происходит под действием одного и того же заданного входного сигнала u = u(t) . Для формирования меры расхождения в виде (5) все компоненты векторов x и u подлежат непосредственному измерению в каждый момент времени.
Вектор ошибки e(t) из (5) удовлетворяет системе уравнений
e = Ae + Axm + B u, (6)
где A = Am - A , B = Bm - B .
Асимптотическое стремление к нулю ошибки e(t) по состоянию объекта и модели будет означать также стремление матриц Am , Bm к матрицам A , B управляемого объекта, то есть к истинным значениям оцениваемых матриц.
Зададим функцию Ляпунова на решениях системы (6) в следующем виде:
V (e, A, B) = eT Pe + tr(AT G-1 A) + tr(BT G-1B),
где матрицы P и Gj, G2 являются постоянными, симметрическими и положительно определенными.
Производная от функции Ляпунова в силу системы (6) имеет вид:
V (e, A, B) = eT ( AT P + PA)e +
(4)
2tr(AT (PexTm + G-1 A )) + 2tr(BT (Peu T + G -1B )).
Здесь матрица P должна удовлетворять матричному уравнению Ляпунова
AT P + PA = -Q,
для произвольной положительно определенной симметрической матрицы Q .
Тогда производная V(e, A, B) является
(6)
отрицательно определенной функцией, если матрицы Am и Bm в эталонной модели (6) выбраны таким образом, чтобы в каждый момент времени выполнялись равенства
'-1~ A T> „..T , /-l-1l'
Pexm + G-1 A = 0, Peu + G-1B = 0.
В этом случае требуемые условия
е(0 — 0 , Ат (*) — А(0 и Вт (*) — В(*)
* — ¥ * — ¥ * — ¥
будут выполняться.
С учетом тождеств
А ° А„ и В ° В„
в силу постоянства компонент матриц А и В объекта (3), задача идентификации сводится к решению нелинейной системы дифференциальных уравнений:
х„, = Ах + В и,
A m =-G,Pexm,
m 1 m5
(7)
B.
= -G 2Peu
В соответствии с таким методом идентификация проводится адаптивно, то есть значения идентифицируемых параметров уточняются на каждом шаге интегрирования указанной системы дифференциальных уравнений.
При наличии внешних возмущений будем считать, что объект описывается уравнениями х = Ах + Ви + Г , где Г е Еп - вектор возмущений, для любых значений * е [0, т ].
В этом случае оценивание параметров в соответствии с системой (7) может не приводить к удовлетворительному результату. Для расширения возможностей по оцениванию параметров модели с помощью метода, можно предложить модификацию соответствующего алгоритма, основанную на адаптивной подстройке под текущее возмущение путем введения в модель и оценки дополнительных параметров.
Введем в модель дополнительный вектор Гт е Еп, элементы которого будем считать не-
: A m x m + B m u + fm
известными:
х т
Функция Ляпунова в этом случае может быть составлена в виде
V (е, А, В) = еТ Ре+и-(АТ о-1 А)+
+ 1Г(ВТ о - В) + £ мX, где М - постоянная, симметрическая и положительно определенная матрица. При этом разность е(*) = хт (* ) - х(*) удовлетворяет системе
e = Ae + Ax + Bu.
(8)
Тогда производная функции V в силу системы (8) принимает вид
V (е, А, В) = еТ (АТ Р + РА )е +
(8)
21х(А т (РехТ, + О-1 А )) + 21х(ВТ (Реи'+О -1В )) + 21-^ (Ре + М-1&т )).
Выбирая вектор Гт таким образом, чтобы выполнялось равенство Ре + М 1т = 0 , приходим к системе дифференциальных уравнений, к
решению которой сводится процесс идентификации в этом случае
х„ = А„х + В и,
т т т
Т
А т = -С:РЄХ т В_ = -^РеиТ
(9)
Г = -МРе.
т
Примеры идентификации модели морского судна
Иллюстрацию возможного преимущества проведения идентификации с подстройкой под возмущающие воздействия посредством введения дополнительных слагаемых в соответствии с приведенными выше рассуждениями, проведем для линейной модели динамики судна в горизонтальной плоскости [6]
¡5 = а11Ь + а12о + Ъ18, о& = а 215 + а 22 о + Ъ28, ф = а>.
Здесь 5 - угол дрейфа судна, о - угловая скорость вращения по курсу, ф - угол курса, 8 -отклонение вертикального руля. В результате решения задачи идентификации тем или иным методом должны быть найдены заранее неизвестные значения коэффициентов а11, а12, а21,
а22 , Ъ1, Ъ2 .
Идентификация проводится на основе измеряемых данных о выходных переменных объекта, полученных при движении судна в специальных режимах при подаче на вход тестовых сигналов с учетом допустимых ограничений на динамические переменные. Для данного примера использованы данные, полученные при проведении экспериментов с нелинейной моделью объекта.
Измеряемыми переменными считаются У1 = 5 и у2 = о, и выходной вектор формируется как у = (у у2 )Т .
В рамках такой постановки задачи неизвестные параметры совпадают с коэффициентами линейной системы второго порядка
5 = а115 + а12о+ Ъ18, о = а215 + а22о + Ъ28,
для полностью измеряемого вектора состояния
у = х = (У1 У2)Т =(5 о)Т .
Истинными значения параметров являются а11 = -0.0936 , а12 = -6.3400 , а21 = -0.0048 , а22 = -0.7170, Ъ1 = -0.1900, Ъ2 = 0.0160.
(10)
Примем в качестве тестового входного сигнала зависимость
8(*) = 0.3вт —*.
15
Дополнительно введем в модель (10) возмущающее воздействие в виде
Г (ї)
БІПІ
где а1, а2 - некоторые постоянные значения.
Проведя дискретизацию модели, можно перейти к описанию динамики объекта в дискретном времени
|Х(^+і) = А й х(ік ) + В й и(ік ), іу(їк ) = Сс х(їк ),
(11)
где Cd = С = Е, А л = Е + Ак, Вл = Вк, к -шаг дискретизации, с матрицами А и В :
Ґ
А =
а
V а21
а
12
Л
а
В=
22 У
Г ь ]
V Ь1 у
Тогда, объединяя элементы матриц А й и Вй в вектор
0 = [а1Ы, а12й , а21й, а22й , ЬЫ , Ь2с1 ] ,
от системы (11) можно легко перейти к модели вида (2).
Для заданных значений а1 = 0.01, а2 = 0 идентификация на промежутке [0,100с] с применением метода наименьших квадратов приводит к оценкам вектора параметров а11 = 0.0827 , а12 = 16.8767 , а21 = -0.0048 , а22 = -0.7170, Ь1 = -0.6074, Ь2 = 0.0160.
В то же время применение того же алгоритма с введением дополнительных слагаемых в соответствии с п.2 позволяет получить достаточно хорошие оценки параметров
а11 = -0.0937
а
12
-6.3424
а
21
-0.0048
а
= -0.7170, Ь1 = -0.1902, Ь2 = 0.0160.
Для иллюстрации работы адаптивного алгоритма (9) для формирования тестового режима введем командный сигнал в виде
11 ф7 (ї) = БІП-----------ї + БІП-----
25 100
ї
Учтем динамику
привода руля 8 = и и на основе сигнала (рг (ї) синтезируем закон управления
и = к1Ь + к 2т+ к3(ф-фг (ї)) + к48 с коэффициентами, гарантирующими асимптотическую устойчивость системы
¡5 = a11b + a12w + Ъ1д, w = a 21b + a 22w + b2d, ф = w,
d = k1b + к 2w + k3 (ф - фг (ґ)) + k4S.
Рис.1 Идентификация параметров Ъ1, Ъ2
Рис.2 Идентификация параметров а11, а12,
а21, а22
Для значений а1 = 0.0001, а2 = 0.00001 адаптивная идентификация параметров модели на промежутке [0,4000с] без дополнительной настройки уровня возмущений иллюстрируется рис.1-2, включение такой настройки даёт процессы, представленные на рис. 3-4.
Таким образом, предлагаемая методика применения идентификации может приводить к улучшению идентификационного процесса.
Рис.3 Идентификация параметров Ъ1, Ъ2 с настройкой возмущения
Рис.4 Идентификация параметров a11, a12,
a21, a22 с настройкой возмущения
Литература
1.Льюинг Л. Идентификация систем. Теория для пользователей: Пер. с англ./ Под ред. Я.З. Цыпкина М.: Наука, 1991. 432 с.
2.Thor I. Fossen. Guidance and Control of Ocean Vehicles. John Wiley & SONS, 1994. 480 c.
3.Petros A. Ioannou. Robust Adaptive Control. Prentice Hall, 1995. 822 c.
4.Hung D. Nguyen. Recursive identification of ship manoeuvring dynamics and hydrodynamics // ANZIAM Journal. 2008. Vol 49 (EMAC2007), PP. 717-732.
5.Веремей Е.И., Корчанов В.М., Коровкин М.В., По-гожев С.В. Компьютерное моделирование систем управления движением морских подвижных объектов. СПб, 2002.
6. Веремей Е.И. Синтез законов многоцелевого управления движением морских объектов // Гироскопия и навигация. 2009. № 4. С. 3 -14.
Санкт-Петербургский государственный университет
PARAMETER IDENTIFICATION OF SHIP DYNAMICAL MODELS N.A. Zhabko
The problem of parameter identification of ship dynamics linear models is considered based on prediction error minimization and model reference adaptive system construction techniques with an external perturbation
Key words: parameter identification, adaptive identification, reference model, external perturbation