Линеаризация динамической системы¶
В этом демонстрационном примере мы идентифицируем LTI-систему (линейную и стационарную) по входам-выходам динамической модели. Это нам позволит линеаризовать выбранный отрезок модели, а затем работать с ним как с обычной линейной стационарной системой.
Описание модели¶
Объектом линеаризации для нас будет служить модель, состоящая из двух RC фильтров, через которые будут проходить прямоугольные импульсы с периодом 0.01 с.
Каждый фильтр состоит из резистора 10 Ом и конденсатора 1e-04 Ф. Входной сигнал u поступает из блока Pulse Generator, а выходной сигнал y является напряжением на последнем конденсаторе каскада, измеренным вольтметром Voltage Sensor.
gr()
modelName = "rc_circuit"
if !(modelName in [m.name for m in engee.get_all_models()]) engee.load( "$(@__DIR__)/$(modelName).engee"); end;
data = engee.run( modelName )
Если мы построим графики переменных u и y, то сможем оценить, насколько цепочка фильтров видоизменяет входной сигнал. Именно это поведение мы попробуем приблизительно описать при помощи линейной модели.
t,u,y = data["u"].time, data["u"].value, data["y"].value
plot( t, [u y], label=["u" "y"] )
Можно ожидать, что апериодического звена может не хватить для достаточно точного приближения поведения системы.
Процесс линеаризации¶
Мы выполним линеаризацию при помощи функции subspaceid из библиотеки ControlSystemIdentification:
using ControlSystems, ControlSystemIdentification
По умолчанию subspaceid возвращает систему в пространстве состояний следующего вида:
$$\begin{aligned} x^+ &= Ax + Bu + Ke\\ y &= Cx + Du + e \end{aligned}$$
В пакете доступно больше десяти других функций для линеаризации систем, но именно эту рекомендуется использовать в первую очередь, если система задана в разомкнутом виде. Для систем в замкнутой форме рекомендуется newpem.
# Получим частоту дискретизации системы из вектора времени
Ts = t[2] - t[1]
# Проведем линеаризацию системы
sys = subspaceid( iddata( y, u, Ts ));
Для наглядности выведем эту систему в виде передаточной функции:
using SymbolicControlSystems
Sym( tf(sys) )
Та же библиотека позволит нам получить формулу системы в синтаксисе $\LaTeX$.
latextf(tf(sys))
График реакции на ступеньку позволит нам оценить качество линеаризации (мы специально выставили амплитуду генератора импульсов равную 1, в ином случае нужно было бы ее нормализовать).
plot( t, y, label="Исходный сигнал" )
plot!( step( sys, 0.03 ), label="Реакция системы на ступеньку" )
Видно, что реакция на единичную ступеньку весьма точно повторяет результаты симуляции – реакцию исходной системы на прямоугольные импульсы. При желании мы могли бы вторым аргументом subspaceid ограничить степень линейной системы и получить более простую, менее точную аппроксимацию.
Изучение устойчивости системы¶
Конечно, наблюдаемая система устойчива. Но в случае более сложной системы мы бы могли делать такие констатации после анализа отклика, АФЧХ или годографов.
Простой график линейной амплитудно-фазово-частотной характеристики системы можно построить при помощи короткой функции:
bodeplot( sys )
Но следующий код позволяет иметь гораздо больше контроля за характером отображения графиков:
mag, phase, w = bode( sys )
plot( [w./(2pi)],
[20*log10.( mag[:]) phase[:]], xscale=:log10, xticks=exp10.(1:8),
linecolor=[1 2], layout=(2,1))
vline!( [1/(2Ts) 1/(2Ts)] )
plot!( title=["ЛАФЧХ системы" ""], xlabel=["" "Частота (Гц)"] )
plot!( ylabel=["Амплитуда (дБ)" "Фаза (град)"], leg=:false )
Красной вертикальной линией на них отмечена частота, равная половине частоты дискретизации. График фазы ограничен значением в 360 градусов, но дальше значения ЛАФЧХ резко падают.
Заключение¶
Этот способ линеаризации позволит построить LTI-модель на основе широкого диапазона динамических систем или просто табличных данных. Полученную упрощенную динамическую систему можно поместить в блок Transfer Function и использовать как суррогатную модель, которая быстро рассчитывается и, например, позволяет в первом приближении изучить устойсивость динамической системы.