Линейная аппроксимация сложных систем идентификацией
В этом примере показано, как получить линейные аппроксимации комплексной, нелинейной системы посредством линейной идентификации модели. Подход основан на выборе входного сигнала, который волнует систему. Линейная аппроксимация получена, подбирая линейную модель к симулированному ответу нелинейной модели для выбранного входного сигнала.
Этот пример использует Simulink®, Control System Toolbox™ и Simulink Control Design™.
Введение
Во многих ситуациях линейная модель получена упрощением более комплексной нелинейной системы при определенных локальных условиях. Например, высококачественная модель динамики самолета может быть описана подробной моделью Simulink. Общий подход, проявленный, чтобы ускорить симуляцию таких систем, изучите их локальное поведение о рабочей точке или спроектируйте компенсаторы, должен линеаризовать их. Если бы мы выполняем линеаризацию исходной модели аналитически о рабочей точке, это, в целом, дало бы к модели порядка целый, или близко к, количество состояний в исходной модели. Этот порядок может быть излишне высоким для класса входных параметров, с которыми это, как предполагается, используется в проекте системы управления или анализе. Следовательно мы можем считать альтернативный подход сосредоточенным вокруг сбора данных ввода - вывода от симуляции системы и использования его, чтобы вывести линейную модель только правильного порядка.
Аналитическая линеаризация модели F14
Рассмотрите модель F14. Это уже - линейная модель, но содержит производные блоки и источники воздействия, которое может влиять на природу его выхода. Мы можем "линеаризовать" его между его одним входным портом и этими двумя выходными портами можно следующим образом:
open_system('idF14Model') IO = getlinio('idF14Model') syslin = linearize('idF14Model',IO)
3x1 vector of Linearization IOs:
--------------------------
1. Linearization input perturbation located at the following signal:
- Block: idF14Model/Pilot
- Port: 1
- Signal Name: Stick Input
2. Linearization output measurement located at the following signal:
- Block: idF14Model/Gain5
- Port: 1
- Signal Name: Angle of Attack
3. Linearization output measurement located at the following signal:
- Block: idF14Model/Pilot G force
- Port: 1
- Signal Name: Pilot G force
syslin =
A =
Transfer Fcn Derivative Transfer Fcn Derivative1
Transfer Fcn -0.6385 0 689.4 0
Derivative 1 -1e+05 0 0
Transfer Fcn -0.00592 0 -0.6571 0
Derivative1 0 0 1 -1e+05
Actuator Mod 0 0 1.424 0
Alpha-sensor 0.001451 0 0 0
Stick Prefil 0 0 0 0
Pitch Rate L 0 0 1 0
Proportional 0 0 -0.8156 0
Actuator Mod Alpha-sensor Stick Prefil Pitch Rate L
Transfer Fcn -1280 0 0 0
Derivative 0 0 0 0
Transfer Fcn -137.7 0 0 0
Derivative1 0 0 0 0
Actuator Mod -20 2.986 -39.32 -1.67
Alpha-sensor 0 -2.526 0 0
Stick Prefil 0 0 -22.52 0
Pitch Rate L 0 0 0 -4.144
Proportional 0 -1.71 22.52 0.9567
Proportional
Transfer Fcn 0
Derivative 0
Transfer Fcn 0
Derivative1 0
Actuator Mod -3.864
Alpha-sensor 0
Stick Prefil 0
Pitch Rate L 0
Proportional 0
B =
Stick Input
Transfer Fcn 0
Derivative 0
Transfer Fcn 0
Derivative1 0
Actuator Mod 0
Alpha-sensor 0
Stick Prefil 1
Pitch Rate L 0
Proportional 0
C =
Transfer Fcn Derivative Transfer Fcn Derivative1
Angle of Att 0.001451 0 0 0
Pilot G forc -3106 3.106e+08 7.083e+04 -7.081e+09
Actuator Mod Alpha-sensor Stick Prefil Pitch Rate L
Angle of Att 0 0 0 0
Pilot G forc 0 0 0 0
Proportional
Angle of Att 0
Pilot G forc 0
D =
Stick Input
Angle of Att 0
Pilot G forc 0
Continuous-time state-space model.
syslin является моделью с 2 выходными параметрами, одним входом и 9 состояниями. Это вызвано тем, что путь к линеаризации от входа "Stick Input" до этих двух выходных параметров имеет 9 состояний в исходной системе. Мы можем проверить это при помощи operpoint:
operpoint('idF14Model')
Operating point for the Model idF14Model.
(Time-Varying Components Evaluated at time t=0)
States:
----------
(1.) idF14Model/Actuator Model
x: 0
(2.) idF14Model/Aircraft Dynamics Model/Transfer Fcn.1
x: 0
(3.) idF14Model/Aircraft Dynamics Model/Transfer Fcn.2
x: 0
(4.) idF14Model/Controller/Alpha-sensor Low-pass Filter
x: 0
(5.) idF14Model/Controller/Pitch Rate Lead Filter
x: 0
(6.) idF14Model/Controller/Proportional plus integral compensator
x: 0
(7.) idF14Model/Controller/Stick Prefilter
x: 0
(8.) idF14Model/Dryden Wind Gust Models/Q-gust model
x: 0
(9.) idF14Model/Dryden Wind Gust Models/W-gust model
x: 0
x: 0
Inputs: None
----------
Это могло заказать уменьшаться при тихом поддержании точности ответам для выбранной прямоугольной волны ("Вход Stick") вход?
Подготовка идентификационных данных
Симулируйте модель и регистрируйте входную прямоугольную волну (u) и выходные параметры "Angle of attack" (y1) и "Сила пилота Г" (y2) в течение 0:30 промежуток второго раза. Эти данные, после интерполяции к однородно расположенному с интервалами временному вектору (шаг расчета 0,0444 секунд), хранятся в "idF14SimData.mat" файле.
load idF14SimData Z = iddata([y1, y2],u,Ts,'Tstart',0); Z.InputName = 'Stick input'; Z.InputUnit = 'rad/s'; Z.OutputName = {'Angle of attack', 'Pilot G force'}; Z.OutputUnit = {'rad', 'g'}; t = Z.SamplingInstants; subplot(311) plot(t,Z.y(:,1)), ylabel('Angle of attack (rad)') title('Logged Input-Output Data') subplot(312) plot(t,Z.y(:,2)), ylabel('Pilot G force (g)') subplot(313) plot(t,Z.u), ylabel('Stick input (rad/s)') axis([0 30 -1.2 1.2]) xlabel('Time (seconds)')
Оценка моделей в пространстве состояний
Оцените модели в пространстве состояний порядков 2 - 4 с помощью ssest команда. Мы конфигурируем оценку, чтобы использовать особое внимание "симуляции" и принять решение не оценить компонент воздействия модели.
opt = ssestOptions('Focus','simulation'); syslin2 = ssest(Z, 2, 'DisturbanceModel', 'none', opt); syslin3 = ssest(Z, 3, 'DisturbanceModel', 'none', opt); syslin4 = ssest(Z, 4, 'DisturbanceModel', 'none', opt);
Сравните производительность линеаризовавшей модели syslin и три идентифицированных модели к данным. Обратите внимание на то, что syslin является моделью SS в то время как syslin2, syslin3 и syslin4 модели IDSS.
syslin.InputName = Z.InputName;
syslin.OutputName = Z.OutputName; % reconcile names to facilitate comparison
clf
compare(Z, syslin, syslin2, syslin3, syslin4)
График показывает что 3-я модель порядка (syslin3) работает вполне прилично линейной аппроксимацией самолета, динамического о его значении по умолчанию (t=0) условия работы. Это соответствует данным немного лучше, чем тот, полученный аналитической линеаризацией (syslin). Если оригинал idF14Model линеен, почему его линеаризация не заканчивается, syslin, произвести 100%-ю подгонку к данным? Существует две причины этого:
Измеренные выходные параметры затронуты порывами ветра, что означает, что регистрируемые выходные параметры не являются просто функцией входа Stick. Существуют воздействия, влияющие на него.
Блок "Pilot G force" использует производные блоки, линеаризация которых зависит от значения постоянной времени "c". "c", должен быть малым (мы используем 1e-5), но не нуль. Ненулевое значение для c вводит ошибку приближения в линеаризации.
Давайте просмотрим параметры модели syslin3 который, кажется, получает ответы вполне прилично:
syslin3
syslin3 =
Continuous-time identified state-space model:
dx/dt = A x(t) + B u(t) + K e(t)
y(t) = C x(t) + D u(t) + e(t)
A =
x1 x2 x3
x1 -1.006 -2.029 0.5842
x2 8.284 -19.39 -5.611
x3 2.784 -12.63 -6.956
B =
Stick input
x1 0.2614
x2 5.512
x3 3.606
C =
x1 x2 x3
Angle of att -8.841 -0.5347 1.402
Pilot G forc -86.42 -15.85 66.12
D =
Stick input
Angle of att 0
Pilot G forc 0
K =
Angle of att Pilot G forc
x1 0 0
x2 0 0
x3 0 0
Parameterization:
FREE form (all coefficients in A, B, C free).
Feedthrough: none
Disturbance component: none
Number of free coefficients: 18
Use "idssdata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using SSEST on time domain data "Z".
Fit to estimation data: [98.4;97.02]%
FPE: 2.367e-05, MSE: 0.1103
Упрощение модели сокращением и оценкой
Мы можем также проявить подход сокращения порядка линеаризовавшей модели syslin и совершенствование параметров упрощенной модели, чтобы лучше всего соответствовать данным Z. Чтобы выяснить хорошее значение для уменьшаемого порядка, мы используем hsvd:
[S, BalData] = hsvd(syslin); clf; bar(S)
Столбчатая диаграмма показывает, что сингулярные значения вполне малы для состояний 4 и выше. Следовательно 3-й порядок может быть оптимальным для сокращения.
sysr = balred(syslin,3,BalData)
opt2 = bodeoptions; opt2.PhaseMatching = 'on';
clf; bodeplot(sysr,syslin,opt2)
sysr =
A =
x1 x2 x3
x1 -1.205 -5.28 45.51
x2 -0.5497 4.546 -27.37
x3 -0.06248 7.327 -27.93
B =
Stick input
x1 -328.5
x2 -1048
x3 682.7
C =
x1 x2 x3
Angle of att -0.0006108 -0.0005806 0.000794
Pilot G forc -0.007392 -0.03346 0.1282
D =
Stick input
Angle of att -0.03784
Pilot G forc -1.617
Continuous-time state-space model.
Диаграмма Боде показывает хорошей точности до 10 рад/с. sysr может эмулировать ответы, столь же хорошие как исходная модель с 9 состояниями как показано compare график:
compare(Z, sysr, syslin)
Давайте совершенствуем параметры sysr улучшить его подгонку к данным. Для этой оценки мы выбираем метод поиска "Levenberg-Marquardt" и изменяем максимальное количество допустимых итераций к 10. Этот выбор был сделан на основе некоторого метода проб и ошибок. Мы также включаем отображение прогресса оценки.
opt.Display = 'on'; opt.SearchMethod = 'lm'; opt.SearchOptions.MaxIterations = 10; sysr2 = ssest(Z, sysr, opt) compare(Z, sysr2)
sysr2 =
Continuous-time identified state-space model:
dx/dt = A x(t) + B u(t) + K e(t)
y(t) = C x(t) + D u(t) + e(t)
A =
x1 x2 x3
x1 -2.227 -5.344 45.5
x2 0.03107 3.495 -27.19
x3 0.6561 6.546 -29.82
B =
Stick input
x1 -328.5
x2 -1048
x3 682.7
C =
x1 x2 x3
Angle of att -0.0006251 -0.0003531 -0.001759
Pilot G forc -0.0112 -0.02443 0.1133
D =
Stick input
Angle of att 0.003051
Pilot G forc 0.6397
K =
Angle of att Pilot G forc
x1 0 0
x2 0 0
x3 0 0
Parameterization:
FREE form (all coefficients in A, B, C free).
Feedthrough: yes
Disturbance component: none
Number of free coefficients: 20
Use "idssdata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using SSEST on time domain data "Z".
Fit to estimation data: [98.78;97.03]%
FPE: 1.434e-05, MSE: 0.1099
Усовершенствованная модель sysr2 соответствует ответу модели F14 вполне хорошо (приблизительно 99% для первого выхода и приблизительно 97% для второго).
Заключения
Мы показали альтернативный подход к аналитической линеаризации для получения линейных аппроксимаций сложных систем. Результаты выведены для конкретного входного сигнала и, строго говоря, применимы к тому входу только. Чтобы улучшить применимость результатов к различным входным профилям, мы могли выполнить несколько симуляций с помощью различных типов входных параметров. Мы могли затем объединить получившиеся наборы данных в один набор данных мультиэксперимента (см. iddata/merge) и используйте его в оценках. В этом примере мы использовали комплекс, но линейную систему для удобства. Реальная выгода этого подхода была бы замечена для нелинейных систем.
Мы также показали подход для сокращения порядка линейной системы при хранении упрощенной модели верной симулированному ответу исходной модели Simulink. Роль упрощенной модели sysr должен был обеспечить исходное предположение для предполагаемой модели sysr2. Подход также подсвечивает то, что любая линейная система, включая те из различного класса, может использоваться в качестве первоначальной модели для оценок.
bdclose('idF14Model')
Дополнительная информация
Для получения дополнительной информации об идентификации динамических систем с System Identification Toolbox посещают страницу информации о продукте System Identification Toolbox.