Имаев Д.Х. Синтез систем управления в среде MATLAB

Подождите немного. Документ загружается.

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

08.0

max0



рад приведены на рис. 3.2.

Система регулирования, синтезированная операторным методом, может

оказаться весьма чувствительной к изменениям параметров объекта и неточной

реализации алгоритма. Для иллюстрации введем коэффициенты передаточной

функции регулятора с клавиатуры, сохраняя четыре значащие цифры

regulator=tf([195.4 1504 -0.6521 -0.06398],[1 16.82 -331.5 -2510])

Передаточную функцию замкнутой системы получим по команде:

>> sysc=feedback(plant,regulator)

2.5 s^5 + 42.05 s^4 - 926.8 s^3 - 7923 s^2 + 32487 s + 245980

-----------------------------------------------------------------------------

s^7 + 16.82 s^6 + 98.2 s^5 + 261 s^4 + 341.4 s^3 + 195.8 s^2 + 63.91 s + 6.27

Видно, что коэффициент знаменателя (характеристического полинома замкнутой

системы) при второй степени равен 195.8 вместо желаемого значения 221. Такое

отличие объясняется необходимостью выполнения операций по вычитанию

больших по модулю примерно равных величин.

0 2 4 6 8 10 12 14 16 18 20

-2

-1

theta

Рис. 3.2. Процессы в системе при максимальном отклонении маятника

Система в данном примере устойчива, однако, ее собственные значения:

-8.0078

-4.1704

-1.9346 + 1.0956i

-1.9346 - 1.0956i

-0.3100 + 0.3910i

-0.3100 - 0.3910i

-0.1526

отличаются от заданных.

Однако это не влияет на размер области устойчивости и характер процесса

при отклонении маятника 0.08 рад (см. рис. 3.2).

3.2. Синтез системы стабилизации маятника

на каретке методами пространства состояний

Задачу будем решать на основе линеаризованной модели объекта, т. е.

система автоматической стабилизации должна будет обеспечивать устойчивость

при малых отклонениях состояния от положения равновесия.

3.2.1 Регулятор состояния

Пусть имеется информация обо всех переменных состояния объекта

BAv



. (3.3)

Алгоритм регулятора состояния запишется так:

   

44332211

4321

vkvkvkvk

kkkkF 













 Kv

. (3.4)

Уравнение замкнутой системы в стандартной форме пространства

состояний получим, если из уравнений объекта (3.3) и регулятора состояния (3.4)

исключим переменную

vBKA

)( 

где

)( BKA 

— матрица замкнутой системы. Отметим, что замкнутая система

получилась автономной — она не имеет входов.

Задача имеет решение, если объект управляем. Для анализа управляемости

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

управляемости

 

BABAABB

>> rank(ctrb(A,B))

ans =

Матрица управляемости имеет полный ранг, что свидетельствует о полной

управляемости объекта. Действуя с силой

на каретку, можно стабилизировать

верхнее положение маятника и привести каретку в исходное положение.

3.2.2. Синтез регулятора методом размещения

собственных значений

Задача синтеза заключается в определении матрицы регулятора состояния

K из условия желаемого расположения собственных значений матрицы системы

(корней характеристического полинома системы).

Существует произвол в выборе желаемых собственных значений системы.

Чем дальше находится от мнимой оси собственное значение, тем быстрее

затухают процессы. Однако требование большего быстродействия означает

необходимость приложения чрезмерных усилий на каретку. Назначим

собственные значения системы с ориентацией на собственные значения объекта,

которые вычислены ранее. Введем желаемые значения

>> p=[-1 -2 -4 -7]';

Матрицу регулятора, обеспечивающего желаемое расположение

собственных значений, найдем по команде:

>> K=place(A,B,p)

-0.5714 -12.3229 -1.6704 -1.0816

Проведем анализ устойчивости системы:

>> eig(A-B*K)

ans =

-7.0000

-4.0000

-2.0000

-1.0000

Замкнутая система имеет желаемое расположение собственных значений.

Регулятор состояния формирует управляющие воздействия на основе

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

измеряется только положение каретки х. Для вычисления остальных переменных

состояния используют так называемый наблюдатель состояния [2, 24]. Если

наблюдатель устойчив, то процессы затухают, в результате чего состояние

модели

приближается к состоянию объекта

Задача синтеза наблюдателя сводится к поиску матрицы L, и имеет

решение, если объект наблюдаем полностью. Наблюдаемость можно проверить

по критерию Калмана, который сводится к проверке ранга матрицы

наблюдаемости















>> rank(obsv(A,C))

ans =

Матрица наблюдаемости имеет полный ранг, следовательно, объект наблюдаем

полностью.

Для синтеза наблюдателя воспользуемся методом размещения собственных

значений – назначим собственные значения наблюдателя несколько дальше от

мнимой оси в левой полуплоскости:

>> po=5*p

po =

-5

-10

-20

-35

Матрица наблюдателя L вычисляется по команде:

>> L=place(A',C',po)'

L =

1.0e+004 *

0.0070

-0.3544

-2.6748

0.1634

Динамический регулятор получим, если объединить регулятор состояния и

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

Поэтому рекомендуется синтезировать редуцированный наблюдатель,

учитывающий, что одна из переменных состояния — положение каретки

—

измеряется непосредственно.

Уравнения динамического регулятора имеют вид:

rrr

BvA



;

rrr

DvC 

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

>> [Ar,Br,Cr,Dr]=reg(A,B,C,D,K,L);

>> regulator=ss(Ar,Br,Cr,Dr);

Передаточная функция регулятора выглядит так

>> [numr,denr]=ss2tf(Ar,Br,Cr,Dr);

>> regulator=tf(numr,denr)

8.655e004 s^3 + 6.58e005 s^2 - 4.543e004 s - 2e004

--------------------------------------------------

s^4 + 84 s^3 + 2677 s^2 - 1.716e005 s - 1.16e006

Знаменатель ПФ — характеристический полином динамического регулятора —

имеет отрицательный коэффициент. Проверим устойчивость регулятора:

>> eig(Ar)

ans =

-55.0088 +52.1090i

-55.0088 -52.1090i

32.2791

-6.2615

Регулятор неустойчив — имеется положительное собственное значение.

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

который заключается в превышении необходимого порядка. Поскольку одна из

переменных состояния, а именно, положение каретки измеряется

непосредственно, то нет необходимости в ее восстановлении. Может быть

синтезирован наблюдатель пониженного порядка [2].

3.2.3. Анализ системы с динамическим регулятором

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

замкнутой системы:

>> sysc=feedback(plant,regulator);

и вычислим собственные значения системы:

>> eig(sysc)

ans =

-35.0000

-20.0000

-10.0000

-1.0000

-2.0000

-4.0000

-5.0000

-7.0000

Замкнутая система с динамическим регулятором устойчива и имеет желаемые

собственные значения.

Анализ линейной модели показывает, что существует область устойчивости

положения равновесия — имеет место устойчивость “в малом”. Однако для

моделей рассматриваемого типа нет аналитических методов определения

размеров области притяжения. Единственным способом важной для практики

оценки является многократное компьютерное моделирование.

Прежде всего, подключим динамический регулятор к нелинейному объекту,

как это показано на рис. 3.3.

Для отредактированной нелинейной системы (см. рис. 3.3) целесообразен

предварительный анализ устойчивости “в малом” положения равновесия. Для

анализа существования области притяжения положения равновесия в

соответствии с первым методом Ляпунова проведём линеаризацию замкнутой

системы и вычислим собственные значения:

>> [Ac,Bc,Cc,Dc]=linmod2('cartpole_closed');

>> eig(Ac)

ans =

-35.0000

-20.0000

-10.0000

-1.0000

-2.0000

-4.0000

-5.0000

-7.0000

Так как получены те же собственные значения, то область устойчивости

существует.

Out1

Scope

regulator

LTI System

Integrator3

Integrator2

Integrator1

Integrator

f(u)

Fcn1

f(u)

Fcn

In1

Рис. 3.3. Система «нелинейный объект + динамический регулятор»

Оценим размеры области притяжения положения равновесия путём

многократных компьютерных экспериментов. Для примера на рис. 3.4.

0 1 2 3 4 5 6 7

-0.8

-0.6

-0.4

-0.2

0.2

0.4

0.6

0.8

cart

pole

Рис. 3.4. Поведение системы при максимальном отклонении маятника на 0.34

рад

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

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

объект.

3.2.4. Аналитическое конструирование регуляторов

Требования устойчивости и качества процессов можно описывать в неявной

форме как экстремали тех или иных функционалов. Наиболее часто, в силу

относительной простоты вычислений при достаточной физичности, применяют

интегральные квадратичные функционалы [24]. В случае, когда объект описан в

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

записывается в виде

I ru dt

 





( )v Qv

(3.5)

где: vk— вектор состояния; uk— скалярное управление; Qk— неотрицательно-

определенная весовая матрица; rk— весовой коэффициент. Безусловная

экстремаль

( )t

функционала (3.5) отвечает желаемому поведению и зависит от

выбора весовых коэффициентов.

Выберем весовую матрицу единичной; коэффициент при управлении также

положим равным единице:

>> Q=eye(4)

1 0 0 0

0 1 0 0

0 0 1 0

0 0 0 1

>> r=1;

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

>> [K,S,E]=lqr(A,B,Q,r)

K =

-1.0000 -16.7054 -2.5344 -1.7914

E =

-14.2914

-4.1053

-1.2343 + 0.3832i

-1.2343 - 0.3832i

Программа выдает коэффициенты регулятора состояния К и собственные

значения замкнутой системы Е.

Примем тот же наблюдатель, что и ранее. Матрицы динамического

регулятора, объединяющего регулятор состояния и наблюдатель, получим по

выполнению команды

>> [Ar,Br,Cr,Dr]=reg(A,B,C,D,K,L);

полученный регулятор назовем так

>> lqreg=ss(Ar,Br,Cr,Dr);

Имитация системы «нелинейный объект + линейный регулятор»

показывает, что процессы (рис. 3.5) затухают при максимальном отклонении

маятника 0.38 рад.

0 1 2 3 4 5 6 7

-0.8

-0.6

-0.4

-0.2

0.2

0.4

0.6

0.8

cart

pole

Рис. 3.5. Процессы в системе линейно-квадратичным регулятором

при максимальном начальном отклонении маятника 0.38 рад

3.3. Синтез системы стабилизации частотным методом

Частотные методы позволяют синтезировать линейные системы с одним

входом и одним выходом (типа SISO) с учетом естественной динамики объектов

или корректируемых систем.

Целесообразно воспользоваться весьма развитыми средствами

MATLAB/SISO Design Tool, позволяющими автоматизировать синтез и

коррекцию систем в частотной области. Однако здесь для наглядности процедуры

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

MATLAB/Control System Toolbox.

Ранее было выявлено, что по паре вход-выход

f

передаточная функция

объекта вырождается до второго порядка -10/(s^2 - 58.8), так как исходная

передаточная функция имеет два нулевых нуля и полюса. Напомним, что это

означает неполную наблюдаемость объекта по этому каналу ― с помощью

обратной связи от измерителя углового положения маятника до привода каретки

нельзя стабилизировать каретку.

Частичная наблюдаемость по одному из каналов может быть и

положительным свойством объекта. Здесь это свойство будет использовано для

декомпозиции синтеза и децентрализации управления.

Вначале найдем обратную связь, стабилизирующую положение маятника ―

от выхода объекта



до входа

. Другими словами, ищется обратная связь,

перемещающая неустойчивый полюс 7.668.

На рис. 3.6. приведена логарифмическая амплитудно-частотная

характеристика канала

f

(кривая 1).

-1

-150

-100

-50

Magnitude (dB)

Bode Diagram

Frequency (rad/sec)

Рис. 3.6 . Логарифмические амплитудно-частотные характеристики

разомкнутого контура стабилизации маятника

Известно, что обратная связь перемещает корни характеристического

полинома системы в разной степени в зависимости от усиления контура на

частотах, равных модулям корней. Например, если усиление контура менее -20

дБ, то корни, модули которых принадлежат этому диапазону частот,

перемещаются мало. Полюсы передаточной функции объекта вообще

неподвижны, если их компенсируют равные им нули, т. е. на комплексной

частоте полюсов фактическое усиление контура равно нулю.

Как следует из рис. 3.6. (кривая 1), усиление объекта по рассматриваемому

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

частотах перемещаемых полюсов усиление превышало 20 дБ. Для повышения

степени подвижности неустойчивого полюса -7.6681 повысим усиление контура

примерно в 200 раз или на 52 дБ (кривая 2). Это расширяет полосу частот контура

и повышает быстродействие системы.

Далее сформируем типовой желаемый вид ЛАЧХ в среднечастотной

области. Как известно, типовая асимптотическая ЛАЧХ в области средних частот

имеет наклон -20 дБ/дек. Чем длиннее этот отрезок, тем больше запас

устойчивости системы.

Асимптотическая ЛАЧХ нескорректированного контура состоит из двух

асимптот ― низкочастотной с наклоном 0 дБ/дек и высокочастотной с наклоном -

40 дБ/дек. Для придания желаемой ЛАЧХ контура типового вида в области

средних частот введем в контур действительный нуль -20 и действительный

полюс -1000, т. е. последовательную коррекцию с передаточной функцией

>> corr=tf([1/20 1],[1/1000 1])

0.05 s + 1

----------- ,

(3.6)

0.001 s + 1

что дает отрезок асимптоты с наклоном -20 дБ/дек в окрестности частоты среза

100

ср



рад/с (кривая 3) на рис. 3.6. При этом крайние частоты отрезка

асимптоты отличаются в 50 раз, что должно обеспечивать большой запас

устойчивости замкнутой системы. Заметим также, что такая активная коррекция

еще более расширяет полосу частот контура стабилизации маятника, т. е.

повышает быстродействие системы.

Ввиду относительной малости постоянной времени знаменателя

передаточной функции (3.6) можно полагать, что для стабилизации маятника

принят ПД-регулятор с передаточной функцией 200(0.05 s + 1)/(0.001 s + 1).

Создадим simulink-модель системы с замкнутым контуром стабилизации

маятника с именем 'linear_pend_1' (рис. 3.7).

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

преобразуем модель к форме пространства состояний

>> [A1,B1,C1,D1]=linmod2('linear_pend_1');

и вычислим собственные значения матрицы

>> eig(A1)

-890.1788

-83.7981

-26.0231

Theta

Out2

Out1

1/20s+1

1/1000s+1

Transfer Fcn

Integrator

200

1/M

m*g

1/l/M

(M+m)*g

In1

Рис. 3.7. Структурная схема линейной системы стабилизации маятника

Цель достигнута ― собственные значения подсистемы являются

отрицательными действительными числами (это следствие большого запаса

устойчивости), и находятся достаточно далеко от границы устойчивости (в силу

большого значения частоты среза). Естественно, что два нулевых собственных

значения неуправляемой части объекта остались на месте.

Приступим ко второму этапу синтеза ― к стабилизации положения каретки.

На этом этапе объектом синтеза оказывается подсистема, синтезированная на