Понижение порядка сложных моделей с помощью инструментов Robust Control Toolbox
Автор: Бильфельд Николай Валентинович, Варламова Светлана Александровна
Рубрика: Управление в технических системах
Статья в выпуске: 4 т.21, 2021 года.
Бесплатный доступ
Город Березники Пермского края расположен на подработанной шахтной территории. Уже несколько лет в городе наблюдаются активные проседания почвы, которые провоцируют разрушения зданий. Поэтому уже несколько лет ведется мониторинг зданий и сооружений города, позволяющий анализировать степень оседания. Для точного анализа ситуации и прогнозирования используются модели достаточно высокого порядка. Работа посвящена возможности моделирования деформации зданий, связанных с проседанием почвы, в результате горных выработок в городе Березники. Целью исследования является рассмотрение возможностей пакета Robust Control Toolbox для понижения порядка сложности моделей на примере восьмиэтажного здания, входящего в сборник эталонных примеров для редукции моделей линейных динамических систем. Материалы и методы. Представлены типичные шаги для решения задачи редукции модели, описаны команды и инструменты, применяемые для решения этой задачи. Определены параметры модели в пространстве состояний, которая насчитывает 48 состояний, являющихся смещениями или скоростями изменения. С помощью сингулярных значений Ганкеля выбраны состояния, которыми можно пренебречь. Выполнено редуцирование модели с использованием адаптивной границы ошибки. Рассмотрено редуцирование с использованием границы мультипликативной ошибки. Выполнено сравнение результатов редуцирования модели всеми описанными способами, обоснован выбор наилучшего способа редуцирования модели. Результаты. Для всех методов выполнен анализ ошибки аппроксимации. Рассчитана максимальная относительная ошибка. Приведен пример расчета порядка модели для заданной величины ошибки в 5 %. Для такой ошибки порядок модели составил 34 состояния, что меньше исходной модели. Для модели с 34 состояниями величина ошибки составляет менее 1 %. В результате построены АФЧХ исходной и редуцированной модели, а также переходные процессы моделей. Графики в частотной области моделей практически совпадают, что говорит об адекватном описании системы. Заключение. В результате было показано, что возможно снизить размер модели на 14 порядков, цель достигнута.
Сложная модель, редуцирование модели, сингулярные значения, алгоритмы редукции, робастные методы
Короткий адрес: https://sciup.org/147236504
IDR: 147236504 | DOI: 10.14529/ctcr210407
Текст научной статьи Понижение порядка сложных моделей с помощью инструментов Robust Control Toolbox
В городе Березники (Пермский край) на протяжении уже нескольких лет проводится мониторинг зданий и сооружений. Это связано с проседанием почвы, так как практически под всем городом находятся не используемые выработки шахт калийного производства.
Для качественного анализа и прогнозирования необходимы достаточно точные модели. Порядок этих моделей может быть очень большим, что усложняет их дальнейшую обработку и увеличивает затрачиваемые ресурсы времени.
Цель данной статьи – рассмотреть набор инструментов Robust Control Toolbox для понижения порядка сложных моделей на конкретном примере.
Так как у нас нет доступа к данным мониторинга, то в качестве примера возьмем модель движения твердого тела здания больницы университета в городе Лос-Анджелесе. Модель взята из рабочей заметки SLICOT 2002-2 «Сборник эталонных примеров для редукции моделей линейных динамических систем, не зависящих от времени», авторы: Y. Chahlaoui и P.V. Dooren .
В здании восемь этажей, на каждом по три степени свободы – два смещения и одно вращение. В данных имеются отношения вход – выход для любого из этих смещений. В результате получается модель в пространстве состояний с 48 состояниями, где каждое из них представляет собой смещение или скорость его изменения.
В toolbox Robust имеется набор инструментов для работы с моделями высоких порядков. Такие модели, как правило, получаются при точной идентификации технологических объектов, когда мы используем функции взвешивания частоты для формирования отклика разомкнутого контура. Часто мы хотим упростить такие модели для целей моделирования или проектирования систем управления. Как правило, достаточно иметь точную модель вблизи частоты среза. Для моделирования достаточно зафиксировать существенную динамику в диапазоне частот сигналов возбуждения. Это означает, что часто можно найти приближения более низкого порядка для моделей высокого порядка. Модуль Robust Control Toolbox предлагает множество алгоритмов сокращения моделей, которые наилучшим образом соответствуют нашим требованиям и характеристикам модели [1].
Процесс редукции модели
Задача редукции модели обычно включает следующие шаги [2, 3].
-
1. Анализ важных характеристик модели по ее откликам во временной или частотной области с помощью команд step или bode [4, 5].
-
2. Определение порядка сокращения, вычислив сингулярные значения Ганкеля ( hankelsv ) исходной модели и построив по ним соответствующие гистограммы, чтобы определить, какие режимы (состояния) можно отбросить, не жертвуя ключевыми характеристиками. Обычно отбрасывают значения, близкие к нулю [6].
-
3. Выбор алгоритма редукции. В наборе инструментов toolbox Robust доступны следующие функции сокращения: balancmr , bstmr , schurmr , hankelmr и ncfmr . Мы можем легко получить доступ к этим алгоритмам через интерфейс верхнего уровня с помощью команды reduce . В методах используются разные оценки «близости» между исходной и редуцированной моделью. Выбор зависит от используемых параметров команды [3].
-
4. Проверка. На данном этапе мы проверяем наши результаты, сравнивая динамику сокращенной модели с исходной моделью. Нам может потребоваться скорректировать параметры редукции (выбор порядка модели, алгоритма, границ ошибок и т. д.), если результаты неудовлетворительны.
В рассмотренном выше здании восемь этажей, на каждом по три степени свободы – два смещения и одно вращение. В данных имеются отношения вход – выход для любого из этих смещений. В результате получается модель в пространстве состояний с 48 состояниями, где каждое состояние представляет собой смещение или скорость его изменения [7–9]. Применим методы редукции к данной модели.
Откроем модель и построим график логарифмической амплитудно-частотной и фазочастотной характеристики (АФЧХ). Для этого выполним команды:
В результате получим график, приведенный на рис. 1.
Как видно из частотной характеристики модели, основная динамика системы лежит в диапазоне частот от 3 до 50 рад/с. Величина падает как в очень низком, так и в высокочастотном диапазоне. Наша цель – найти модель низкого порядка, которая сохраняет информационное содержание в этом частотном диапазоне до приемлемого уровня точности.
Чтобы понять, какие состояния модели можно безопасно отбросить, вычислим сингулярные значения Ганкеля для нашей модели [10–12]. Для этого выполним команду 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.

Рис. 1. График АФЧХ исходной системы Fig. 1. APАС plot of the original system

Рис. 2. Гистограмма сингулярных значений Ганкеля для исходной модели
Fig. 2. Hankel singular values for the original model
Из графика сингулярных значений Ганкеля видно, что в этой системе есть четыре доминирующих режима. Однако вклад остальных мод все же значителен. Мы проведем линию через 10 состояний и отбросим оставшиеся, чтобы найти сокращенную модель 10-го порядка ( Gr ), которая наилучшим образом приближает исходную систему ( G ).
Функция reduce является шлюзом ко всем процедурам редукции моделей, доступным в наборе инструментов [13, 14]. Мы будем использовать параметр усечения баланса квадратного корня по умолчанию ' balancmr ' для функции reduce как первый шаг. В этом методе для уменьшения порядка используется «аддитивная» граница ошибки, что означает, что он пытается поддерживать одинаково малую абсолютную ошибку аппроксимации по частотам.
Данный метод используется в функции по умолчанию.
Вычислим сокращенную модель 10-го порядка. Для этого выполним команду
[Gr_add,info_add] = reduce(G,10);
Для сравнения исходной ( G ) и редуцированной модели ( 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.

Рис. 3. АФЧХ исходной и редуцируемой модели Fig. 3. APFC original and reduced model
Как видно из графиков на рис. 3, упрощенная модель довольно хорошо улавливает резонансы ниже 30 рад/с, но согласование в области низких частот (< 2 рад/с) плохое. Кроме того, упрощенная модель не полностью отражает динамику в диапазоне частот 30–50 рад/с. Возможное объяснение больших ошибок на низких частотах – относительно низкий коэффициент усиления модели на этих частотах. Следовательно, даже большие ошибки на этих частотах мало влияют на общую ошибку [15, 16].
Чтобы обойти эту проблему, мы можем попробовать метод мультипликативной ошибки, такой как bstmr . Этот алгоритм делает упор на относительные ошибки, а не на абсолютные [17, 18]. Поскольку относительные сравнения не работают, когда коэффициенты усиления близки к нулю, нам нужно добавить порог минимального усиления, например, добавив сквозное усиление.
Под сквозным усилением понимается непосредственная связь входа с выходом. Такая связь существует, если значение матрицы d пространства состояний не равно нулю. Предполагая, что нас не беспокоят ошибки при усилении ниже –100 дБ, мы можем установить коэффициент усиления сквозной связи, равным 1e-5.
Для этого выполним команду
GG = G; GG.d = 1e-5;
Снова вычислим сингулярные значения Ганкеля и построим гистограмму:
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.
Multiplicative-Error Singular Values of the Model (GG)

Рис. 4. Гистограмма сингулярных значений Ганкеля для редуцированной модели
Fig. 4. Histogram of Hankel singular values for the reduced model
Из диаграммы видно, что нам необходима модель 26-го порядка, но для сравнения с предыдущим результатом давайте придерживаться редукции на 10-й порядок. Используем опцию алгоритма bstmr для сокращения модели:
[Gr_mult,info_mult] = reduce(GG,10,'algorithm','bst');
Построим графики АФЧХ для исходной ( G ) и редуцированной ( Gr_add и 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.

Рис. 5. Графики АФЧХ для моделей G, Gr_add, Gr_mult
Fig. 5. APFC plots for models 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.

Рис. 6. Переходные процессы для моделей G, Gr_add, Gr_mult Fig. 6. Transients for models G, Gr_add, Gr_mult
Проверка результатов и корректировка модели
Все алгоритмы обеспечивают границы ошибки аппроксимации. Для методов аддитивной ошибки, таких как balancmr , ошибка аппроксимации измеряется пиковым (максимальным) коэффициентом усиления модели ошибки | G-Greduced на всех частотах. Максимальное значение можно получить, если вычислить бесконечную (Чебышевскую) норму вектора с помощью функции norm(A,inf). Теоретическая граница хранится в поле Error Bound структуры INFO, возвращаемой функцией reduce .
Для анализа выполним команды:
E1 = norm(G-Gr_add,inf)
E2 = info_add.ErrorBound
[E3,f] = norm(GG\(GG-Gr_mult),inf)
E4 = info_mult.ErrorBound
В результате получим:
E 1 = 6.004322380100927e–04 – абсолютная ошибка модели Gr_add ;
E 2 = 0.004718864240517 – теоретическое значение границы модели Gr_add ;
E 3 = 0.594921503759114 – максимальная относительная ошибка модели Gr_mult и соответствующая ей частота f = 17.2;
E 4 = 4.732954442311137e+02 – теоретическое значение границы модели Gr_multi .
Построим график относительных ошибок модели Gr_mult в частотной области. Для этого выполним команды:
figure;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.

Рис. 7. График относительной ошибки модели Gr_mult Fig. 7. Graph of the relative error of the Gr_mult model
Из приведенного выше графика и сделанных расчетов можно увидеть, что максимальная относительная ошибка достигает 59 % при частоте 17,2 рад/с. В нашем случае это больше допустимого.
Функция reduce позволяет выбрать наименьший порядок, совместимый с желаемым уровнем точности. Чтобы повысить точность модели Gr_mult , нам придется увеличить порядок.
Зададимся максимальной ошибкой, равной 5 %, и выполним редуцирование с заданной точностью.
[Gred,info] = reduce(GG,'ErrorType','mult','MaxError',0.05);
Определим порядок полученного объекта:
P=size(Gred.a)
В результате получим P = 34. Следовательно, минимальный порядок модели, обеспечивающий 5%-ную ошибку, составляет 34.
Вычислим относительную ошибку и теоретическую границу модели Gred :
E5 = norm(GG\(GG-Gred),inf)
E6 = info.ErrorBound
В результате получим
E 5 = 0.0068 – максимальная относительная ошибка модели Gred .
E 6 = 0.0422 – теоретическое значение границы модели Gred .
Повышение точности достигнуто за счет увеличения порядка модели (с 10 до 34). Обратите внимание, что фактическая максимальная ошибка составляет 0,6 %, что намного меньше целевого показателя в 5 %. Это несоответствие является результатом функции bstmr , которая использует границу ошибки, а не фактическую ошибку, чтобы выбрать необходимый порядок.
И в заключение построим АФЧХ и переходные процессы исходной и редуцированной модели:
figure; bode(G,Gred,{1e-2,1e4}), grid on legend('G - 48-state original','Gred - 34-state reduced')
figure; step(G,'b',Gred,'r--',15),grid %step response until 15 seconds legend('G: 48-state original ','Gred: 34-state (bstmr)')
text(4,-5e-4,'Maximum relative error <0.05').

Рис. 8. АФЧХ исходной и редуцированной модели Fig. 8. APFC of the original and reduced model

Рис. 9. Переходные процессы исходной и редуцированной модели Fig. 9. Transient processes of the original and reduced model
Анализируя рис. 8 и 9, мы видим, что графики в частотной области, где проявляется основная динамика объекта, практически совпадают. Ошибка аппроксимации составляет менее 5 %. В результате для адекватного описания системы мы смогли понизить ее размер на 14 порядков, что положительным образом должно сказаться на дальнейших расчетах и работе с данной моделью.
Список литературы Понижение порядка сложных моделей с помощью инструментов Robust Control Toolbox
- Перельмутер, В.М. Пакеты расширения MATLAB. Control System Toolbox и Robust Control Toolbox: практ. пособие /В.М. Перельмутер. - М.: СОЛОН-ПРЕСС, 2008. - 224 с.
- Романова, И.К. Современные методы редукции систем и их применение к задачам анализа и синтеза систем управления / И.К. Романова // Вестник МГТУ им. Н.Э. Баумана. Спец. вып. «Специальная робототехника и мехатроника». - 2011. - С. 142-152.
- Конструирование объекта управления. Часть I. Проблема редуцирования модели, размещения управляющих органов, измерительных устройств и оценки потенциальной робастности / Д.С. Бирюков, Н.А. Дударенко, О.В. Слита, А.В. Ушаков //Мехатроника, автоматизация, управление. - 2013. - № 6 - С. 2-6.
- Vidyasagar, M. Control System Synthesis: A Factorization Approach /М. Vidyasagar. - London: The MIT Press, 1985.
- An information theortic approach to model reduction based on frequency-domain cross-gramian information / J. Fu, C. Zhong, Y. Ding, J. Zhou // Proceedings of 8th World Congress on Intelligent Control and Automation. - 2010. - P. 3679-3683.
- Качественный анализ динамики позиционного регулирования температуры процесса восстановления титана / Ю.П. Кирин, А.В. Затонский, В.Ф. Беккер, Н.В. Бильфельд // Приборы и системы. Управление, контроль, диагностика. - 2008. - № 10. - С. 54-56.
- Minimum information loss method on cross-gramian matrix for model reduction / J. Fu, H. Zhang, Y. Sun, C. Zhong // 129 proceedings of 7th World Congress on Intelligent Control and Automation. -2008. - P. 7339-7343.
- Дядик, В.Ф. Теория автоматического управления /В.Ф. Дядик, С.А. Байдали, Н.С. Крини-цын. - Томск: Изд-во Томского политехн. ун-та, 2011. - 196 с.
- Бильфельд, Н.В. Многокритериальное исследование систем управления / Н.В. Бильфельд. -Пермь: Изд-во Перм. нац. исслед. политехн. ун-та, 2015. - 440 с.
- Matrix Eigensystem Routines - EISPACK Guide / B.T. Smith, J.M. Boyle, J.J. Dongarraet et al. //Lecture Notes in Computer Science. - Berlin, 1976. - Vol. 6. DOI: 10.1007/3-540-07546-1
- Matrix Eigensystem Routines - EISPACK Guide Extension / B.S. Garbow, J.M. Boyle, J.J. Dongarr, C.B. Moler // Lecture Notes in Computer Science. - Berlin, 1977. - Vol. 51. DOI: 10.1007/3-540-08254-9
- Поршнев, С.В. Об особенностях сингулярных чисел и сингулярных векторов выборочной корреляционной матрицы в методе SSA / С.В. Поршнев, Р. Фуад // Научно-технический вестник Поволжья. - 2012. - № 3. - С. 146-150.
- Freitas, F. Gramian-based reduction method applied to large sparse power system descriptor models / F. Freitas, J. Rommes, N. Martins // Proceedings of IEEE Power and Energy Society General Meeting. - 2009. - P. 1. DOI: 10.1109/PES.2009.5275568
- Симахин, В.А. Робастные непараметрические оценки / В.А. Симахин. - Саарбрюкен: LAP LAMBERT Academic Publishing GmbH & Co, 2011. - 292 с.
- Никифоров, В.О. Управление в условиях неопределенности: чувствительность, адаптация, робастность /В.О. Никифоров, А.В. Ушаков. - СПб.: СПбГИТМО(ТУ), 2002.
- Бильфельд, Н.В. Современные средства моделирования динамики системы автоматизации /Н.В. Бильфельд, Ю.И. Володина. - Пермь: Изд-во Перм. нац. исслед. политехн. ун-та, 2020. -143 с.
- Влияние нестационарности объекта управления на параметры установившихся автоколебаний / М.Н. Ерыпалова, В.Ф. Беккер, А.В. Затонский, Ю.П. Кирин // Известия высших учебных заведений. Поволжский регион. Технические науки. - 2008 - № 4 (8) - С. 50-57.
- Safonov, M.G. Model Reduction for Robust Control: A Schur Relative Error Method / M.G. Safonov, R.Y. Chiang // International J. of Adaptive Control and Signal Processing. - 1988. -Vol. 2. - P. 259-272.