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

Собственно сабж



Подвижная масса на роликах находится на безмассовой платформе.
Картинко:


Планарную неподвижную систему координат сонаправим с $u$ и $y$.
Поскольку единственная масса в системе -- $m$, то для нее можно записать второй закон Ньютона:
$$ m a = \sum_i F_i $$
где $a$ -- ускорение, $F_i$ -- силы, действующие на подвижный груз.
$$ a = \frac{d^2 y}{d t^2} $$
$$ m \frac{d^2 y}{d t^2} = - k (y - u)$$
$$ \frac{d^2 y}{d t^2} = - k (y - u) /m$$
Теперь представим систему в пространстве состояний. Проблема в том, что уравнения объекта управления должны быть записаны в т.н. нормальной форме Коши -- т.е. системой уравнений первого порядка, а в данном случае уравнение второго порядка. Нужно просто ввести фиктивную переменную $z$, что $\dot y = z$. Тогда представление в пространстве состояний:

$$\dot y = z$$
$$\dot z = - k (y - u) /m$$

Для моделирования такой системы в Simulink можно воспользоваться блоком State Space (Simulink/Continous). Для использования этого блока необходимо представить уравнения объекта управления в матричном виде. И, кроме того, блок State Space требует, чтобы число входов, число выходов и число переменных состояний системы совпадало, однако в рассматриваемой системе один вход и выход, но два состояния -- решается с помощью фиктивных входов и выходов.

Вектором состояний является $X = (y, z)^T$, вход $U = (u, .)^T$, выход $Y = (y, .)^T$ (точка означает фиктивную переменную). Тогда в матричном виде динамика системы запишется:

$$\begin{pmatrix} \dot y \cr \dot z \end{pmatrix} = \begin{pmatrix} 0 && 1 \cr -k/m && 0 \end{pmatrix} \begin{pmatrix} y \cr z \end{pmatrix} + \begin{pmatrix} 0 && 0 \cr k/m && 0 \end{pmatrix}

\begin{pmatrix} u \cr . \end{pmatrix}$$

Или в матричном виде окончательно:
$$\dot X = A X + B U$$
$$Y = C X + D U$$

где

$$A = \begin{pmatrix} 0 && 1 \cr -k/m && 0 \end{pmatrix}$$
$$B = \begin{pmatrix} 0 && 0 \cr k/m && 0 \end{pmatrix}$$
$$C = \begin{pmatrix} 1 && 0 \cr 0 && 0 \end{pmatrix}$$
$$D = 0 = \begin{pmatrix} 0 && 0 \cr 0 && 0 \end{pmatrix}$$

Положив для определенности $m=k=1$, запишем матрицы в свойства блока Simulink/Continous:


Модель системы без обратной связи:


Поскольку блок State Space принимает на вход вектор и выдает также вектор, то для подключения в модели необходимо использовать элементы Mux и Demux для входа и выхода соответственно.

Видно, что после подачи единичного импульса на вход, в системе начинаются незатухающие колебания, что согласуется с физическим смыслом подвижного груза на пружине без учета трения и рассеивания энергии (математический маятник):



Осуществим автоматизированный синтез регулятора с помощью функции lqr в пакете Control System Toolbox.

Во-первых, для управления системой необходимо, чтобы все состояния были наблюдаемы. Для этого выведем вторым выходом переменную состояния $z$ и тогда новая матрица $C$ запишется в виде единичной диагональной:

$$C = \begin{pmatrix} 1 && 0 \cr 0 && 1 \end{pmatrix}$$

Функция lqr синтезирует с использованием численной оптимизации стабилизирующий регулятор таким образом, чтобы минимизировать квадратичную целевую функцию:

$$ J = \int_0^{\infty} x^T(t)Qx(t) + u^T(t)Ru(t) dt \to min$$

Матрицы $Q$ и $R$ предназначены просто для выбора весов (важности) минимизации одного критерия перед другим. Для начала можно выбрать единичные диагональные матрицы $Q = R = C = 1$,

Следующий код на языке MATLAB запускает функцию lqr

% ввод системы в пространстве состояний
n = 2; % число состояний
A = [0 1; -1 0]; 
B = [0 0; 1 0];
C = [1 0; 0 1];
D = [0 0; 0 0];

% синтез регулятора
Q = diag(ones(1,n)); % веса выходов
R = diag(ones(1,n)); % веса входов
K = lqr(A,B,Q,R) % синтез регулятора

После выполнения последней строчки кода получаем в консоли MATLAB

K =

    0.4142    1.3522
         0         0

Промоделируем систему с обратной связью. Для этого добавим в систему блок Gain со свойствами:


Система с обратной связью:


Реакция на единичный импульс:


Заметно, что установившееся значение выхода системы не равно входному (0.7 вместо 1). Эта ошибка является пропорциональной и означает, что коэффициент усиления системы с обратной связью отличается от 1. Включив усилитель Gain с коэффициентом 1.412, получаем в результате систему с единичным коэффициентом усиления:


После проделанного система повторяет значение на входе в установившемся режиме.

Реакция на единичный импульс амплитудой 2:


Дополнительные ссылки:
http://matlab.exponenta.ru/controlsystem/book1/5.php
http://www.engin.umich.edu/group/ctm/examples/pend/invSS.html

Файлы модели с приведенным решением для Simulink:
- без обратной связи
- с обратной связью