ПИД-регулятор: теория, настройка и реальные примеры

Что такое ПИД и почему он везде

ПИД-регулятор (Пропорционально-Интегрально-Дифференциальный) — самый распространённый алгоритм автоматического управления в промышленности. По различным оценкам, более 90% всех промышленных регуляторов используют ПИД или его вариации.

Термостат в вашей духовке. Круиз-контроль в автомобиле. Система стабилизации квадрокоптера. Регулятор давления в насосной станции. Система автопилота самолёта. Везде — ПИД.

Почему? Потому что он:

• Прост в понимании — три компоненты, три параметра

• Работает для подавляющего большинства объектов управления

• Хорошо изучен — 80 лет теории и практики

• Легко реализуется — 10–20 строк кода

Математика: просто о сложном

Пусть у нас есть:

• SP (Setpoint) — желаемое значение (уставка)

• PV (Process Variable) — измеренное значение

• e(t) = SP – PV — ошибка регулирования

• u(t) — управляющий сигнал (выход регулятора)

Формула ПИД:

u(t) = Kp × e(t)  +  Ki × ∫e(t)dt  +  Kd × de(t)/dt
          │                │                    │
    Пропорциональная  Интегральная        Дифференциальная
    составляющая      составляющая         составляющая

Разберём каждую составляющую на примере: нагреватель, цель — держать 80°C.

P — Пропорциональная составляющая

Управляющий сигнал пропорционален текущей ошибке. Температура 60°C, уставка 80°C, ошибка 20°C → включить нагрев на 20×Kp%.

Проблема чистого P-регулятора: статическая ошибка (offset). При стабильном нагреве нагрев должен компенсировать теплопотери, значит ошибка никогда не станет нулём — иначе нагрев выключится и температура упадёт. Система "зависает" с постоянной небольшой ошибкой.

I — Интегральная составляющая

Суммирует ошибку со временем. Даже маленькая постоянная ошибка, накапливаясь, создаёт всё больший управляющий сигнал — и в конечном счёте устраняет статическую ошибку.

Проблема: интегральное насыщение (windup). Если выход ограничен (0–100%), а ошибка долго накапливается (например, нагрев не справляется при открытом окне), интеграл "разгоняется" до огромного значения. Потом, когда ошибка уменьшается, нужно много времени чтобы интеграл "разрядился" — система перерегулирует.

D — Дифференциальная составляющая

Реагирует на скорость изменения ошибки. Если ошибка быстро уменьшается (температура быстро растёт к уставке) — D-составляющая "притормаживает" для предотвращения перерегулирования.

Проблема: шум. Производная усиливает высокочастотный шум измерений. Любой дёрг датчика превращается в огромный кратковременный скачок управляющего сигнала. Поэтому D всегда применяется с фильтром.

Дискретная реализация (для цифровых систем)

В реальных системах регулятор выполняется с дискретным шагом Ts (время выборки). Непрерывная формула преобразуется в разностное уравнение:

P = Kp × e[k]
I = I[k-1] + Ki × e[k] × Ts
D = Kd × (e[k] - e[k-1]) / Ts

u[k] = P + I + D

Но правильная реализация сложнее. Полная профессиональная реализация:

class PIDController:
    """
    Промышленный ПИД-регулятор с:
    - Anti-windup (ограничение интеграла)
    - Производная по измерению (не по ошибке)
    - Фильтр производной
    - Ограничение выхода
    - Bump-less transfer (безударное переключение авто/ручной)
    """
    
    def __init__(self,
                 kp: float,
                 ki: float,
                 kd: float,
                 ts: float,              # Период дискретизации, секунды
                 out_min: float = 0.0,
                 out_max: float = 100.0,
                 filter_coeff: float = 0.1):  # Коэффициент фильтра D
        
        self.kp = kp
        self.ki = ki
        self.kd = kd
        self.ts = ts
        self.out_min = out_min
        self.out_max = out_max
        self.filter_coeff = filter_coeff  # Nf: чем меньше, тем сильнее фильтрация
        
        # Состояние
        self.setpoint  = 0.0
        self.integral  = 0.0
        self.prev_measurement = None
        self.filtered_deriv   = 0.0
        self.output    = 0.0
        
        # Режим
        self.auto_mode = True
        
    def set_setpoint(self, sp: float):
        self.setpoint = sp
    
    def set_manual(self, output: float):
        """Ручной режим — задать выход напрямую"""
        self.auto_mode = False
        self.output = max(self.out_min, min(self.out_max, output))
        # При переходе обратно в авто — интеграл синхронизируется
        
    def set_auto(self):
        """Переход в автоматический режим (bump-less)"""
        if not self.auto_mode:
            # Инициализируем интеграл текущим выходом
            # чтобы не было скачка при переключении
            self.integral = self.output
            self.auto_mode = True
    
    def compute(self, measurement: float) -> float:
        """
        Вычислить управляющий сигнал.
        measurement — текущее значение регулируемой величины.
        Вызывать строго с периодом ts!
        """
        if not self.auto_mode:
            return self.output
        
        if self.prev_measurement is None:
            self.prev_measurement = measurement
        
        error = self.setpoint - measurement
        
        # === Пропорциональная составляющая ===
        p_term = self.kp * error
        
        # === Интегральная составляющая с anti-windup ===
        # Предварительно вычисляем интеграл
        integral_new = self.integral + self.ki * error * self.ts
        
        # === Дифференциальная составляющая ===
        # Производная по ИЗМЕРЕНИЮ (не по ошибке) — избегаем derivative kick
        # при изменении уставки
        raw_deriv = -(measurement - self.prev_measurement) / self.ts
        
        # Фильтр производной (экспоненциальный)
        self.filtered_deriv = (self.filter_coeff * raw_deriv + 
                               (1 - self.filter_coeff) * self.filtered_deriv)
        
        d_term = self.kd * self.filtered_deriv
        
        # === Предварительный выход ===
        output_raw = p_term + integral_new + d_term
        
        # === Ограничение выхода ===
        self.output = max(self.out_min, min(self.out_max, output_raw))
        
        # === Anti-windup: обновляем интеграл только если не насыщены ===
        # Или используем "back-calculation anti-windup"
        if self.output == output_raw:
            # Нет насыщения — обновляем интеграл нормально
            self.integral = integral_new
        else:
            # Насыщение — не накапливаем интеграл дальше
            # Back-calculation: уменьшаем интеграл на величину насыщения
            saturation = output_raw - self.output
            self.integral = integral_new - saturation
        
        self.prev_measurement = measurement
        
        return self.output
    
    def get_components(self) -> dict:
        """Для отладки и мониторинга"""
        error = self.setpoint - (self.prev_measurement or 0)
        return {
            'setpoint':   self.setpoint,
            'error':      error,
            'p_term':     self.kp * error,
            'i_term':     self.integral,
            'd_term':     self.kd * self.filtered_deriv,
            'output':     self.output,
        }

Методы настройки коэффициентов

Метод 1: Инженерная настройка (руками)

Алгоритм по шагам:

1. Ki = 0, Kd = 0. Постепенно увеличиваем Kp до появления устойчивых колебаний (нарастающих) — это Kc. Уменьшаем Kp до 0.5 × Kc.

2. Медленно увеличиваем Ki. Он должен устранить остаточную ошибку. Если появились колебания — Ki велик.

3. Если нужно улучшить реакцию и уменьшить перерегулирование — добавляем Kd. Часто тепловые объекты обходятся без D.

Метод 2: Циглера-Николса (классика)

1. Установить Ki = 0, Kd = 0

2. Увеличивать Kp до возникновения незатухающих колебаний

3. Записать: Ku — критический коэффициент, Tu — период колебаний

Результат: Обычно агрессивный, дающий ~25% перерегулирования. Хорош как стартовая точка.

Метод 3: Ступенчатый отклик (Step Response)

Более безопасный метод — не нужно доводить до автоколебаний:

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

2. Скачком изменить управляющий сигнал на 10–20%

3. Записать переходный процесс

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

def fit_first_order_plus_delay(t, data, step_size):
    """
    Аппроксимация объекта управления моделью первого порядка с запаздыванием:
    G(s) = K × exp(-L×s) / (T×s + 1)
    """
    # Нормализация данных
    y_norm = (data - data[0]) / step_size
    
    # Начальные параметры: K, T, L
    K_init = y_norm[-1]         # Статический коэффициент усиления
    T_init = t[np.argmin(np.abs(y_norm - 0.63 * K_init))]
    L_init = t[np.argmin(np.abs(y_norm - 0.05 * K_init))]
    
    def model(t, K, T, L):
        result = np.zeros_like(t)
        for i, ti in enumerate(t):
            if ti > L:
                result[i] = K * (1 - np.exp(-(ti - L) / T))
        return result
    
    try:
        popt, _ = curve_fit(model, t, y_norm,
                            p0=[K_init, T_init, max(L_init, 0.001)],
                            bounds=([0, 0, 0], [np.inf, np.inf, np.inf]))
        K, T, L = popt
        return K, T, L
    except RuntimeError:
        return K_init, T_init, max(L_init, 0.001)

def imc_tuning(K, T, L, lambda_factor=1.0):
    """
    IMC (Internal Model Control) настройка — лучший метод для объектов
    с запаздыванием.
    lambda_factor: чем больше — тем медленнее, но стабильнее (рекомендуется L...3L)
    """
    lambda_c = lambda_factor * max(L, 0.1 * T)
    
    Kp = (2 * T + L) / (2 * K * (lambda_c + L))
    Ti = T + L / 2  # Постоянная интегрирования
    Td = T * L / (2 * T + L)  # Постоянная дифференцирования
    
    Ki = Kp / Ti
    Kd = Kp * Td
    
    return {
        'Kp': round(Kp, 4),
        'Ki': round(Ki, 4),
        'Kd': round(Kd, 4),
        'Ti': round(Ti, 4),
        'Td': round(Td, 4),
    }

# Пример использования:
# K=2.0 (усиление объекта), T=30с (постоянная времени), L=5с (запаздывание)
params = imc_tuning(K=2.0, T=30.0, L=5.0, lambda_factor=1.0)
print(f"Kp={params['Kp']}, Ki={params['Ki']}, Kd={params['Kd']}")

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

Перед применением в реальной системе — всегда симулируйте:

import numpy as np
import matplotlib.pyplot as plt

def simulate_pid(kp, ki, kd, setpoint, sim_time=200, ts=0.1,
                 K=2.0, T=30.0, L=5.0, noise_std=0.1):
    """
    Симуляция ПИД с объектом первого порядка + запаздывание.
    Возвращает массивы времени, выходного значения, управляющего сигнала.
    """
    
    steps = int(sim_time / ts)
    delay_steps = int(L / ts)
    
    t       = np.zeros(steps)
    pv      = np.zeros(steps)
    u       = np.zeros(steps)
    sp      = np.zeros(steps)
    
    pid = PIDController(kp, ki, kd, ts, out_min=0.0, out_max=100.0)
    
    # Ступенчатое изменение уставки
    sp[:steps//4] = 0
    sp[steps//4:] = setpoint
    
    # Буфер запаздывания
    u_delayed = np.zeros(delay_steps + 1)
    
    for i in range(1, steps):
        t[i] = i * ts
        sp_now = sp[i]
        pid.set_setpoint(sp_now)
        
        # Добавляем шум измерения
        noise = np.random.normal(0, noise_std)
        
        # Управляющий сигнал
        u[i] = pid.compute(pv[i-1] + noise)
        
        # Обновляем буфер запаздывания
        u_delayed = np.roll(u_delayed, 1)
        u_delayed[0] = u[i]
        
        # Объект управления: первый порядок с запаздыванием
        # dy/dt = (K × u_delayed - y) / T
        pv[i] = pv[i-1] + ts * (K * u_delayed[-1] - pv[i-1]) / T
    
    return t, pv, u, sp

# Сравнение разных настроек
fig, axes = plt.subplots(2, 1, figsize=(12, 8))

configs = [
    {'kp': 0.5, 'ki': 0.01, 'kd': 2.0, 'label': 'IMC (lambda=1)', 'color': 'blue'},
    {'kp': 1.5, 'ki': 0.05, 'kd': 5.0, 'label': 'Агрессивный', 'color': 'red'},
    {'kp': 0.2, 'ki': 0.005, 'kd': 0.5, 'label': 'Консервативный', 'color': 'green'},
]

for cfg in configs:
    t, pv, u, sp = simulate_pid(
        kp=cfg['kp'], ki=cfg['ki'], kd=cfg['kd'],
        setpoint=60.0
    )
    axes[0].plot(t, pv, label=cfg['label'], color=cfg['color'])
    axes[1].plot(t, u, color=cfg['color'], alpha=0.7)

axes[0].plot(t, sp, 'k--', label='Уставка', linewidth=2)
axes[0].set_ylabel('Температура, °C')
axes[0].set_title('Сравнение настроек ПИД-регулятора')
axes[0].legend()
axes[0].grid(True)

axes[1].set_xlabel('Время, с')
axes[1].set_ylabel('Управляющий сигнал, %')
axes[1].grid(True)

plt.tight_layout()
plt.savefig('pid_comparison.png', dpi=150)

Специальные случаи и решения

Каскадный ПИД (Cascade Control)

Используется когда одного ПИД недостаточно. Классический пример: температурный реактор с нагревателем и теплоносителем.

Внешний ПИД (медленный):
    SP температуры → [ПИД_темп] → SP расхода теплоносителя

Внутренний ПИД (быстрый):
    SP расхода → [ПИД_расход] → Управление клапаном

Внутренний контур в 5–10 раз быстрее внешнего. Преимущество: внутренний контур компенсирует возмущения по теплоносителю до того, как они повлияют на температуру.

class CascadePID:
    """Каскадный регулятор с двумя ПИД"""
    
    def __init__(self, outer_pid: PIDController, inner_pid: PIDController):
        self.outer = outer_pid
        self.inner = inner_pid
    
    def compute(self, outer_measurement: float, 
                       inner_measurement: float) -> float:
        """
        outer_measurement — медленная переменная (температура)
        inner_measurement — быстрая переменная (расход)
        """
        # Внешний контур вырабатывает уставку для внутреннего
        inner_setpoint = self.outer.compute(outer_measurement)
        
        # Ограничиваем уставку внутреннего контура
        inner_setpoint = max(0.0, min(100.0, inner_setpoint))
        self.inner.set_setpoint(inner_setpoint)
        
        # Внутренний контур управляет исполнительным устройством
        return self.inner.compute(inner_measurement)

ПИД с ограничением скорости изменения выхода (Rate Limiting)

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

def rate_limited_output(new_output: float, prev_output: float,
                        max_rate: float, ts: float) -> float:
    """
    Ограничивает скорость изменения выходного сигнала.
    max_rate — максимальное изменение в секунду (например, 10%/сек)
    """
    max_change = max_rate * ts
    change = new_output - prev_output
    change = max(-max_change, min(max_change, change))
    return prev_output + change

ПИД для нестационарных объектов (Gain Scheduling)

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

class GainSchedulingPID:
    """ПИД с переключением коэффициентов в зависимости от уставки"""
    
    def __init__(self):
        # Коэффициенты для разных диапазонов температур
        self.schedules = [
            # (max_sp, kp, ki, kd)
            (50.0,  2.0, 0.05, 1.0),   # Низкая температура
            (100.0, 1.5, 0.03, 0.8),   # Средняя
            (200.0, 1.0, 0.02, 0.5),   # Высокая (объект другой!)
        ]
        
        self.pid = PIDController(kp=2.0, ki=0.05, kd=1.0, ts=1.0)
    
    def update_gains(self, setpoint: float):
        for max_sp, kp, ki, kd in self.schedules:
            if setpoint <= max_sp:
                self.pid.kp = kp
                self.pid.ki = ki
                self.pid.kd = kd
                break
    
    def compute(self, setpoint: float, measurement: float) -> float:
        self.update_gains(setpoint)
        self.pid.set_setpoint(setpoint)
        return self.pid.compute(measurement)

Реализация на C для встраиваемых систем

#include <stdint.h>
#include <float.h>

typedef struct {
    float kp;
    float ki;
    float kd;
    float ts;           // Период дискретизации, секунды
    float out_min;
    float out_max;
    float filter_coeff; // Коэффициент фильтра производной (0.01..0.2)
    
    // Состояние
    float setpoint;
    float integral;
    float prev_measurement;
    float filtered_deriv;
    float output;
    
    uint8_t first_run;
} PIDState;

void PID_Init(PIDState *pid, float kp, float ki, float kd, float ts,
              float out_min, float out_max)
{
    pid->kp = kp;
    pid->ki = ki;
    pid->kd = kd;
    pid->ts = ts;
    pid->out_min = out_min;
    pid->out_max = out_max;
    pid->filter_coeff = 0.1f;
    
    pid->setpoint          = 0.0f;
    pid->integral          = 0.0f;
    pid->prev_measurement  = 0.0f;
    pid->filtered_deriv    = 0.0f;
    pid->output            = 0.0f;
    pid->first_run         = 1;
}

void PID_SetSetpoint(PIDState *pid, float sp)
{
    pid->setpoint = sp;
}

void PID_Reset(PIDState *pid)
{
    pid->integral         = 0.0f;
    pid->filtered_deriv   = 0.0f;
    pid->first_run        = 1;
}

float PID_Compute(PIDState *pid, float measurement)
{
    float error, p_term, d_term, output_raw;
    float raw_deriv, integral_new;
    
    if (pid->first_run) {
        pid->prev_measurement = measurement;
        pid->first_run = 0;
    }
    
    error = pid->setpoint - measurement;
    
    // Пропорциональная
    p_term = pid->kp * error;
    
    // Интегральная (предварительный расчёт)
    integral_new = pid->integral + pid->ki * error * pid->ts;
    
    // Дифференциальная по измерению
    raw_deriv = -(measurement - pid->prev_measurement) / pid->ts;
    
    // Фильтр производной (IIR первого порядка)
    pid->filtered_deriv = pid->filter_coeff * raw_deriv +
                          (1.0f - pid->filter_coeff) * pid->filtered_deriv;
    d_term = pid->kd * pid->filtered_deriv;
    
    // Суммируем
    output_raw = p_term + integral_new + d_term;
    
    // Ограничение выхода
    if (output_raw > pid->out_max)
        pid->output = pid->out_max;
    else if (output_raw < pid->out_min)
        pid->output = pid->out_min;
    else
        pid->output = output_raw;
    
    // Anti-windup: обновляем интеграл только если нет насыщения
    if (pid->output == output_raw) {
        pid->integral = integral_new;
    } else {
        // Back-calculation
        float saturation = output_raw - pid->output;
        pid->integral = integral_new - saturation;
    }
    
    pid->prev_measurement = measurement;
    
    return pid->output;
}

// Использование на STM32/Arduino (в прерывании таймера):
PIDState temp_pid;

void PID_ISR(void)  // Вызывается каждые 100 мс
{
    float measured_temp = read_temperature_sensor();
    float heater_power = PID_Compute(&temp_pid, measured_temp);
    
    // Перевести в ШИМ 0-255
    uint8_t pwm = (uint8_t)(heater_power * 255.0f / 100.0f);
    set_heater_pwm(pwm);
}

Диагностика плохой работы ПИД

Симптом: Медленный выход на уставку, большая статическая ошибка

• Мало Kp и/или Ki

• Проверьте нет ли механических ограничений (клапан не открывается полностью)

• Проверьте масштабирование — правильно ли в % интерпретируется выход?

Симптом: Сильные устойчивые колебания

• Kp слишком велик

• Ki слишком велик (самая частая причина!)

• Период дискретизации слишком большой относительно динамики объекта

Симптом: Большое перерегулирование при изменении уставки

• Велико Kp, увеличьте Kd

• Используйте "setpoint weighting": P-составляющую считайте как Kp × (w × SP – PV), где w < 1

Симптом: Шумный выход, дёргает исполнительный механизм

• Шум датчика, увеличьте filter_coeff для фильтрации D

• Уменьшите Kd

• Добавьте мёртвую зону: не реагировать если |error| < deadband

Заключение

ПИД-регулятор — это не просто формула, это искусство. Одни и те же Kp, Ki, Kd дадут отличные результаты на одном объекте и катастрофу на другом. Понимание физики объекта управления важнее знания формул.

Практический совет: начинайте с IMC-метода настройки — он даёт предсказуемые результаты и физический смысл у параметра lambda (желаемая скорость реакции). Используйте симуляцию перед внедрением. Добавляйте anti-windup всегда — без него реальная система будет вести себя непредсказуемо после насыщения.

И помните: 80% задач управления решается PI-регулятором (без D). D добавляйте только когда PI действительно недостаточно и когда шум измерений не является проблемой.