Урок 3 — Анализ устойчивости систем управления¶

1. Понятие устойчивости¶

Устойчивость – это одно из важнейших свойств любой динамической системы, включая летательные аппараты. Говоря простыми словами, устойчивость – это способность системы возвращаться в исходное состояние равновесия после прекращения действия возмущающих сил.

Представьте себе шарик в разных положениях:

Устойчивое равновесие: Шарик на дне ямы. Если его толкнуть, он покачается и вернется обратно на дно.
Неустойчивое равновесие: Шарик на вершине холма. Малейший толчок – и он скатится вниз, уже не возвращаясь.
Безразличное равновесие: Шарик на ровной горизонтальной поверхности. Если его толкнуть, он остановится в новом положении и останется там.

Для летательного аппарата устойчивость означает, что после воздействия порыва ветра или кратковременного отклонения рулей он самостоятельно вернется к исходному режиму полета (например, к прямолинейному горизонтальному полету).

2. Виды устойчивости¶

Различают два основных вида устойчивости:

Статическая устойчивость: Способность системы создавать моменты и силы, стремящиеся вернуть ее в исходное положение равновесия, сразу после возникновения возмущения. Это начальная тенденция к возврату.
Динамическая устойчивость: Характер переходного процесса возврата системы в равновесие. Динамически устойчивая система не просто стремится вернуться, но и делает это без нарастающих колебаний. Процесс затухает со временем.

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

3. Устойчивость летательного аппарата¶

Движение летательного аппарата для удобства анализа разделяют на два вида:

Продольное движение: Движение в вертикальной плоскости (изменение высоты, угла тангажа, скорости).
Боковое движение: Движение в горизонтальной плоскости (изменение курса, угла крена, угла скольжения).

Соответственно, и устойчивость рассматривают отдельно для каждого вида движения [1]:

Продольная устойчивость: Способность ЛА сохранять заданный угол атаки и скорость полета.
Путевая (курсовая) устойчивость: Способность ЛА сохранять заданное направление полета (курс).
Поперечная устойчивость: Способность ЛА восстанавливать исходный угол крена.

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

4. Критерии устойчивости линейных систем¶

Для анализа устойчивости линейных систем, описанных в пространстве состояний, ключевую роль играют собственные числа (eigenvalues) матрицы системы A.

Система ẋ = Ax является устойчивой, если все действительные части всех собственных чисел матрицы A отрицательны.

Re(λᵢ) < 0 для всех i = 1, ..., n

где λᵢ – собственные числа матрицы A.

Если хотя бы одно собственное число имеет положительную действительную часть, система будет неустойчивой. Если есть собственные числа с нулевой действительной частью (и нет с положительной), система находится на границе устойчивости.

4.1. Алгебраические критерии (Рауса-Гурвица)¶

Эти критерии позволяют определить, есть ли у характеристического уравнения системы (которое находится как det(A - λI) = 0) корни с положительной действительной частью, не вычисляя сами корни. Это удобно для аналитического исследования.

4.2. Частотные критерии (Найквиста, Михайлова)¶

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

5. Ссылки на материалы¶

Алексеенков В. Устойчивость летательных аппаратов (презентация). – Авиационный колледж. – URL: http://taviak.ru/distance/wp-content/uploads/2013/PLA/Ustoi_chivost_letatel_nykh_apparatov(Alekseenkov).pdf
Продольная устойчивость летательного аппарата // Большая российская энциклопедия. – URL: https://bigenc.ru/c/prodol-naia-ustoichivost-letatel-nogo-apparata-5e7446
Устойчивость и управляемость самолета // Vzletim.ru. – URL: https://vzletim.ru/upload/iblock/133/aerodynamics09.pdf

6. Практический пример на Python¶

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

import numpy as np
import matplotlib.pyplot as plt
import control as ctrl
from scipy.linalg import eigvals
import sympy as sp

def analyze_stability(A, system_name="Система"):
    """Анализ устойчивости системы по собственным числам матрицы A"""

    print(f"\n{'='*50}")
    print(f"АНАЛИЗ УСТОЙЧИВОСТИ: {system_name}")
    print(f"{'='*50}")

    print(f"Матрица A:")
    print(A)

    # Вычисление собственных чисел
    eigenvalues = eigvals(A)
    print(f"\nСобственные числа (полюсы):")

    stable = True
    for i, lam in enumerate(eigenvalues):
        real_part = np.real(lam)
        imag_part = np.imag(lam)

        if abs(imag_part) < 1e-10:  # действительное число
            print(f"λ_{i+1} = {real_part:.4f}")
        else:
            print(f"λ_{i+1} = {real_part:.4f} + {imag_part:.4f}j")

        if real_part >= 0:
            stable = False

    print(f"\nВсе действительные части отрицательны: {stable}")

    if stable:
        print("✓ СИСТЕМА УСТОЙЧИВА")
    else:
        print("✗ СИСТЕМА НЕУСТОЙЧИВА")

    return eigenvalues, stable

def plot_poles(eigenvalues, title="Расположение полюсов"):
    """Построение карты полюсов на комплексной плоскости"""

    plt.figure(figsize=(8, 6))

    for lam in eigenvalues:
        plt.plot(np.real(lam), np.imag(lam), 'rx', markersize=10, markeredgewidth=2)

    # Добавляем ось мнимых чисел (граница устойчивости)
    y_max = max(abs(np.imag(eigenvalues))) * 1.2 if len(eigenvalues) > 0 else 1
    plt.axvline(x=0, color='k', linestyle='--', alpha=0.5, label='Граница устойчивости')

    # Заштриховываем область устойчивости
    x_min = min(np.real(eigenvalues)) * 1.2 if len(eigenvalues) > 0 else -1
    plt.axvspan(x_min, 0, alpha=0.2, color='green', label='Область устойчивости')

    plt.grid(True, alpha=0.3)
    plt.xlabel('Действительная часть')
    plt.ylabel('Мнимая часть')
    plt.title(title)
    plt.legend()
    plt.axis('equal')

    return plt.gcf()

# Пример 1: Устойчивая система (масса-пружина-демпфер)
A1 = np.array([[0, 1],
               [-4, -1]])

eigenvals1, stable1 = analyze_stability(A1, "Масса-пружина-демпфер")

# Пример 2: Неустойчивая система (обращенный маятник)
A2 = np.array([[0, 1],
               [1, 0]])  # без демпфирования и с "отрицательной пружиной"

eigenvals2, stable2 = analyze_stability(A2, "Обращенный маятник")

# Пример 3: Система на границе устойчивости
A3 = np.array([[0, 1],
               [-1, 0]])  # консервативный осциллятор

eigenvals3, stable3 = analyze_stability(A3, "Консервативный осциллятор")

# Пример 4: Система с комплексными полюсами
A4 = np.array([[0, 1, 0],
               [-2, -0.5, 0],
               [0, 0, -3]])

eigenvals4, stable4 = analyze_stability(A4, "Система третьего порядка")

# Построение карт полюсов
plt.figure(figsize=(15, 10))

plt.subplot(2, 2, 1)
for lam in eigenvals1:
    plt.plot(np.real(lam), np.imag(lam), 'rx', markersize=10, markeredgewidth=2)
plt.axvline(x=0, color='k', linestyle='--', alpha=0.5)
plt.axvspan(-5, 0, alpha=0.2, color='green')
plt.grid(True, alpha=0.3)
plt.title('Устойчивая система')
plt.xlabel('Re(λ)')
plt.ylabel('Im(λ)')

plt.subplot(2, 2, 2)
for lam in eigenvals2:
    plt.plot(np.real(lam), np.imag(lam), 'rx', markersize=10, markeredgewidth=2)
plt.axvline(x=0, color='k', linestyle='--', alpha=0.5)
plt.axvspan(-2, 0, alpha=0.2, color='green')
plt.axvspan(0, 2, alpha=0.2, color='red')
plt.grid(True, alpha=0.3)
plt.title('Неустойчивая система')
plt.xlabel('Re(λ)')
plt.ylabel('Im(λ)')

plt.subplot(2, 2, 3)
for lam in eigenvals3:
    plt.plot(np.real(lam), np.imag(lam), 'rx', markersize=10, markeredgewidth=2)
plt.axvline(x=0, color='k', linestyle='--', alpha=0.5)
plt.axvspan(-2, 0, alpha=0.2, color='green')
plt.grid(True, alpha=0.3)
plt.title('Граница устойчивости')
plt.xlabel('Re(λ)')
plt.ylabel('Im(λ)')

plt.subplot(2, 2, 4)
for lam in eigenvals4:
    plt.plot(np.real(lam), np.imag(lam), 'rx', markersize=10, markeredgewidth=2)
plt.axvline(x=0, color='k', linestyle='--', alpha=0.5)
plt.axvspan(-4, 0, alpha=0.2, color='green')
plt.grid(True, alpha=0.3)
plt.title('Система 3-го порядка')
plt.xlabel('Re(λ)')
plt.ylabel('Im(λ)')

plt.tight_layout()
plt.show()

# Демонстрация переходных процессов
print(f"\n{'='*60}")
print("СРАВНЕНИЕ ПЕРЕХОДНЫХ ПРОЦЕССОВ")
print(f"{'='*60}")

# Создание систем
sys1 = ctrl.ss(A1, np.array([[0], [1]]), np.array([[1, 0]]), 0)
sys2 = ctrl.ss(A2, np.array([[0], [1]]), np.array([[1, 0]]), 0)
sys3 = ctrl.ss(A3, np.array([[0], [1]]), np.array([[1, 0]]), 0)

# Переходные характеристики
t = np.linspace(0, 10, 1000)

plt.figure(figsize=(15, 5))

# Устойчивая система
plt.subplot(1, 3, 1)
try:
    t1, y1 = ctrl.step_response(sys1, t)
    plt.plot(t1, y1)
    plt.title('Устойчивая система\n(затухающий процесс)')
    plt.xlabel('Время (с)')
    plt.ylabel('Выход')
    plt.grid(True)
except:
    plt.text(0.5, 0.5, 'Ошибка моделирования', ha='center', va='center')

# Неустойчивая система
plt.subplot(1, 3, 2)
try:
    t2, y2 = ctrl.step_response(sys2, t)
    plt.plot(t2, y2)
    plt.title('Неустойчивая система\n(расходящийся процесс)')
    plt.xlabel('Время (с)')
    plt.ylabel('Выход')
    plt.grid(True)
    plt.ylim([-10, 10])  # ограничиваем масштаб для наглядности
except:
    plt.text(0.5, 0.5, 'Система неустойчива\n(расходящийся процесс)', ha='center', va='center')

# Система на границе устойчивости
plt.subplot(1, 3, 3)
try:
    t3, y3 = ctrl.step_response(sys3, t)
    plt.plot(t3, y3)
    plt.title('Граница устойчивости\n(незатухающие колебания)')
    plt.xlabel('Время (с)')
    plt.ylabel('Выход')
    plt.grid(True)
except:
    plt.text(0.5, 0.5, 'Ошибка моделирования', ha='center', va='center')

plt.tight_layout()
plt.show()

# Критерий Рауса-Гурвица (символьный анализ)
print(f"\n{'='*60}")
print("КРИТЕРИЙ РАУСА-ГУРВИЦА (символьный анализ)")
print(f"{'='*60}")

def routh_hurwitz_table(coeffs):
    """Построение таблицы Рауса-Гурвица"""
    n = len(coeffs)
    table = np.zeros((n, (n+1)//2))

    # Заполнение первых двух строк
    table[0, :] = coeffs[::2]  # четные коэффициенты
    if n > 1:
        table[1, :len(coeffs[1::2])] = coeffs[1::2]  # нечетные коэффициенты

    # Заполнение остальных строк
    for i in range(2, n):
        for j in range(table.shape[1]-1):
            if table[i-1, 0] != 0:
                table[i, j] = (table[i-1, 0] * table[i-2, j+1] - 
                              table[i-2, 0] * table[i-1, j+1]) / table[i-1, 0]

    return table

# Пример: характеристический полином s³ + 2s² + 3s + 4
coeffs = [1, 2, 3, 4]  # коэффициенты от старшей степени к младшей
print("Характеристический полином: s³ + 2s² + 3s + 4")
print("Таблица Рауса-Гурвица:")
routh_table = routh_hurwitz_table(coeffs)
print(routh_table)

# Проверка знаков в первом столбце
first_column = routh_table[:, 0]
sign_changes = sum(1 for i in range(len(first_column)-1) 
                  if first_column[i] * first_column[i+1] < 0)
print(f"\nКоличество перемен знака в первом столбце: {sign_changes}")
print(f"Количество полюсов в правой полуплоскости: {sign_changes}")

if sign_changes == 0:
    print("✓ Система устойчива по критерию Рауса-Гурвица")
else:
    print("✗ Система неустойчива по критерию Рауса-Гурвица")

Этот пример демонстрирует:

Анализ устойчивости через вычисление собственных чисел матрицы системы
Визуализацию полюсов на комплексной плоскости с выделением областей устойчивости
Сравнение переходных процессов для устойчивых, неустойчивых систем и систем на границе устойчивости
Применение критерия Рауса-Гурвица для символьного анализа устойчивости
Практические примеры различных типов динамических систем