Практика

При работе с реальными объектами управления, как правило, сначала производится построение уравнений динамики объекта управления, синтез регулятора и последующий анализ в среде математического моделирования. Только после этого полученный закон управления и коэффициенты обратной связи, удовлетворяющие желаемым требованиям качества системы, переносят на аппаратную платформу. Такой порядок особенно важен для устройств, которые при неверной настройке могут выйти из строя или причинить вред человеку. К таким устройствам можно отнести летательные аппараты, автомобили или шагающие роботы. Именно поэтому сначала производится моделирование и уже потом тестирование на реальном объекте. Ход практической работы будет строиться по этому же принципу.

Предупреждение

Напомним, что перед выполнением дальнейшей части работы вам следует установить и настроить ПО по инструкции из предыдущего раздела (Windows, Ubuntu).

Математическое моделирование

Откройте файл Python/model.ipynb. В нём, в качестве примера, реализована модель коллекторного двигателя. Его динамика описывается системой уравнений:

$$\begin{aligned} \ddot{\theta} &= -\frac{b}{J}\dot{\theta} + \frac{K}{J}i, \\ \dot{i} &= -\frac{K}{L}\dot{\theta} - \frac{R}{L}i + \frac{1}{L}V. \end{aligned} \tag{22}$$

где $\theta$ — угол поворота выходного вала двигателя, $i$ — ток в обмотках статора, $V$ — напряжение, $b$ — коэффициент трения, $J$ — момент инерции ротора, $K$ — постоянная противо-ЭДС, $L$ — индуктивность обмоток, $R$ — сопротивление обмоток. В матричной форме эти уравнения выглядят как:

$$\frac{d}{dt} \begin{bmatrix} \theta \\ \dot{\theta} \\ i \end{bmatrix} = \begin{bmatrix} 0 & 1 & 0 \\ 0 & -\frac{b}{J} & \frac{K}{J} \\ 0 & -\frac{K}{L} & -\frac{R}{L} \end{bmatrix} \begin{bmatrix} \theta \\ \dot{\theta} \\ i \end{bmatrix} + \begin{bmatrix} 0 \\ 0 \\ \frac{1}{L} \end{bmatrix} V. \tag{23}$$

Параметры маятника

В блоке Pendulum Physical Parameters приведены все необходимые физические параметры исследуемого объекта управления. Обычно эти параметры берутся из документации к составным устройствам или получаются экспериментальным путём. В данном случае в блоке приведены параметры коллекторного двигателя. Отредактируйте этот блок и впишите все необходимые параметры маятника, используя табл. 3. Используйте те же названия переменных, что и в графе Обозначение.

Таблица 3. Физические параметры маятника

Параметр	Обозначение	Ед. изм.	Значение	Описание
$l_1$	L1	м	0.044	Расстояние от точки $O$ до центра масс стержня
$l_2$	L2	м	0.088	Расстояние от точки $O$ до маховика
$R$	R	м	0.0655	Радиус маховика
$m_1$	m1	кг	0.1975	Масса стержня
$m_2$	m2	кг	0.2545	Масса маховика
$I_1$	I1	кг·м²	m1*L1**2/12	Момент инерции стержня
$I_2$	I2	кг·м²	m2*R**2/2	Момент инерции маховика
$g$	g	м/с²	9.81	Ускорение свободного падения
$\mu_p$	b_p	Н·м·с/рад	0.007	Коэффициент вязкого трения в шарнире $O$
$\mu_m$	b_m	Н·м·с/рад	0.0	Коэффициент вязкого трения между статором и ротором двигателя
$L_m$	Lm	Гн	0.01	Индуктивность обмотки двигателя
$R_m$	Rm	Ом	22	Сопротивление обмотки двигателя
$K_e$	Ke	рад/с/В	0.29	Скоростная постоянная
$K_t$	Kt	Н·м/А	0.27	Моментная постоянная
$V_{\max}$	Vmax	В	12	Максимальное напряжение двигателя

Нелинейная модель

Блок Nonlinear model описывает систему уравнений $(22)$. Преобразуйте этот блок так, чтобы он описывал нелинейную модель обратного маятника $(10)$. Используйте следующий шаблон кода:

# constants
a = ...
b = ...

# function that describes nonlinear dynamics
def nonlinear_dynamics(theta, d_theta, d_phi, voltage):
    # Input:
    #   theta - pendulum angle [rad]
    #   d_theta - pendulum angular velocity [rad/s]
    #   d_phi - wheel angular velocity [rad/s]
    #   voltage - motor input voltage [volt]
    _dd_theta = ...
    _dd_phi = ...
    _theta = theta + d_theta * dt
    _d_theta = d_theta + _dd_theta * dt
    _d_phi = d_phi + _dd_phi * dt
    return np.array([_theta, _d_theta, _d_phi])

Чтобы удостовериться в корректности работы вашего кода, используйте блок Nonlinear Model Simulation, который производит симуляцию и вывод графика переходных процессов. При правильно составленной модели вы должны увидеть такой график, как на рис. 11.

Рисунок 11. Результаты моделирования нелинейной модели.

Вопрос

Объясните, почему объект ведёт себя именно так. Попробуйте изменить параметры x0 и u, чтобы лучше понять его поведение.

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

В блоке State Space Model описывается линеаризованная модель двигателя в пространстве состояний $(23)$. Преобразуйте этот блок так, чтобы он описывал модель маятника $(12)$. Используйте приведённый ниже шаблон:

# constants
a = ...
b = ...

# elements of matrix A and vector B
a21 = ...
a22 = ...
a23 = ...
a31 = ...
a32 = ...
a33 = ...
b2 = ...
b3 = ...

# State transition matrix A and control vector B
A = np.array([
    [0,   1,   0],
    [a21, a22, a23],
    [a31, a32, a33],
])
B = np.array([
    [0],
    [b2],
    [b3],
])

# Output matrices and vectors
C = np.array([
    [1, 0, 0],  # pend angle
    [0, 1, 0],  # pend angular vel
    [0, 0, 1],  # wheel speed
])
D = np.array([[0], [0], [0]])

Чтобы удостовериться в корректности работы вашего кода, используйте блок State Space Model Simulation, который производит симуляцию и вывод графика переходных процессов. При правильно составленной модели вы получите график, как на рис. 12.

Рисунок 12. Результаты моделирования линейной модели.

Объясните, почему объект ведёт себя именно так. Попробуйте изменить параметры x0 и u, чтобы лучше понять его поведение.

Дискретизация, управляемость и стабилизация

• Блок Convert to Discrete Time Domain производит преобразование непрерывной модели в дискретную, используя средства библиотеки Control Systems Library. Пользуясь уравнениями (14) и (15), напишите собственную реализацию такого преобразования.
• В блоке Controllability Analysis произведите анализ управляемости системы. Переходите к следующим пунктам только если система управляема. Если окажется, что система не управляема, то скорее всего вы допустили ошибку во время выполнения предыдущих пунктов.
• Блок Stabilization in upright position via LQR Control производит решение линейно-квадратичной задачи управления. Получение матрицы коэффициентов обратной связи $K$ производится с помощью функции dlqr() библиотеки Control Systems Library. В библиотеке scipy существует функция solve_discrete_are(Ad, Bd, Q, R), решающая дискретное уравнение Риккати. На вход она принимает значения матриц Ad, Bd, Q и R в указанном порядке, а на выходе даёт решение $P$. Она уже импортирована в проект в блоке 1. Напишите собственную реализацию получения матрицы $K$, пользуясь предложенной функцией и уравнением $(18)$.

подсказка

Матричное умножение в Python делается при помощи оператора @.

• В блоке Stability Analysis проанализируйте устойчивость замкнутой системы, пользуясь критерием устойчивости дискретных систем из теоретической части. Чтобы получить собственные числа матрицы, используйте функцию библиотеки numpy:

eigvals = np.linalg.eigvals(M)

• На вход, вместо M, подайте матрицу $A_d - B_dK$, на выходе вы получите все собственные числа.
• В блоках Closed-Loop System Simulation производится симуляция замкнутой системы, используя линейные и нелинейные уравнения. Сравните полученные результаты. Поэкспериментируйте с весовыми коэффициентами матриц Q и R, посмотрите, как меняются переходные процессы. Выберите наилучшие, по вашему мнению, коэффициенты.
• В блоке Swing-Up Control производится управление раскачиванием маятника с последующим переходом в режим стабилизации. Коэффициент обратной связи $k$ и порог переключения режимов управления $\varepsilon$ заданы случайным образом. Вам необходимо подобрать такие значения этих переменных, чтобы маятник, раскачавшись, смог перейти в режим стабилизации. Значение $k$ подбирайте, опираясь на формируемый в блоке график Energy. Старайтесь подобрать его таким, чтобы ошибка между $E$ и $E^*$ была минимальной. Значение $\varepsilon$ выбирайте таким, чтобы избежать преждевременного переключения в режим стабилизации, когда маятник ещё не успел достичь области притяжения. В нашем случае максимально возможная ошибка между $E$ и $E^*$ составляет приблизительно 0.7 Дж. Соответственно, если вы поставите, например, $\varepsilon = 0.4$, то контроллер перейдёт в режим стабилизации, не успев раскачаться.

Эксперимент

В процессе выполнения вышеописанных пунктов вы должны были получить коэффициенты матрицы $K$ для режима стабилизации, коэффициент обратной связи $k$ для управления раскачиванием и порог переключения режимов управления $\varepsilon$. Эти значения понадобятся для управления аппаратной платформой маятника. Чтобы применить их, проделайте следующие действия:

1. Откройте файл Firmware/wheel_pendulum_practicum/wheel_pendulum_practicum.ino с помощью Arduino IDE. В этом файле реализовано взаимодействие с датчиками, управление двигателем и работа в режимах раскачивания, торможения и стабилизации.
2. Объявление коэффициентов обратной связи $K$ и $k$ вместе с порогом переключения $\varepsilon$ реализовано в строках 17-21. Функция control() вызывается каждые 10 мс. В ней реализованы контуры стабилизации, раскачивания и торможения. Изучите эту функцию и разберитесь, при каких условиях включается тот или иной режим управления.
3. Введите полученные вами в ходе работы значения коэффициентов вместо указанных в строках 17-21.
4. Подключите программатор к плате управления маятником, а затем к компьютеру через порт USB.
5. Скомпилируйте и загрузите прошивку кнопкой Upload в Arduino IDE.

Рисунок 13. Микроконтроллерный модуль VB32G4.

6. После выполнения загрузки прошивки отсоедините программатор от платы управления и нажмите кнопку RST (рис. 13). Двигатель начнёт вращаться. В этот момент система инициализирует все необходимые параметры управления двигателем. Дождитесь окончания вращения двигателя.
7. Нажмите кнопку USR. Должно включиться управление маятником.

Если через секунду после нажатия на кнопку USR двигатель не пришёл в движение, сделайте следующее:

• качните маятник рукой;
• отключите программатор, выключите и снова включите плату управления.

Если вы включили управление маятником, он раскачался и встал в режим стабилизации, значит вы правильно настроили коэффициенты обратной связи. Часто при разработке алгоритмов управления получается так, что коэффициенты, отлично подходящие для модели, не всегда стабилизируют реальный объект. Это связано с расхождением между моделью и реальностью. Чаще всего это происходит из-за погрешности вычислений физических параметров объекта управления и неучтённых параметров в уравнениях движения, например силы трения, упругости или ограничения по скорости двигателя. Линеаризация также влияет на расхождение, так как происходит пренебрежение некоторыми параметрами и переменными. В таком случае сначала стоит попробовать подобрать другие коэффициенты. Если это не поможет, необходимо произвести уточнение модели. Соответственно, в случае проблем со стабилизацией попробуйте подобрать другие коэффициенты, меняя весовые коэффициенты матриц Q и R. Также в случае проблем с раскачиванием или порогом переключения попробуйте изменить коэффициенты $k$ и $\varepsilon$. Если ничего не помогло, ещё раз убедитесь, что вы правильно составили математическую модель в файле model.ipynb.