第四章 模型預測控制器
4.2 架構實現方法
4.2 架構實現方法
由於 MPC 為 model-based 控制器，所以在此我們使用 state-space model 去表 示一 single-input single-output 系統，如下式所示：
𝑥
𝑚
(𝑘 + 1) = 𝐴𝑚
𝑥𝑚
(𝑘) + 𝐵𝑚
𝑢(𝑘) (4-1) 𝑦(𝑘) = 𝐶𝑚
𝑥𝑚
(𝑘) (4-2)其中u 為輸入，y 為輸出，𝑥
𝑚
為 state variable (dimension : 𝑛1
)。對(4-1)式兩邊進行微分可得:
𝑥
𝑚
(𝑘 + 1) − 𝑥𝑚
(𝑘) = 𝐴𝑚
[𝑥𝑚
(𝑘) − 𝑥𝑚
(𝑘 − 1)] + 𝐵𝑚
[𝑢(𝑘) − 𝑢(𝑘 − 1)](4-3) 在此定義∆𝑥
𝑚
(𝑘) = 𝑥𝑚
(𝑘) − 𝑥𝑚
(𝑘 − 1) ; ∆u(k) = 𝑢(𝑘) − 𝑢(𝑘 − 1)上式(4-3)可以表示為:
∆𝑥
𝑚
(𝑘 + 1) = 𝐴𝑚
∆𝑥𝑚
(𝑘) + 𝐵𝑚
∆𝑢(𝑘) (4-4) 此時 state-space model 的表示式之輸入為∆u(k)由於 MPC 在決定新的 state 時必須將系統輸出納入考慮，所以定義一新的 state 變數 𝑥(𝑘) = [∆𝑥
𝑚
(𝑘)𝑇
𝑦(𝑘)]𝑇
，其中上標 T 代表矩陣的轉置(transpose)運算 符號。接著將式子(4-2)兩邊做相同的微分運算:
y(k + 1) − y(k) = 𝐶
𝑚
(𝑥𝑚
(𝑘 + 1) − 𝑥𝑚
(𝑘)) = 𝐶𝑚
∆𝑥𝑚
(𝑘 + 1)= 𝐶
𝑚
𝐴𝑚
∆𝑥𝑚
(𝑘) + 𝐶𝑚
𝐵𝑚
∆𝑢(𝑘) ( 4 - 5 ) 將(4-4)與(4-5)整理在一起，可以得到系統新的 state-space model 表示式:[∆𝑥
𝑚
(𝑘 + 1)𝑦(𝑘 + 1) ] = [ 𝐴
𝑚
𝑜𝑚 𝑇
𝐶
𝑚
𝐴𝑚
1] [∆𝑥𝑚
(𝑘)𝑦(𝑘) ] + [ 𝐵
𝑚
𝐶
𝑚
𝐵𝑚
] ∆𝑢(𝑘) = 𝑨𝑥(𝑘) + 𝑩∆𝑢(𝑘)43
y(k) = [𝑜
𝑚
1] [∆𝑥𝑚
(𝑘) 𝑦(𝑘) ]= 𝑪𝑥(𝑘) ( 4 - 6 ) 其中𝑜
𝑚
= [0 0 0 . . . 0⏞𝑛
1]
將系統的 state-space model 轉換為符合 MPC 需求的表示方法後，現在要做 的是預測出系統未來的輸出響應，要預測出未來的系統響應，首先必須先決定未 來的控制訊號，這些控制訊號即為可調整的變數。由於 MPC 在每個時間點都會 對未來控制訊號與系統響應進行估測，在離散系統中即為每個取樣時間發生的時 刻。在此定義目前取樣的時間點為𝑘
𝑖
, 𝑘𝑖
> 0，用來進行最佳化運算的響應預測 窗口長度為𝑁𝑝
，並假設在此刻系統的狀態變數𝑥(𝑘𝑖
)是可以被量測到的。而我們 期望得到的最佳控制訊號串：∆𝑢(𝑘
𝑖
), ∆𝑢(𝑘𝑖
+ 1), … , ∆𝑢(𝑘𝑖
+ 𝑁𝑐
− 1) 𝑁𝑐
為期望預測的控制訊號長度(control horizon)。而上述提到之𝑁
𝑝
，為由目前狀態𝑥(𝑘𝑖
)預測未來狀態的訊號長度(prediction horizon)，如下式表示 :𝑥(𝑘
𝑖
+ 1|𝑘𝑖
), 𝑥(𝑘𝑖
+ 2|𝑘𝑖
), 𝑥(𝑘𝑖
+ 3|𝑘𝑖
), … , 𝑥(𝑘𝑖
+ m|𝑘𝑖
), … , 𝑥(𝑘𝑖
+ 𝑁𝑝
|𝑘𝑖
) 一般而言，𝑁𝑐
會選擇小於或等於𝑁𝑝
。從系統的 state-space model (A,B,C)，可以將未來的狀態一一算出：
x(𝑘
𝑖
+ 1|𝑘𝑖
) = 𝐴𝑥(𝑘𝑖
) + 𝐵∆𝑢(𝑘𝑖
) (4-7) x(𝑘𝑖
+ 2|𝑘𝑖
) = 𝐴𝑥(𝑘𝑖
+ 1|𝑘𝑖
) + 𝐵∆𝑢(𝑘𝑖
+ 1)= 𝐴
2
𝑥(𝑘𝑖
) + 𝐴𝐵∆𝑢(𝑘𝑖
) + 𝐵∆𝑢(𝑘𝑖
+ 1) x(𝑘𝑖
+ 3|𝑘𝑖
) = 𝐴𝑥(𝑘𝑖
+ 2|𝑘𝑖
) + 𝐵∆𝑢(𝑘𝑖
+ 2)= 𝐴
3
𝑥(𝑘𝑖
) + 𝐴2
𝐵∆𝑢(𝑘𝑖
) + 𝐴𝐵∆𝑢(𝑘𝑖
+ 1) + 𝐵∆𝑢(𝑘𝑖
+ 2)⋮
44
x(𝑘
𝑖
+ 𝑁𝑝
|𝑘𝑖
) = 𝐴𝑁
𝑝𝑥(𝑘𝑖
) + 𝐴𝑁
𝑝−1
𝐵∆𝑢(𝑘𝑖
) + 𝐴𝑁
𝑝−2
𝐵∆𝑢(𝑘𝑖
+ 1) + ⋯ + 𝐴𝑁
𝑝−𝑁
𝑐𝐵∆𝑢(𝑘𝑖
+ 𝑁𝑐
− 1)(4-8) 有了預測的 state 後，我們就可以得到預測的系統輸出響應：
y(𝑘
𝑖
+ 1|𝑘𝑖
) = 𝐶𝐴𝑥(𝑘𝑖
) + 𝐶𝐵∆𝑢(𝑘𝑖
)y(𝑘
𝑖
+ 2|𝑘𝑖
) = 𝐶𝐴2
𝑥(𝑘𝑖
) + 𝐶𝐴𝐵∆𝑢(𝑘𝑖
) + 𝐶𝐵∆𝑢(𝑘𝑖
+ 1)y(𝑘
𝑖
+ 3|𝑘𝑖
) = 𝐶𝐴3
𝑥(𝑘𝑖
) + 𝐶𝐴2
𝐵∆𝑢(𝑘𝑖
) + 𝐶𝐴𝐵∆𝑢(𝑘𝑖
+ 1) + 𝐶𝐵∆𝑢(𝑘𝑖
+ 2)⋮
y(𝑘
𝑖
+ 𝑁𝑝
|𝑘𝑖
) = C𝐴𝑁
𝑝𝑥(𝑘𝑖
) + 𝐶𝐴𝑁
𝑝−1
𝐵∆𝑢(𝑘𝑖
) + 𝐶𝐴𝑁
𝑝−2
𝐵∆𝑢(𝑘𝑖
+ 1) + ⋯ + C𝐴𝑁
𝑝−𝑁
𝑐𝐵∆𝑢(𝑘𝑖
+ 𝑁𝑐
− 1)(4-9)
值得注意的是，式(4-9)中，預測訊號的輸入訊號皆為𝑥(𝑘
𝑖
)與∆𝑢(𝑘𝑖
)，所以在此 我們再定義Y 與∆U：Y = [𝑦(𝑘
𝑖
+ 1|𝑘𝑖
) 𝑦(𝑘𝑖
+ 2|𝑘𝑖
) 𝑦(𝑘𝑖
+ 3|𝑘𝑖
) ⋯ 𝑦(𝑘𝑖
+ 𝑁𝑝
|𝑘𝑖
)]𝑇
∆U = [∆𝑢(𝑘
𝑖
) ∆𝑢(𝑘𝑖
+ 1) ∆𝑢(𝑘𝑖
+ 2) ⋯ ∆𝑢(𝑘𝑖
+ 𝑁𝑐
− 1)]𝑇
此時便可將(4-9)整理如下：Y = 𝐅 𝑥(k
i
) + 𝚽 ΔU (4-10) 其中:F = [
CA C𝐴
2
C𝐴3
⋮ C𝐴
𝑁
𝑝]Φ = [
CB 0 0 ⋯ 0 𝐶𝐴𝐵 𝐶𝐵 0 ⋯ 0 𝐶𝐴
2
𝐵 𝐶𝐴𝐵 𝐶𝐵 ⋯ 0⋮
𝐶𝐴
𝑁
𝑝−1
𝐵 𝐶𝐴𝑁
𝑝−2
𝐵 𝐶𝐴𝑁
𝑝−3
𝐵 ⋯ 𝐶𝐴𝑁
𝑝−𝑁
𝑐𝐵]45
接下來的目標為如何決定最佳控制訊號的方法。考慮一給定輸入訊號，並且 在預測範圍中為保持不變的常數 r(𝑘
𝑖
)，而此輸入訊號與系統預測的響應可構成 一差值函數(error function)，藉由將此函式最小化(minimizing)，可計算出控制訊 號列 ∆U，此時計算出的結果即為在此定義的最佳控制訊號列。令 𝑅
𝑠 𝑇
為包含著輸入訊號資訊的資料向量 :𝑅
𝑠 𝑇
= [1 1 ⋯ 1]⏞ 𝑟(𝑘𝑖
)𝑁
𝑝定義系統之成本函數(cost function) J：
J = (𝑅
𝑠
− 𝑌)𝑇
(𝑅𝑠
− 𝑌) + ∆𝑈𝑇
𝑅̅∆𝑈 (4-11) 其中𝑅̅ = 𝑟𝑤
𝐼𝑁
𝑐×𝑁
𝑐(𝑟𝑤
≥ 0)顯而易見的，成本函數 J 的前半部(𝑅
𝑠
− 𝑌)𝑇
(𝑅𝑠
− 𝑌)，為輸入訊號與輸出訊 號的差值；後半部∆𝑈𝑇
𝑅̅∆𝑈，則反映了是否將∆U的大小納入考慮，藉由調整𝑟𝑤
， 可改變∆U在成本函數中佔的比重。舉例來說，如果𝑟𝑤
= 0，則我們要考慮的就是 如何算出控制訊號，讓系統輸出與輸入的差值在最短時間內降至最低；然而當 𝑟𝑤
> 0時，代表著必須將控制訊號∆U的大小考慮進來，此時的狀況就變成不僅只 要使輸出與輸入響應差值降低，還必須考慮到∆U的變化，感覺就像是我們必須 謹慎地將許多因素考慮進來，以達到期望的目的。關於𝑟𝑤
的影響在接下來的模擬 與實作中會再次驗證。接著繼續介紹求出期望控制訊號列∆U之公式推導，將(4-10)代入(4-11)可得:
J = (𝑅
𝑠
− 𝐹𝑥(𝑘𝑖
))𝑇
(𝑅𝑠
− 𝐹𝑥(𝑘𝑖
)) − 2∆𝑈𝑇
Φ𝑇
(𝑅𝑠
− 𝐹𝑥(𝑘𝑖
)) + ∆𝑈𝑇
(Φ𝑇
Φ + 𝑅̅)∆𝑈 (4-12) 由成本函數 J 對∆𝑈 做一階微分 :𝜕𝐽
𝜕∆𝑈
= −2Φ𝑇
(𝑅𝑠
− 𝐹𝑥(𝑘𝑖
)) + ∆𝑈𝑇
(Φ𝑇
Φ + 𝑅̅)∆𝑈 (4-13) 由於我們期望此成本函數降至最低，所以令：46
𝜕𝐽
𝜕∆𝑈= 0
由此條件與(4-13)式，最後我們可以得到最佳的控制訊號列∆𝑈為：
∆U = (Φ
𝑇
Φ + 𝑅̅)−1
Φ𝑇
(𝑅̅𝑠
𝑟(𝑘𝑖
) − 𝐹𝑥(𝑘𝑖
)) (4-14)其中𝑅̅
𝑠
= [1 1 1 ⋯ 1]⏞𝑇
𝑁
𝑝在離散系統中的每個取樣時間時，皆會算出一組∆U，其中包含著𝑁
𝑐
個控制 訊號的變化量(∆𝑢(𝑘𝑖
) ∆𝑢(𝑘𝑖
+ 1) ∆𝑢(𝑘𝑖
+ 2) ⋯ ∆𝑢(𝑘𝑖
+ 𝑁𝑐
− 1))，但我們只將 第一個∆𝑢(𝑘𝑖
)拿來使用，誠如上述計畫工作的案例所敘，每個時間點皆執行第一 個小時的工作後，檢視成果與重新訂定計畫，此方法稱為移動時域控制(receding horizon control)。而由(3-14)式可觀察到，(Φ𝑇
Φ + 𝑅̅)−1
Φ𝑇
𝑅̅𝑠
𝑟(𝑘𝑖
) 項中，因為Φ為 只與系統 model 特性(C,A,B)有關之矩陣，所以如果輸入訊號為固定時，此項即 為非時變矩陣；−(Φ𝑇
Φ + 𝑅̅)−1
𝐹𝑥(𝑘𝑖
) 項則因為 F 為只與系統 model (C,A)有關 之矩陣，所以−(Φ𝑇
Φ + 𝑅̅)−1
𝐹亦為非時變矩陣，非時變矩陣的好處在於我們可 以先將其結果求出，如此在實作時就不必在每個取樣時間都重新計算一次。基於移動時域控制法則(RHC)，我們只執行第一個控制訊號，所以我們將
∆U 的第一個元素取出 :
∆𝑢(𝑘
𝑖
) = [1 0 0 ⋯ 0]⏞𝑁
𝑝(Φ
𝑇
Φ + 𝑅̅)−1
Φ𝑇
(𝑅̅𝑠
𝑟(𝑘𝑖
) − 𝐹𝑥(𝑘𝑖
))= 𝐾
𝑦
𝑟(𝑘𝑖
) − 𝐾𝑚𝑝𝑐
𝑥(𝑘𝑖
) (4-15) 其中𝐾𝑦
為(Φ𝑇
Φ + 𝑅̅)−1
Φ𝑇
𝑅̅𝑠
的第一個元素，𝐾𝑚𝑝𝑐
為(Φ𝑇
Φ + 𝑅̅)−1
𝐹矩陣的第一列。值得注意的是，式(3-15)其實為一閉迴路控制系統的標準表示式，新的控制訊號 由輸入訊號減去回授訊號得到，而回授訊號之增益為𝐾
𝑚𝑝𝑐
。將(4-15)代回(4-7)中，可得：
𝑥(𝑘
𝑖
+ 1|𝑘𝑖
) = 𝐴𝑥(𝑘𝑖
) − 𝐵𝐾𝑚𝑝𝑐
𝑥(𝑘𝑖
) + 𝐵𝐾𝑦
𝑟(𝑘𝑖
)= (𝐴 − 𝐵𝐾
𝑚𝑝𝑐
)𝑥(𝑘𝑖
) + 𝐵𝐾𝑦
𝑟(𝑘𝑖
)(4-16)
47
把這個結果演繹到一般的情況中，將𝑘
𝑖
換成 k 可得：𝑥(𝑘 + 1) = (𝐴 − 𝐵𝐾
𝑚𝑝𝑐
)𝑥(k) + 𝐵𝐾𝑦
𝑟(k) (4-17) 而由於矩陣 F 的最後一欄與𝑅̅𝑠
相同，皆為[1 1 1 ⋯ 1]𝑇
，因此𝐾𝑦
與𝐾𝑚𝑝𝑐
的最後 一個元素相同，藉此，我們可以將𝐾𝑚𝑝𝑐
重新表示為:𝐾
𝑚𝑝𝑐
= [𝐾𝑥
𝐾𝑦
] 𝐾𝑥
為∆𝑥𝑚
(k)之回授增益。最後，由(4-17)，閉迴路控制器架構如圖 4-2 所示：
圖 4-2 離散系統中之 MPC 控制架構
其中𝑞
−1
表示延遲一個取樣時間；1−𝑞 1
−1表示離散系統中之積分控制項，相反的，1 − 𝑞
−1
表示離散系統中之微分控制項。實現的流程細節如圖 4-3 所示：
48
圖 4-3 MPC 控制訊號計算方法實現流程圖
① : 首先我們必須得到系統初始之控制訊號 u(k)，與狀態量測值
② : 有了目前的狀態後，我們可以得到y(k)，由(3-2)反推可得到𝑥
𝑚
(𝑘)③ : 將𝑥
𝑚
(𝑘)與控制訊號u(k)，由式(3-1)可算出新的系統狀態𝑥𝑚
(𝑘 + 1)④ : 由前一個 state 與新算出的𝑥
𝑚
(𝑘 + 1)，可得∆𝑥𝑚
(𝑘 + 1)與y(k + 1)⑤ : 由新的系統狀態，可對新的控制訊號∆U進行一定長度(𝑁
𝑝
)的預估計算⑥ : 但是我們只拿出第一個控制訊號變量∆u(k + 1|k)
⑦ : 計算出新的控制訊號後，重複此流程