[MPC] MPC를 위한 두가지 QP 모델 - 2

MPC(model predictive control) 문제를 최적화 문제인 QP(quadratic program)문제로 변환할 때 널리 사용되는 일반적인 방법( https://pasus.tistory.com/229 )은 모델이 조밀해져서 문제의 구조가 손실되는 단점이 있다.

 

고속 MPC(fast MPC)에서는 QP 문제로 변환할 시 적절한 변수 재정렬을 사용하여 희소(sparse)구조를 최대로 이용할 수 있도록 특별한 방법으로 변환하며, 최적화 기법을 적용할 시 warm start, fixed iteration, early termination등의 휴리스틱 기법을 이용하여 MPC 계산량을 대폭 줄이는 방법을 사용한다. fast MPC알고리즘은 스탠퍼드의 Boyd 교수와 그 제자의 논문인 'Fast Model Predictive Control Using Online Optimization'에 자세히 설명되어 있다. 여기서는 fast MPC에서 사용된 QP 변환 모델을 유도해 보고자 한다.

MPC는 다음과 같은 제약조건을 갖는 선형 시스템에서

 

$$\begin{array}{r}\text{(1)}& & {\mathbf{x}}_{t+1}=A{\mathbf{x}}_t+B{\mathbf{u}}_t\\ & {\mathbf{y}}_t=C{\mathbf{x}}_t\\ \\ \text{(2)}& & {\mathbf{u}}_{min}\le {\mathbf{u}}_{t+i}\le {\mathbf{u}}_{max},i=0,...,N-1\\ & {\mathbf{y}}_{min}\le {\mathbf{y}}_{t+i}\le {\mathbf{y}}_{max},i=1,...,N\end{array}$$

 

매 시간 스텝마다 다음 목적함수가 일정 성능 예측구간 $[t,t+N]$ 에서 최적화되도록 제어입력의 시퀀스 ${\mathbf{u}}_{t:t+N-1}^{\ast }=({\mathbf{u}}_t^{\ast },{\mathbf{u}}_{t+1}^{\ast },...,{\mathbf{u}}_{t+N-1}^{\ast })$ 를 계산하는 문제다. 여기서 ${\mathbf{x}}_t\in {\mathbb{R}}^n$, ${\mathbf{u}}_t\in {\mathbb{R}}^p$ 이다.

 

$$\begin{array}{}\text{(3)}& \underset{{\mathbf{u}}_{t:t+N-1}}{min}J={\mathbf{x}}_{t+N}^T{Q}_f{\mathbf{x}}_{t+N}+\sum _{i=0}^{N-1}{\mathbf{x}}_{t+i}^{T}Q{\mathbf{x}}_{t+i}+{\mathbf{u}}_{t+i}^{T}R{\mathbf{u}}_{t+i}\end{array}$$

 

식 (1)에 의하면 시간 스텝(time step) $t$ 부터 $t+N$ 까지 미래의 상태를 다음과 같이 계산할 수 있다.

 

$$\begin{array}{r}\text{(4)}& & -B{\mathbf{u}}_t+{\mathbf{x}}_{t+1}=A{\mathbf{x}}_t\\ & -A{\mathbf{x}}_{t+1}-B{\mathbf{u}}_{t+1}+{\mathbf{x}}_{t+2}=0\\ \\ & -A{\mathbf{x}}_{t+2}-B{\mathbf{u}}_{t+2}+{\mathbf{x}}_{t+3}=0\\ \\ & \cdots \\ \\ & -A{\mathbf{x}}_{t+N-1}-B{\mathbf{u}}_{t+N-1}+{\mathbf{x}}_{t+N}=0\end{array}$$

 

식 (2)의 제약조건 수식은 다음과 같이 재구성할 수 있다.

 

$$\begin{array}{r}\text{(5)}& & \left[\begin{array}{c}I\\ -I\end{array}\right]{\mathbf{u}}_{t+i}\le \left[\begin{array}{c}{\mathbf{u}}_{max}\\ -{\mathbf{u}}_{min}\end{array}\right],i=0,...,N-1\\ & \left[\begin{array}{c}C\\ -C\end{array}\right]{\mathbf{x}}_{t+i}\le \left[\begin{array}{c}{\mathbf{y}}_{max}\\ -{\mathbf{y}}_{min}\end{array}\right],i=1,...,N\end{array}$$

 

식 (3)의 목적함수도 다음과 같이 행렬 형식으로 표현할 수 있다.

 

$$\begin{array}{r}\text{(6)}& J& ={\mathbf{x}}_t^{T}Q{\mathbf{x}}_t+\left[\begin{array}{ccccc}{\mathbf{x}}_{t+1}& {\mathbf{x}}_{t+2}& \cdots & {\mathbf{x}}_{t+N-1}& {\mathbf{x}}_{t+N}\end{array}\right]\left[\begin{array}{ccccc}Q& & & & \\ & Q& & & \\ & & \ddots & & \\ & & & Q& \\ & & & & Q_f\end{array}\right]\left[\begin{array}{c}{\mathbf{x}}_{t+1}\\ {\mathbf{x}}_{t+2}\\ {\mathbf{x}}_{t+3}\\ ⋮\\ {\mathbf{x}}_{t+N}\end{array}\right]\\ & +\left[\begin{array}{ccccc}{\mathbf{u}}_t& {\mathbf{u}}_{t+1}& \cdots & {\mathbf{u}}_{t+N-2}& {\mathbf{u}}_{t+N-1}\end{array}\right]\left[\begin{array}{ccccc}R& & & & \\ & R& & & \\ & & \ddots & & \\ & & & R& \\ & & & & R\end{array}\right]\left[\begin{array}{c}{\mathbf{u}}_t\\ {\mathbf{u}}_{t+1}\\ {\mathbf{u}}_{t+2}\\ ⋮\\ {\mathbf{u}}_{t+N-1}\end{array}\right]\end{array}$$

 

이제 최적화 변수 $\mathbf{z}$ 를 다음과 같이 정의한다.

 

$$\begin{array}{}\text{(7)}& \mathbf{z}=\left[\begin{array}{c}{\mathbf{u}}_t\\ {\mathbf{x}}_{t+1}\\ {\mathbf{u}}_{t+1}\\ ⋮\\ {\mathbf{u}}_{t+N-1}\\ {\mathbf{x}}_{t+N}\end{array}\right]\in {\mathbb{R}}^{N(n+p)}\end{array}$$

 

그러면 식 (6)의 목적함수를 다음과 같이 최적화 변수의 2차함수로 바꿀 수 있다.

 

$$\begin{array}{}\text{(8)}& \tilde{J}={\mathbf{z}}^{T}H\mathbf{z}\end{array}$$

 

여기서 ${\mathbf{x}}_t^{T}Q{\mathbf{x}}_t$ 는 최적화 대상이 아니기 때문에 제외하였으며,

 

$$H=\left[\begin{array}{ccccccc}R& & & & & & \\ & Q& & & & & \\ & & R& & & & \\ & & & \ddots & & & \\ & & & & Q& & \\ & & & & & R& \\ & & & & & & Q_f\end{array}\right]\in {\mathbb{R}}^{N(n+p)\times N(n+p)}$$

 

이다. 식 (5)의 제약조건도 다음과 같이 최적화 변수의 부등식으로 바꿀 수 있다.

 

$$\begin{array}{}\text{(9)}& P\mathbf{z}\le \mathbf{h}\end{array}$$

 

여기서

 

$$P=\left[\begin{array}{cccc}I& & & \\ & C& & \\ & & \ddots & \\ & & & C\\ -I& & & \\ & -C& & \\ & & \ddots & \\ & & & -C\end{array}\right]\in {\mathbb{R}}^{2N(n+p)\times N(n+p)},\mathbf{h}=\left[\begin{array}{c}{\mathbf{u}}_{max}\\ {\mathbf{x}}_{max}\\ ⋮\\ {\mathbf{x}}_{max}\\ -{\mathbf{u}}_{min}\\ -{\mathbf{x}}_{min}\\ ⋮\\ -{\mathbf{x}}_{min}\end{array}\right]$$

 

이다.

 

 

시스템 방정식 (4)는 다음과 같이 최적화 변수의 등식으로 바꿀 수 있다.

 

$$\begin{array}{}\text{(10)}& D\mathbf{z}=\mathbf{b}\end{array}$$

 

여기서

 

$$\begin{array}{rl}& D=\left[\begin{array}{cccccccccc}-B& I& 0& 0& 0& 0& & & & 0\\ 0& -A& -B& I& 0& 0& & & & \\ & & & -A& -B& I& & & & & \\ & & & & & & \ddots & & & \\ & & & & & & & I& 0& 0\\ & & & & & & & -A& -B& I\end{array}\right]\in {\mathbb{R}}^{Nn\times N(n+p)},\\ \\ & \mathbf{b}=\left[\begin{array}{c}A{\mathbf{x}}_t\\ 0\\ 0\\ ⋮\\ 0\\ 0\end{array}\right]\in {\mathbb{R}}^{Nn}\end{array}$$

 

이다. 식 (8)-(10)에 의하면 식 (1)-(3)의 MPC 최적화 문제는 다음과 같이 QP 최적화 문제로 변환된다.

 

$$\begin{array}{r}\text{(11)}& & \underset{\mathbf{z}}{min}\tilde{J}={\mathbf{z}}^{T}H\mathbf{z}\\ & \text{subject to }P\mathbf{z}\le \mathbf{h},D\mathbf{z}=\mathbf{b}\end{array}$$

 

위 최적화 문제의 해는 ${\mathbf{x}}_t$ 값에 대한 ${\mathbf{z}}^{\ast }=({\mathbf{u}}_t^{\ast },{\mathbf{x}}_{t+1}^{\ast },{\mathbf{u}}_{t+1}^{\ast },...,{\mathbf{x}}_{t+N}^{\ast })$ 으로서 제어입력과 상태변수의 최적 시퀀스이다.

 

 

식 (11)의 QP 문제의 행렬과 벡터의 구조를 보면 널리 사용되는 일반적인 방법과는 달리 MPC 모델의 구조가 훼손되지 않았고 희소(sparse)행렬임을 알 수 있다. 따라서 범용 QP 솔버 대신 이를 최대로 이용할 수 있는 특화된 QP 솔버를 만들 수 있다.