動(dòng)態(tài)系統(tǒng)建模被各領(lǐng)域廣泛應(yīng)用，例如電動(dòng)汽車，能源系統(tǒng)，航空航天。我們本文提到動(dòng)態(tài)系統(tǒng)主要是被控對(duì)象，對(duì)被控對(duì)象進(jìn)行建模是因?yàn)槲覀兿Ｍ私膺@個(gè)系統(tǒng)（被控對(duì)象）的物理特性以及接受一些外部輸入（力，扭矩，電流等等）時(shí)會(huì)有什么樣的動(dòng)態(tài)響應(yīng)，基于此從而可以更好的給出控制輸入得到我們期望的系統(tǒng)的輸出，以及理解系統(tǒng)的退化或最大化提升系統(tǒng)效率。

這些動(dòng)態(tài)系統(tǒng)的行為是由多物理場(chǎng)復(fù)雜的交互作用決定的，因此系統(tǒng)行為和系統(tǒng)響應(yīng)建模通常需要復(fù)雜的第一原理支撐，仿真時(shí)也需要大量的計(jì)算(例如有限元模型)。

這也是本文的出發(fā)點(diǎn)，提供數(shù)據(jù)驅(qū)動(dòng)（主要介紹深度學(xué)習(xí)和系統(tǒng)辨識(shí)）的模型降階（Reduced Order Modeling）提速的方法，通過數(shù)據(jù)得到具有一定保真度的數(shù)據(jù)模型，在捕捉到系統(tǒng)動(dòng)態(tài)特性的同時(shí)也提升仿真速度。

本文中將涉及多個(gè) demo，數(shù)據(jù)以及腳本文件，若您感興趣進(jìn)一步獲取這些鏈接，可以在文末填寫反饋問卷，獲取這些鏈接。

動(dòng)態(tài)系統(tǒng)

動(dòng)態(tài)系統(tǒng)包含狀態(tài)空間 S, 時(shí)間集 T 和一個(gè)映射（規(guī)則）來描述狀態(tài)隨時(shí)間的演變規(guī)則 R:S×T→S。例如給定一個(gè)時(shí)刻 t 的狀態(tài) st ,通過這個(gè)規(guī)則可以計(jì)算后面一個(gè)或幾個(gè)時(shí)刻狀態(tài) st+1=R(st)，st+2=R(R(st)) 等等。動(dòng)態(tài)系統(tǒng)通?？梢杂秒S時(shí)間變化的方程或方程組來描述。尤其對(duì)于連續(xù)時(shí)間系統(tǒng)，可以通過微分方程來表示。

我們先看一個(gè)簡單的常微分系統(tǒng)（ODE）, [鏈接1]

其中 y(t) 是系統(tǒng)狀態(tài)。例如：一個(gè)簡單二自由度線性系統(tǒng)，

其中 A 是一個(gè) 2x2 的矩陣。初始條件 x0 = [2; 0], 可以通過求解 ODE 方程得到的相應(yīng)的二維輸出 x(t)，包含兩個(gè)狀態(tài)，時(shí)序和動(dòng)態(tài)圖如下：

圖表 1 系統(tǒng)輸出x(t):

x0 = [2; 0];

A = [-0.1 -1; 1 -0.1];

trueModel = @(t,y) A*y; % 定義系統(tǒng)函數(shù)，此處就是一個(gè)狀態(tài)空間方程

numTimeSteps = 2000;

T = 15;

odeOptions = odeset(RelTol=1.e-7);

t = linspace(0, T, numTimeSteps);

[~, xTrain] = ode45(trueModel, t, x0, odeOptions);

后面（在介紹 Neural ODE 部分）我們會(huì)嘗試?yán)眠@個(gè)系統(tǒng)的數(shù)據(jù) xTrain 進(jìn)行深度學(xué)習(xí)模型的訓(xùn)練來得到這個(gè)系統(tǒng)的數(shù)據(jù)代理模型（Surrogate Model），這種思路也可以同樣用于復(fù)雜系統(tǒng)。

既然數(shù)據(jù)驅(qū)動(dòng)，有很多機(jī)器學(xué)習(xí)和深度學(xué)習(xí)算法可以用，那動(dòng)態(tài)系統(tǒng)建模有什么特殊性呢，不是都適用嗎？

目前工程中已經(jīng)用到很多穩(wěn)態(tài)（靜態(tài)）模型。例如在發(fā)動(dòng)機(jī)排放標(biāo)定，通過 DoE 試驗(yàn)時(shí)我們會(huì)將發(fā)動(dòng)機(jī)維持在不同的穩(wěn)態(tài)工況（轉(zhuǎn)速恒定，扭矩恒定等等），通過試驗(yàn)數(shù)據(jù)建模得到用于標(biāo)定的穩(wěn)態(tài)數(shù)據(jù)模型[鏈接2]。

穩(wěn)態(tài)工況下，對(duì)于方程（1）這樣一個(gè)簡單系統(tǒng)，其中 y(t)' 可以看作 0，到達(dá)平衡點(diǎn)，于是 y(t) 和 u(t) 關(guān)系恒定，不再在時(shí)間維度上與歷史狀態(tài) y(t-1)，y(t-2) 等等有關(guān)，因此穩(wěn)態(tài)模型針對(duì)穩(wěn)態(tài)工況是非常準(zhǔn)的。

而在瞬態(tài)工況下通常 y(t)' 非零，因此方程（1）在求解系統(tǒng)輸出 y(t) 時(shí)不僅由當(dāng)前時(shí)刻的輸入 u(t) 決定，還取決于 y(t-1)，y(t-2)，u(t-1), u(t-2) 等等，這就是動(dòng)態(tài)系統(tǒng)的特殊性，當(dāng)前輸出不僅依賴于當(dāng)前的輸入，還依賴于系統(tǒng)過去的行為（歷史輸入和歷史輸出）。我們?cè)谙乱黄P(guān)于系統(tǒng)辨識(shí)的文章中會(huì)基于示例詳細(xì)說明這一點(diǎn)。

不考慮動(dòng)態(tài)系統(tǒng)，單純從系統(tǒng)建模來說，通常有如下的兩個(gè)方向：基于第一原理的和基于數(shù)據(jù)驅(qū)動(dòng)的。

第一原理建模是領(lǐng)域工程師都比較熟悉的，例如可以使用 M 腳本語言，Simulink 或 Simscape 建模語言從物理原理進(jìn)行系統(tǒng)模型的創(chuàng)建。

那什么情況下搭建系統(tǒng)會(huì)考慮使用或部分使用數(shù)據(jù)代理模型（Surrogate Model）？

• 物理系統(tǒng)原理比較復(fù)雜或者不夠清晰，無法構(gòu)建第一原理模型

• 數(shù)據(jù)獲取相對(duì)簡單

• 第一原理物理模型過于復(fù)雜，求解耗時(shí)，在控制開發(fā)時(shí)有時(shí)效要求，需要加速仿真計(jì)算

可以看到右半部分主要是基于數(shù)據(jù)驅(qū)動(dòng)的建模手段，其中針對(duì)動(dòng)態(tài)系統(tǒng)的建模主要是系統(tǒng)辨識(shí)和神經(jīng)網(wǎng)絡(luò)/深度學(xué)習(xí)。也是本系列兩個(gè)方向。

本文介紹神經(jīng)網(wǎng)絡(luò)的幾種用于動(dòng)態(tài)系統(tǒng)建模的模型，下一篇文章會(huì)介紹系統(tǒng)辨識(shí)的幾種模型。

前饋神經(jīng)網(wǎng)絡(luò)（Feedforward Neural Network）

前神經(jīng)網(wǎng)絡(luò)大家都不陌生，各神經(jīng)元分層排列。如圖2

圖2 簡單前饋神經(jīng)網(wǎng)絡(luò)示意圖

每層神經(jīng)元只接受上層輸出，結(jié)果只傳給下一層，沒有反饋。稍微復(fù)雜點(diǎn)的如圖3。

圖3 squeeze net MATLAB 示例

前饋神經(jīng)網(wǎng)絡(luò)是相對(duì)于循環(huán)神經(jīng)網(wǎng)絡(luò)（Recurrent Neural Network）而言，后者具有反饋，后面我們也會(huì)介紹。

那么問題來了，前饋神經(jīng)網(wǎng)絡(luò)因?yàn)闆]有反饋，如何表達(dá)動(dòng)態(tài)系統(tǒng)（Dynamic System）的時(shí)間狀態(tài)依賴呢？

通常主要手段使用不同尺度的滑動(dòng)窗口來構(gòu)建衍生特征從而表征系統(tǒng)在時(shí)間上的動(dòng)態(tài)。

電池 SoC 預(yù)測(cè)示例

我們以電池管理系統(tǒng)為例，通過使用深度學(xué)習(xí)來估計(jì)電池的荷電狀態(tài) SoC（State of Charge）。

本示例主要介紹 SoC 數(shù)據(jù)驅(qū)動(dòng)的建模方法，在即使不清楚電池電化學(xué)模型以及物理非線性特性的情況下，依然可以進(jìn)行 SoC 估計(jì)。我們可以通過實(shí)驗(yàn)室中收集到的實(shí)測(cè)數(shù)據(jù)進(jìn)行一個(gè)前饋神經(jīng)網(wǎng)絡(luò)代理模型訓(xùn)練實(shí)現(xiàn) SoC 估計(jì)[1]。

腳本和數(shù)據(jù)[鏈接3]

圖表4 數(shù)據(jù)集預(yù)覽和模型的5個(gè)輸入1個(gè)輸出

數(shù)據(jù)準(zhǔn)備

我們利用實(shí)測(cè)數(shù)據(jù)通過預(yù)處理得到訓(xùn)練數(shù)據(jù)集 cdsTrain，我們看其中一條數(shù)據(jù)記錄：

preview(cdsTrain)

ans = 1×2 cell array

{5×669956 double}    {1×669956 double}

數(shù)據(jù)包含五個(gè)特征，分別是電壓 V、電流 I、溫度 T、滑動(dòng)平均電壓 V_avg、滑動(dòng)平均電流 I _avg,其中后面兩個(gè)衍生特征用于表征動(dòng)態(tài)信息。

搭建模型

接下來我們構(gòu)建神經(jīng)網(wǎng)絡(luò)模型，結(jié)構(gòu)比較簡單，三層全連接網(wǎng)絡(luò)，一個(gè)輸出的 Feedforward Neural Network模型。

layers =[sequenceInputLayer(numFeatures,"Normalization","zerocenter")

    fullyConnectedLayer(numHiddenUnits)

    tanhLayer                            % HyperbolicTangent

    fullyConnectedLayer(numHiddenUnits)

    leakyReluLayer(0.3)                  % 激活函數(shù)

    fullyConnectedLayer(numResponses)

    clippedReluLayer(1)                  % 激活函數(shù)

regressionLayer];

設(shè)置訓(xùn)練選項(xiàng)并訓(xùn)練

options =trainingOptions('adam', ...  % Adam optimizer

'MaxEpochs',Epochs,...

'ExecutionEnvironment','cpu',...%可以選擇GPU

    'InitialLearnRate',InitialLearnRate, ...

    'LearnRateSchedule','piecewise', ...

    'LearnRateDropPeriod',LearnRateDropPeriod,...

    'LearnRateDropFactor',LearnRateDropFactor,...

    'ValidationData', {X,Y}, ...

    'ValidationFrequency',validationFrequency,...

    'MiniBatchSize',miniBatchSize, ...

    'CheckpointPath', NET_Path);

圖表 5 訓(xùn)練過程 Loss 變化

導(dǎo)入測(cè)試數(shù)據(jù)驗(yàn)證模型在測(cè)試集上的準(zhǔn)確度

圖表 6 深度學(xué)習(xí)模型的預(yù)測(cè)值與實(shí)測(cè)值比較