Embedded HandbookПриложения
Все разделы › XIV. Машинное обучение на микроконтроллере
STM32F407, STM32F103 · уровень продвинутый · 12 мин чтения
Открыть в интерактивном справочнике

121. Классика против ML: фильтр Калмана и оценка заряда батареи

Заряд батареи нельзя измерить — его можно только оценить по току и напряжению. Задача кажется созданной для нейросети: есть входы, есть ответ, данных сколько угодно со стенда. Раздел ставит честное сравнение: счёт кулонов, таблица по напряжению, расширенный фильтр Калмана на физической модели батареи и нейросеть, обученная на тысячах секунд разряда, — на привычной нагрузке, на незнакомой и на постаревшей батарее. Итог части XIV: когда физика известна, её модель с обратной связью по измерению обыгрывает обученную — и особенно там, где мир ушёл от обучающих данных.

Опыты — в образе stm-ml:1; фильтр на C сверен с Python и собран для F407 и F103. Фильтр Калмана вводится здесь с нуля; модели и их обучение — разделы 116–118.

Как это называется в документации

Оригинал Перевод
state of charge (SoC) заряд: доля от полной ёмкости
coulomb counting счёт кулонов: интеграл тока
open-circuit voltage (OCV) напряжение холостого хода
equivalent circuit model схема замещения батареи
extended Kalman filter (EKF) расширенный фильтр Калмана
state, covariance состояние и его ковариация
process / measurement noise шум процесса / шум измерения
innovation невязка: измерение минус прогноз
Jacobian матрица производных (якобиан)
observability наблюдаемость: видно ли состояние по измерениям
sensor offset (bias) смещение нуля датчика
state of health, ageing состояние здоровья, старение
joint (augmented) estimation совместная оценка состояния и параметра
battery management system (BMS) система управления батареей

121.1 Батарея как схема

Схема замещения литиевой ячейки: источник напряжения холостого хода, зависящего от заряда (таблица из 11 точек от 3,00 до 4,20 В), последовательное сопротивление R0 = 0,05 Ом и RC-цепочка R1 = 0,03 Ом, C1 = 1000 Ф — медленная часть просадки с постоянной 30 с. Ёмкость 2,5 А·ч. Шаг 1 с. Датчик тока врёт на 0,05 А постоянно и шумит на 0,02 А, датчик напряжения шумит на 5 мВ — обычные цифры для шунта с усилителем и АЦП 12 бит.

Модель — для обычной ячейки на 4,2 В при комнатной температуре. На холоде R0 растёт в разы, поэтому в изделии параметры берут таблицей по температуре. У LiFePO4 напряжение в середине заряда почти не меняется (3,2–3,3 В): наклон мал, и фильтр там опирается в основном на счёт кулонов — поправку дают лишь края шкалы.

Три нагрузки проверки: тот же характер, что при обучении нейросети (броски 0,5–3 А разной длины); другой — ровные 1,2 А с короткими бросками до 4 А; тот же другой, но на постаревшей батарее — ёмкость на 12 % меньше, R0 в 1,6 раза больше. Все оценщики, кроме счёта кулонов, стартуют с заряда 0,5 при истинном 0,9: прибор перезапустился под нагрузкой и начального заряда не знает.

121.2 Счёт кулонов и таблица

Счёт кулонов складывает ток: заряд убывает на I·Δt / Q каждый шаг. Он гладкий и точный на коротком отрезке, но всё копит: смещение 0,05 А за полтора часа — это 3 % заряда, а неверная ёмкость у старой батареи добавляет ещё. И ему нужен начальный заряд — здесь его дали точно.

Таблица по напряжению — обратная задача: заряд, при котором напряжение холостого хода равно измеренному. Под нагрузкой измеренное напряжение просело на I·R0 и на RC-цепочку, и ошибка доходит до десятков процентов. Таблица верна только после долгого покоя — так её и используют: для начального заряда после ночи.

121.3 Расширенный фильтр Калмана

Фильтр Калмана соединяет оба способа: прогноз по току, как счёт кулонов, и поправка по напряжению, как таблица, — с весом, который он сам считает из того, насколько доверяет прогнозу и измерению. Состояние — заряд и напряжение RC-цепочки; напряжение на выводах зависит от заряда нелинейно, поэтому фильтр расширенный: в поправке участвует наклон таблицы в текущей точке.

Шаг фильтра на C — прогноз:

e->x[0] -= i * DT_Q;
e->x[1] = A_RC * e->x[1] + R1 * (1.0f - A_RC) * i;

Поправка по невязке — разнице измеренного и предсказанного напряжения:

yhat = ocv(e->x[0], &slope) - e->x[1] - e->x[2] * i;
S = H[0] * PH[0] + H[1] * PH[1] + H[2] * PH[2] + R_V;
y = v - yhat;
K[r] = PH[r] / S;
e->x[r] += K[r] * y;

Код взят из фильтра, который работает на приборе (121.8), и состояние в нём уже из трёх чисел: третье, x[2], — сопротивление R0 (121.5). В простом фильтре на два состояния вместо e->x[2] стоит постоянная R0 = 0,05 Ом, а H, PH и K — из двух элементов.

Кроме оценки, фильтр ведёт ковариацию P — насколько он в оценке уверен. Прогноз её растит, P = F·P·Fᵀ + Q; матрица перехода F диагональна — 1, A_RC, 1, — поэтому хватает двух циклов:

for (r = 0; r < 3; r++)
    for (c = 0; c < 3; c++)
        e->P[r][c] *= F[r] * F[c];
for (r = 0; r < 3; r++)
    e->P[r][r] += QD[r];

H — чувствительность напряжения на выводах к каждому состоянию: к заряду — наклон таблицы, к напряжению RC-цепочки — −1, к R0 — −i; PH = P·Hᵀ:

H[0] = slope;
H[1] = -1.0f;
H[2] = -i;
for (r = 0; r < 3; r++)
    PH[r] = e->P[r][0] * H[0] + e->P[r][1] * H[1] + e->P[r][2] * H[2];

После поправки ковариация сжимается, P = (I − K·H)·P:

for (r = 0; r < 3; r++)                 /* P = (I − K·H)·P, P симметрична: H·P = PHᵀ */
    for (c = 0; c < 3; c++)
        e->P[r][c] -= K[r] * PH[c];

Две настройки решают всё: Q — насколько фильтр не доверяет прогнозу, R_V — насколько не доверяет измерению. В опыте Q по диагонали — 10⁻⁷ для заряда, 10⁻⁶ В² для RC-цепочки и 10⁻¹⁰ Ом² для R0, а R_V = 10⁻⁴ В², то есть (10 мВ)²: вдвое больше шума датчика напряжения 5 мВ, с запасом на неточность модели. Больше Q — фильтр охотнее верит напряжению, больше R_V — крепче держится за счёт кулонов.

С неверного старта 0,5 фильтр сошёлся к истинному заряду с точностью 2 % за 50 секунд: несколько бросков тока — и напряжение выдало, где батарея на самом деле.

Первые 120 секунд EKF со старта 0,5 при истинном заряде 0,9: оценка с полосой ±2σ, ток с броском 4 А, ошибка меньше 2 % на 50-й секунде

121.4 Смещение датчика в состоянии — и наблюдаемость

Смещение тока копит ошибку у счёта кулонов; естественная мысль — добавить его в состояние фильтра, пусть оценит. Результат обратный ожидаемому: на новой батарее фильтр со смещением хуже простого, а само смещение он оценил как 0,085 А вместо 0,05. Смещение тока плохо наблюдаемо: его след в напряжении — ошибка I·R0 в 2,5 мВ, меньше шума датчика, а в заряде оно неотличимо от чуть другой ёмкости. Простой фильтр справляется со смещением и так: поправка по напряжению не даёт ошибке копиться. Лишнее состояние, которое измерения не видят, только добавляет шум.

121.5 R0 в состоянии: фильтр следит за старением

Старение меняет сразу два параметра: ёмкость и сопротивление. Фильтр с паспортными R0 и ёмкостью объясняет лишнюю просадку напряжения у старой батареи меньшим зарядом — и ошибается на 5–10 %. А вот сопротивление наблюдаемо отлично: каждый бросок тока показывает его прямо, как ступеньку напряжения. Третье состояние — R0, медленно блуждающее:

x = np.array([s0, 0.0, R0])
H = np.array([docv(x[0]), -1.0, -i])

Фильтр оценил R0 у старой батареи в 0,0787 Ом при истинных 0,08 — и ошибка заряда упала до процента. Ёмкость в состояние можно добавить так же, но она видна только за целый цикл разряда; в промышленных системах её подстраивают отдельно, по полным циклам.

121.6 Нейросеть

Нейросеть получает то же, что фильтр: напряжение и ток, плюс их средние за минуту — чтобы видеть медленную часть просадки. Два скрытых слоя по 16 нейронов, 369 весов, 48 252 обучающих отсчётов из 12 разрядов от разного начального заряда с той же схемой замещения. Преимущество у неё одно, но важное: она не накапливает состояние и не требует начального заряда — ответ по последней минуте.

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

121.7 Сравнение

Постаревшая батарея: ошибка заряда по времени — счёт кулонов уходит линейно, EKF с паспортными параметрами и нейросеть занижают заряд, EKF с R0 в состоянии держится у нуля

Среднеквадратичная / наибольшая ошибка заряда, %, после первых 10 минут:

Способ Та же нагрузка Другая нагрузка Постаревшая батарея
Счёт кулонов 1,8 / 3,0 1,8 / 3,0 4,5 / 7,3
Таблица по напряжению 17,3 / 37,9 16,6 / 48,4 21,7 / 61,3
EKF 0,5 / 0,9 0,6 / 0,9 5,6 / 10,6
EKF со смещением тока 1,1 / 1,8 0,9 / 1,5 2,8 / 5,7
EKF с R0 в состоянии 0,2 / 0,5 0,3 / 0,7 0,7 / 1,1
Нейросеть 4 → 16 → 16 → 1 0,9 / 3,4 1,0 / 4,0 6,7 / 20,9

Счёту кулонов дан точный начальный заряд — у остальных его нет. Фильтр с R0 в состоянии лучший во всех трёх случаях. Нейросеть проигрывает простому фильтру и фильтру с R0 везде, а на старой батарее ошибается сильнее даже счёта кулонов.

Ошибка оценки заряда на постаревшей батарее, среднеквадратичная и наибольшая в процентах: счёт кулонов 4,5 и 7,3, таблица по напряжению 21,7 и 61,3, EKF 5,6 и 10,6, EKF со смещением тока 2,8 и 5,7, EKF с R0 в состоянии 0,7 и 1,1 — лучший, нейросеть 6,7 и 20,9 — хуже счёта кулонов

121.8 На приборе

Фильтр на C — 12 чисел float: три числа состояния и матрица 3 × 3, состояние — 48 байт; таблица напряжения — 44 байта. armclang: 558 байт кода на F407 и 908 на F103, где float программный. Прогон на 4700 секундах старой батареи совпал с numpy до миллионной доли. Шаг раз в секунду — несколько сотен операций с плавающей точкой, на F103 без FPU это доли миллисекунды.

У нейросети 369 весов — 1476 байт во float, в int8 по разделу 119 около половины килобайта. Меньше фильтра она не выходит, а считать её приходится так же раз в секунду.

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

121.9 Когда что выбирать

121.10 Грабли

Фильтр сходится к неверному заряду. Неверные параметры модели: просадку от старения фильтр списывает на заряд. Сопротивление — в состояние, ёмкость — подстраивать по полным циклам.

Добавили состояние — стало хуже. Состояние плохо наблюдаемо: смещение тока оценено как 0,085 А вместо 0,05. Проверять, видит ли его измерение, и не добавлять лишнего.

Фильтр разошёлся через сутки работы. Матрица P потеряла симметрию или положительность из-за округления float. Симметризовать P, ограничивать снизу диагональ, форма Джозефа для обновления.

Счёт кулонов уплыл на проценты. Смещение нуля датчика тока. Калибровать ноль при нулевом токе, поправлять по напряжению после покоя.

Нейросеть уверенно врёт на новой батарее. Условия вне обучающих данных. Проверять на условиях, которых не было в обучении, до выпуска.

121.11 Итоги

Часть XIV закрыта: от выбора подхода (116) и признаков (117) через обучение (118), квантование (119) и модели нормы (120) к границе, где обученная модель уступает физической.

Полные исходники

ekf121.c — расширенный фильтр Калмана заряда батареи с состоянием [заряд, U_RC, R0].

/* Раздел 121: расширенный фильтр Калмана заряда батареи — состояние [заряд, U_RC, R0].
 * С -DML_HOST main() прогоняет измерения постаревшей батареи и сверяет с numpy. */
#include "model121.h"

typedef struct {
    float x[3];                 /* заряд 0…1, напряжение RC-цепочки, В, сопротивление R0, Ом */
    float P[3][3];
} ekf_t;

static const float QD[3] = { 1e-7f, 1e-6f, 1e-10f };

/* Напряжение холостого хода и его наклон по таблице из 11 точек через 0,1. */
static float ocv(float s, float *slope)
{
    int k = (int)(s * 10.0f);
    float t;
    if (k < 0)
        k = 0;
    if (k > N_PTS - 2)
        k = N_PTS - 2;
    *slope = (OCV_T[k + 1] - OCV_T[k]) * 10.0f;
    t = s * 10.0f - (float)k;
    if (s <= 0.0f)
        return OCV_T[0];
    if (s >= 1.0f)
        return OCV_T[N_PTS - 1];
    return OCV_T[k] + (OCV_T[k + 1] - OCV_T[k]) * t;
}

void ekf_init(ekf_t *e, float soc0)
{
    int i, j;
    for (i = 0; i < 3; i++)
        for (j = 0; j < 3; j++)
            e->P[i][j] = 0.0f;
    e->x[0] = soc0;
    e->x[1] = 0.0f;
    e->x[2] = R0;
    e->P[0][0] = 0.1f;
    e->P[1][1] = 0.01f;
    e->P[2][2] = 1e-4f;
}

/* Один шаг 1 с: прогноз по току, коррекция по напряжению. */
void ekf_step(ekf_t *e, float v, float i)
{
    float F[3] = { 1.0f, A_RC, 1.0f };     /* матрица перехода диагональна */
    float H[3], PH[3], S, K[3], y, slope, yhat;
    int r, c;

    e->x[0] -= i * DT_Q;
    e->x[1] = A_RC * e->x[1] + R1 * (1.0f - A_RC) * i;
    for (r = 0; r < 3; r++)
        for (c = 0; c < 3; c++)
            e->P[r][c] *= F[r] * F[c];
    for (r = 0; r < 3; r++)
        e->P[r][r] += QD[r];

    yhat = ocv(e->x[0], &slope) - e->x[1] - e->x[2] * i;
    H[0] = slope;
    H[1] = -1.0f;
    H[2] = -i;
    for (r = 0; r < 3; r++)
        PH[r] = e->P[r][0] * H[0] + e->P[r][1] * H[1] + e->P[r][2] * H[2];
    S = H[0] * PH[0] + H[1] * PH[1] + H[2] * PH[2] + R_V;
    y = v - yhat;
    for (r = 0; r < 3; r++) {
        K[r] = PH[r] / S;
        e->x[r] += K[r] * y;
    }
    for (r = 0; r < 3; r++)                 /* P = (I − K·H)·P, P симметрична: H·P = PHᵀ */
        for (c = 0; c < 3; c++)
            e->P[r][c] -= K[r] * PH[c];
}

#ifdef ML_HOST
#include <math.h>
#include <stdio.h>
int main(void)
{
    ekf_t e;
    int n;
    float worst = 0.0f, worst_r = 0.0f;
    ekf_init(&e, S0);
    for (n = 0; n < N_SAMPLES; n++) {
        ekf_step(&e, MEAS_V[n], MEAS_I[n]);
        if (fabsf(e.x[0] - REF_SOC[n]) > worst)
            worst = fabsf(e.x[0] - REF_SOC[n]);
        if (fabsf(e.x[2] - REF_R0[n]) > worst_r)
            worst_r = fabsf(e.x[2] - REF_R0[n]);
    }
    printf("steps=%d\nworst_soc_ppm=%d\nworst_r0_uohm=%d\nsoc_end=%.4f\nr0_end=%.4f\n", N_SAMPLES,
           (int)(worst * 1e6f + 0.5f), (int)(worst_r * 1e6f + 0.5f), e.x[0], e.x[2]);
    return 0;
}
#else
volatile float g_v, g_i, g_soc;
void SystemInit(void) {}
int main(void)
{
    ekf_t e;
    ekf_init(&e, S0);
    for (;;) {
        ekf_step(&e, g_v, g_i);
        g_soc = e.x[0];
    }
}
#endif

Дальше: раздел 121а «OpenMV: камера с машинным зрением на STM32H7 — DCMI, MicroPython и поиск объектов» — открытая камера на STM32H743: кадр через DCMI и MDMA, пороги в LAB и поиск пятен find_blobs, нейросети через TensorFlow Lite Micro.

Открыть в интерактивном справочнике