Введение

Надежная панель инструментов управления предлагает инструменты для работы с крупными моделями, такие как:

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

Процесс сокращения модели

Задача сокращения модели обычно включает в себя следующие шаги:

Мы можем легко получить доступ к этим алгоритмам через интерфейс верхнего уровня reduce. Способы используют различные показатели «близости» между исходной и уменьшенной моделями. Выбор зависит от приложения. Давайте попробуем каждому из них исследовать свои относительные достоинства.

Пример: Модель движения жесткого тела здания

В этом примере мы применяем методы сокращения к модели здания Университетской больницы Лос-Анджелеса. Модель взята из рабочей записки SLICOT 2002-2, «Сборник эталонных примеров для модельного сокращения линейных временных инвариантных динамических систем», Y. Chahlaoui и P.V. Dooren. Имеет восемь этажей, каждый с тремя степенями свободы - два смещения и одно вращение. Мы представляем соотношение «вход-выход» для любого из этих смещений с использованием модели с 48 состояниями, где каждое состояние представляет смещение или его скорость изменения (скорость).

Давайте загрузим данные для примера:

load buildingData.mat

Изучение динамики завода

Начнем с анализа частотной характеристики модели:

bode(G)
grid on

Рис. 1: Схема Bode для анализа частотной характеристики

Как видно из частотной характеристики модели, существенная динамика системы лежит в диапазоне частот от 3 до 50 радиан/с. Величина падает как в очень низком, так и в высокочастотном диапазоне. Наша цель - найти модель низкого порядка, которая сохранит содержание информации в этом частотном диапазоне до приемлемого уровня точности.

Вычисление сингулярных значений Hankel

Чтобы понять, какие состояния модели можно безопасно отбросить, посмотрите на сингулярные значения Ханкеля модели:

hsv_add = hankelsv(G);
bar(hsv_add)
title('Hankel Singular Values of the Model (G)');
xlabel('Number of States')
ylabel('Singular Values (\sigma_i)')
line([10.5 10.5],[0 1.5e-3],'Color','r','linestyle','--','linewidth',1)
text(6,1.6e-3,'10 dominant states.')

Рис. 2: Сингулярные значения Hankel модели (G).

График сингулярных значений Ханкеля предполагает, что в этой системе существует четыре доминирующих режима. Однако вклад остальных режимов по-прежнему значителен. Мы проведем линию в 10 состояниях и отбросим оставшиеся, чтобы найти уменьшенную модель 10-го порядка Gr это лучше всего соответствует исходной системе G.

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

Функция reduce является шлюзом для всех подпрограмм сокращения модели, доступных на панели инструментов. Мы будем использовать опцию усечения баланса квадратного корня по умолчанию ('balancmr') reduce в качестве первого шага. Этот метод использует «аддитивную» ошибку, связанную с уменьшением, что означает, что он пытается поддерживать абсолютную аппроксимационную ошибку равномерно малой по частотам.

% Compute 10th-order reduced model (reduce uses balancmr method by default)
[Gr_add,info_add] = reduce(G,10);

% Now compare the original model G to the reduced model Gr_add
bode(G,'b',Gr_add,'r')
grid on
title('Comparing Original (G) to the Reduced model (Gr\_add)')
legend('G - 48-state original ','Gr\_add - 10-state reduced','location','northeast')

Рис. 3: Сравнение оригинала (G) с уменьшенной моделью (Gr_add)

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

Как видно из диаграммы Боде на фиг.3, уменьшенная модель достаточно хорошо фиксирует резонансы ниже 30 рад/с, но совпадение в низкочастотной области (< 2 рад/с) является плохим. Также уменьшенная модель не полностью фиксирует динамику в диапазоне частот 30-50 рад/с. Возможным объяснением больших ошибок на низких частотах является относительно низкий коэффициент усиления модели на этих частотах. Следовательно, даже большие ошибки на этих частотах мало способствуют общей ошибке .

Чтобы обойти эту проблему, мы можем попробовать метод мультипликативной ошибки, такой как bstmr. Этот алгоритм подчеркивает относительные ошибки, а не абсолютные. Поскольку относительные сравнения не работают, когда коэффициенты усиления близки к нулю, нам необходимо добавить порог минимального коэффициента усиления, например, путем добавления коэффициента усиления прохождения D к нашей исходной модели. Предполагая, что мы не обеспокоены ошибками при выигрыше ниже -100 дБ, мы можем установить прохождение 1e-5.

GG = G;
GG.D = 1e-5;

Теперь рассмотрим сингулярные значения для мультипликативных (относительных) ошибок (используя опцию «мульт» hankelsv)

hsv_mult = hankelsv(GG,'mult');
bar(hsv_mult)
title('Multiplicative-Error Singular Values of the Model (G)');
xlabel('Number of States')
ylabel('Singular Values (\sigma_i)')

Рис. 4: Сингулярные значения мультипликативной ошибки модели (G)

Модель 26-го порядка выглядит многообещающе, но ради сравнения с предыдущим результатом давайте придерживаться сокращения 10-го порядка.

% Use bstmr algorithm option for model reduction 
[Gr_mult,info_mult] = reduce(GG,10,'algorithm','bst');

%now compare the original model G to the reduced model Gr_mult
bode(G,Gr_add,Gr_mult,{1e-2,1e4}), grid on
title('Comparing Original (G) to the Reduced models (Gr\_add and Gr\_mult)')
legend('G - 48-state original ','Gr\_add (balancmr)','Gr\_mult (bstmr)','location','northeast')

Рис. 5: Сравнение оригинала (G) с уменьшенными моделями (Gr_add и Gr_mult)

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

step(G,Gr_add,Gr_mult,15) %step response until 15 seconds
legend('G: 48-state original ','Gr\_add: 10-state (balancmr)','Gr\_mult: 10-state (bstmr)')

Рис. 6: Пошаговые ответы трех моделей

Проверка результатов

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

где сумма превышает все отброшенные сингулярные значения Ханкеля G (позиции 11-48 из hsv_add). Мы можем проверить, что эта граница удовлетворена, сравнивая две стороны этого неравенства:

norm(G-Gr_add,inf) % actual error

ans = 6.0251e-04

% theoretical bound (stored in the "ErrorBound" field of the "INFO" 
% struct returned by |reduce|)
info_add.ErrorBound

ans = 0.0047

Для методов мультипликативной ошибки, таких как bstmr, аппроксимационная ошибка измеряется пиковым усилением частоты модели относительной ошибки G\(G-Greduced). Граница ошибки выглядит следующим образом

где сумма превышает отброшенные мультипликативные сингулярные значения Ханкеля, вычисленные hankelsv(G,'mult'). Снова мы можем сравнить эти границы для уменьшенной модели Gr_mult

norm(GG\(GG-Gr_mult),inf) % actual error

ans = 0.5949

% Theoretical bound 
info_mult.ErrorBound

ans = 546.1730

Постройте график относительной ошибки для подтверждения

bodemag(GG\(GG-Gr_mult),{1e-2,1e3})
grid on
text(0.1,-50,'Peak Gain: -4.6 dB (59%) at 17.2 rad/s')
title('Relative error between original model (G) and reduced model (Gr\_mult)')

Рис. 7: Относительная ошибка между исходной моделью (G) и уменьшенной моделью (Gr_mult)

На приведенном выше графике относительных ошибок имеется до 59% относительных ошибок при 17,2 рад/с, что может быть больше, чем мы готовы принять.

Выбор минимального порядка, совместимого с требуемым уровнем точности

Для повышения точности Gr_mult, нам нужно увеличить заказ. Для достижения по крайней мере 5% относительной точности, какой самый низкий порядок мы можем получить? Функция reduce может автоматически выбрать модель самого низкого порядка, совместимую с нашим требуемым уровнем точности.

% Specify a maximum of 5% approximation error
[Gred,info] = reduce(GG,'ErrorType','mult','MaxError',0.05);
size(Gred)

State-space model with 1 outputs, 1 inputs, and 35 states.

Алгоритм выбрал уменьшенную модель с 34 состояниями Gred. Сравните фактическую ошибку с теоретической границей:

norm(GG\(GG-Gred),inf)

ans = 0.0068

info.ErrorBound

ans = 0.0342

Рассмотрим относительную величину ошибки как функцию частоты. Более высокая точность достигается за счет большего порядка моделей (от 10 до 34). Следует отметить, что фактическая максимальная погрешность составляет 0,6%, что намного меньше целевого значения 5%. Это расхождение является результатом функции bstmr использование границы ошибки, а не фактической ошибки для выбора заказа.

bodemag(GG\(GG-Gred),{1,1e3})
grid on
text(5,-75,'Peak Gain: -43.3 dB (0.6%) at 73.8 rad/s')
title('Relative error between original model (G) and reduced model (Gred)')

Рис. 8: Относительная ошибка между исходной моделью (G) и уменьшенной моделью (Grd)

Сравнение ответов Боде

bode(G,Gred,{1e-2,1e4})
grid on
legend('G - 48-state original','Gred - 34-state reduced')

Рис. 9: Блок-схема исходной модели с 48 состояниями и уменьшенной модели с 34 состояниями

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

step(G,'b',Gred,'r--',15) %step response until 15 seconds
legend('G: 48-state original ','Gred: 34-state (bstmr)')
text(5,-4e-4,'Maximum relative error <0.05')

Рисунок 10: График ответа шага исходной модели с 48 состояниями и уменьшенной модели с 34 состояниями