Получение данных моделирования из модели Simulink

Рассмотрим систему, описанную следующей моделью Simulink:

open_system('iddemsl1')
set_param('iddemsl1/Random Number','seed','0')

Красная часть является системой, синяя часть является контроллером, а опорный сигнал является стреловидной синусоидой (сигнал щебета). Значение шага расчета данных устанавливается на 0,5 секунды.

Эта система может быть представлена с помощью idpoly структура:

m0 = idpoly(1,0.1,1,1,[1 0.5],'Ts',0,'InputDelay',1,'NoiseVariance',0.01)

Давайте симулируем модель iddemsl1 и сохраните данные в iddata объект:

sim('iddemsl1')
dat1e = iddata(y,u,0.5); % The IDDATA object

Давайте сделаем вторую симуляцию режима в целях валидации.

set_param('iddemsl1/Random Number','seed','13')
sim('iddemsl1')
dat1v = iddata(y,u,0.5);

Давайте посмотрим на данные оценки, полученные во время первой симуляции:

plot(dat1e)

Оценка дискретных моделей с помощью данных моделирования

Начнем с оценки дискретной модели по умолчанию, чтобы получить некоторое предварительное представление о характеристиках данных:

m1 = n4sid(dat1e, 'best')     % A default order model

Проверьте, насколько хорошо модель воспроизводит данные валидации

compare(dat1v,m1)

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

ImpModel = impulseest(dat1e,[],'negative');
clf
h = impulseplot(ImpModel);
showConfidence(h,3)

ImpModel - модель конечная импульсная характеристика, чей порядок (количество коэффициентов) определяется автоматически. Мы также принимаем решение анализировать эффекты обратной связи путем вычисления импульсной характеристики для отрицательных задержек. Влияния от негативных лагов не все незначительны. Это связано с регулятором (обратной связью на выходе). Это означает, что оценка импульсной характеристики не может использоваться, чтобы определить временную задержку. Вместо этого создайте несколько моделей ARX низкого порядка с различными задержками и узнайте лучшую подгонку:

V = arxstruc(dat1e,dat1v,struc(1:2,1:2,1:10));
nn = selstruc(V,0) %delay is the third element of nn

Задержка определяется как 3 лага. (Это правильно: время ожидания 1 секунду дает две задержки, а ZOH-блокирует еще одну.) Соответствующая ARX-модель также может быть вычислена следующим образом:

m2 = arx(dat1e,nn)
compare(dat1v,m1,m2);

Уточнение оценки

Эти две модели m1 и m2 вести себя аналогично в симуляции. Давайте теперь попробуем подстроить порядки и задержки. Исправьте задержку до 2 (которая в сочетании с отсутствием сквозного соединения дает чистую задержку в 3 выборки) и найдите модель пространства состояний по умолчанию с этой задержкой:

m3 = n4sid(dat1e,'best','InputDelay',2,'Feedthrough',false);
% Refinement for prediction error minimization using pem (could also use
% |ssest|)
m3 = pem(dat1e, m3);

Посмотрим на оцененную системную матрицу

m3.a  % the A-matrix of the resulting model

Динамика третьего порядка выбирается автоматически, что вместе с 2 «дополнительными» задержками дает модель пространства состояний 5-го порядка.

Всегда желательно не слепо полагаться на автоматический выбор порядка. На них влияют случайные ошибки. Хороший способ - взглянуть на нули и полюсы модели, наряду с доверительными областями:

clf
h = iopzplot(m3);
showConfidence(h,2) % Confidence region corresponding to 2 standard deviations

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

m4 = ssest(dat1e,1,'Feedthrough',false,'InputDelay',2,'Ts',dat1e.Ts);
compare(dat1v,m4,m3,m1)

The compare график показывает, что простая модель первого порядка m4 дает очень хорошее описание данных. Таким образом, мы выберем эту модель в качестве нашего конечного результата.

Преобразование дискретной модели в непрерывное время (LTI)

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

mc = d2c(m4);
idtf(mc)

Было получено хорошее описание системы, как показано выше.

Непосредственная оценка модели непрерывного времени

Модель непрерывного времени также может быть оценена непосредственно. Дискретная модель m4 имеет 2 выборки входной задержки, которая представляет 1 секунду задержки. Используем ssest команда для этой оценки:

m5 = ssest(dat1e,1,'Feedthrough',false,'InputDelay',1);
present(m5)

Анализ неопределенности

Параметры модели m5 показывают высокие уровни неопределенности, хотя модель подходит для данных 87%. Это связано с тем, что модель использует больше параметров, чем абсолютно требуется, что приводит к потере уникальности в оценках параметров. Чтобы просмотреть истинный эффект неопределенности в модели, существует два возможных подхода:

Позвольте использовать оба подхода. Сначала мы оценим модель в канонической форме.

m5Canon = ssest(dat1e,1,'Feedthrough',false,'InputDelay',1,'Form','canonical');
present(m5Canon)

m5Canon использует каноническую параметризацию модели. Он подходит для данных оценки так же хорошо, как и модель m5. Он показывает небольшие неопределенности в значениях его параметров, что свидетельствует о его надежности. Однако, как мы видели m5большая неопределенность не обязательно означает «плохую» модель. Чтобы установить качество этих моделей, позвольте использовать просмотр их откликов во временных и частотных диапазонах с доверительными регионами, соответствующими 3 стандартным отклонениям. Мы также строим график исходной системы m0 для сравнения.

The диаграммы Боде.

clf
opt = bodeoptions;
opt.FreqScale = 'linear';
h = bodeplot(m0,m5,m5Canon,opt);
showConfidence(h,3)
legend show

График шага.

clf
showConfidence(stepplot(m0,m5,m5Canon),3)
legend show

Границы неопределенности для этих двух моделей фактически идентичны. Мы можем аналогично сгенерировать карту полюс-ноль (iopzplot) и Годограф Найквиста (nyquistplot) с доверительными областями для этих моделей.

idtf(m5)

Учет межпозиционного поведения при оценке в непрерывном времени

При сравнении моделей непрерывного времени, вычисленных из выборочных данных, важно учитывать поведение интерсампов входного сигнала. В данном примере вход в систему был частичным постоянным, из-за схемы Нуля order-Hold (zoh) в контроллере. Теперь удалите эту схему и рассмотрим действительно непрерывную систему. Входной и выходной сигналы все еще дискретизированы на 2 Гц, и все остальное то же самое:

open_system('iddemsl3')
sim('iddemsl3')
dat2e = iddata(y,u,0.5);

Модели в дискретном времени все еще хорошо справятся с этими данными, поскольку, когда они скорректированы с учетом измерений, они будут включать свойства дискретизации и интерсамплируемого входного поведения (для токового входа). Однако при построении моделей непрерывного времени знание свойств intersamples очень важно. Сначала создайте модель, так же как и для случая ZOH:

m6 = ssest(dat2e,1,'Feedthrough',false,'InputDelay',1,'Form','canonical');
idtf(m6)

Сравним предполагаемую модель (m6) против истинной модели (m0):

step(m6,m0)  % Compare with true system

Соглашение сейчас не так хорошо. Однако мы можем включать в объект данных информацию о входе. В качестве приближения позвольте описать его как кусочно-линейный (задержка первого порядка, FOH) между моментами дискретизации. Эта информация затем используется оценщиком для правильной выборки:

dat2e.Intersample = 'foh';
m7 = ssest(dat2e,1,'Feedthrough',false,'InputDelay',1,'Form','canonical');  % new estimation with correct intersample behavior
idtf(m7)

Рассмотрим переходную характеристику сравнение еще раз:

step(m7,m0)  % Compare with true system

Эта модель (m7) дает намного лучший результат, чем m6. На этом этот пример завершается.

bdclose('iddemsl1');
bdclose('iddemsl3');