[Continuous-Time] 밸런싱 변환을 이용한 모델 차원 축소

다음과 같은 선형 시불변(LTI) 시스템이 있다.

 

\[ \begin{align} \dot{\mathbf{x}} &= A \mathbf{x}+B \mathbf{u} \tag{1} \\ \\ \mathbf{y} &= C\mathbf{x}+D \mathbf{u} \end{align} \]

 

여기서 \(\mathbf{x}(t) \in \mathbb{R}^n\) 는 상태변수, \(\mathbf{u}(t) \in \mathbb{R}^p\) 는 제어입력, \(\mathbf{y}(t) \in \mathbb{R}^q\) 는 출력이다.

이 LTI 시스템은 제어가능(controllable)하고 관측가능(observable)하며 안정(stable)하다고 가정한다.

그러면 이 시스템의 무한 제어가능성 그래미안(infinite controllability gramian) \(W_c\) 와 무한 관측가능성 그래미안(infinite observability gramian) \(W_o\) 는 다음 식의 해다.

 

\[ \begin{align} & AW_c+W_c A^T+BB^T=0 \tag{2} \\ \\ & A^T W_o+W_o A+C^T C=0 \end{align} \]

 

여기서 \(W_c \gt 0\), \(W_o \gt 0\) 이다. 또한 앞선 논의에 의하면 상태변수 \(\mathbf{x}\) 까지 도달하기 위한 최소 제어 에너지와 상태변수 \(\mathbf{x}\) 가 산출하는 최대 관측 에너지는 다음과 같이 계산할 수 있다.

 

\[ \begin{align} & \min \lVert \mathbf{u}(t) \rVert^2 = \mathbf{x}^T W_c^{-1} \mathbf{x} \tag{3} \\ \\ & \max \lVert \mathbf{y}(t) \rVert^2 = \mathbf{x}^T W_o \mathbf{x} \end{align} \]

 

\( \min \lVert \mathbf{u}(t) \rVert^2\) 은 상태변수 \(\mathbf{x}\) 까지 도달하기 위해 필요한 최소의 제어 에너지이므로 이 값이 클 수록 상태 상태변수 \(\mathbf{x} \) 까지 도달하기가 어려워진다. 비슷하게 \( \max \lVert \mathbf{y}(t) \rVert^2 \) 이 작을 수록 상태변수 \(\mathbf{x}\) 는 관측하기 어려워진다.

모델 차원 축소(model reduction)에서 절삭(truncation) 방법은 상태공간(state-space) 모델에서 중요하지 않은 상태를 제거하거나 잘라내는 것이다. 이러한 관점에서 볼 때 입력 \(\mathbf{u}(t)\) 에서 출력 \(\mathbf{y}(t)\) 로의 에너지 전달에 가장 적게 관여하는 상태변수 벡터 \(\mathbf{x}(t)\) 의 구성 요소를 삭제하는 것도 한 방법이다.

제어 에너지를 이용한 모델 차원 축소 아이디어는 상태벡터의 구성 요소 중에서 상대적으로 제어 에너지를 크게 만드는 성분을 잘라내는 것이 될 수 있다. 동일한 제어 에너지로 다른 상태벡터 성분보다도 상대적으로 덜 움직인다면 입력과 출력 관계에서 덜 영향을 주기 때문이다. 마찬가지로 관측 에너지를 이용한 모델 차원 축소 아이디어는 상태벡터의 구성 요소 중에서 관측 에너지에 미치는 영향이 작은 성분을 제거하는 것이다.

 

 

그러면 이제 상태벡터의 구성 요소 중에서 어느 성분이 제어 에너지가 더 큰지 알아보자.

무한 제어가능성 그래미안 \(W_c\) 는 정정행렬이기 때문에 고유값과 고유벡터는 실수값을 가지며, 고유벡터들은 서로 직각이다. 또한 고유값은 모두 \(0\) 보다 크다 (https://pasus.tistory.com/14). 따라서 행렬 \(W_c\) 는 다음과 같이 분해할 수 있다.

 

\[ \begin{align} W_c=E_c \Lambda_c E_c^T \tag{4} \end{align} \]

 

여기서 \(\Lambda_c\) 는 \(W_c\) 의 고유값 \(\lambda_{ci}\) 를 성분으로 하는 대각행렬이고 \(E_c\) 는 그에 대응하는 고유벡터 \(\mathbf{v}_{ci}\) 를 열로 하는 고유벡터 행렬이다. 위 식에 의하면 \(W_c\) 의 역행렬의 고유값 분해는 다음과 같다.

 

\[ \begin{align} W_c^{-1}=E_c \Lambda_c^{-1} E_c^T \tag{5} \end{align} \]

 

따라서 \(W_c\) 와 \(W_c^{-1}\) 의 고유벡터는 동일하며 고유값만 서로 역의 관계를 갖는다. 식 (5)를 식 (3)에 대입하면 다음과 같이 된다.

 

\[ \begin{align} \min \lVert \mathbf{u}(t) \rVert^2 &= \mathbf{x}^T E_c \Lambda_c^{-1} E_c^T \mathbf{x} \tag{6} \\ \\ &= \mathbf{z}_c^T \Lambda_c^{-1} \mathbf{z}_c \\ \\ &= \frac{z_{c1}^2}{\lambda_{c1}} + \frac{z_{c2}^2}{\lambda_{c2}}+ \cdots + \frac{z_{cn}^2}{\lambda_{cn}} \end{align} \]

 

여기서 \(\mathbf{z}_c=E_c^T \mathbf{x}\) 이고 이를 구성 요소별로 쓰면 \(z_{ci}= \mathbf{v}_{ci}^T \mathbf{x}\) 가 되므로 \(z_{ci}\) 는 상태벡터 \(\mathbf{x}\) 의 \(i \)-번째 고유벡터 \(\mathbf{v}_{ci}\) 방향 성분에 해당한다.

만약 \(W_c\) 의 고유값을 큰값부터 작은값 순으로 정렬한다면, 즉 \(\lambda_{c1} \ge \lambda_{c2} \ge \cdots \ge \lambda_{cn} \gt 0\), 큰 제어 에너지를 요구하는 상태 또는 도달하기 어려운 상태는 상태벡터 구성 요소중에서 작은 고유값에 해당하는 고유벡터로 스팬(span)되는 부분공간에 위치한 성분이 된다. 따라서 이러한 성분을 잘라내는 내는 것이 제어 에너지 관점에서의 모델 축소 방법이다.

 

 

이번에는 상태벡터의 구성 요소 중에서 어느 성분이 관측 에너지가 더 작은지 알아보자. 마찬가지로 행렬 \(W_o\) 는 다음과 같이 분해할 수 있다.

 

\[ \begin{align} W_o=E_o \Lambda_o E_o^T \tag{7} \end{align} \]

 

여기서 \(\Lambda_o\) 는 \(W_o\) 의 고유값 \(\lambda_{oi}\) 를 성분으로 하는 대각행렬이고 \(E_o\) 는 그에 대응하는 고유벡터 \(\mathbf{v}_{oi}\) 를 열로 하는 고유벡터 행렬이다. 식 (7)을 식 (3)에 대입하면 다음과 같이 된다.

 

\[ \begin{align} \max \lVert \mathbf{y}(t) \rVert^2 &= \mathbf{x}^T E_o \Lambda_o E_o^T \mathbf{x} \tag{8} \\ \\ &= \mathbf{z}_o^T \Lambda_o \mathbf{z}_o \\ \\ &= \lambda_{o1} z_{o1}^2 + \lambda_{o2} z_{o2}^2+ \cdots + \lambda_{on} z_{on}^2 \end{align} \]

 

여기서 \(\mathbf{z}_o=E_o^T \mathbf{x}\) 이고 이를 구성 요소별로 쓰면 \(z_{oi}=\mathbf{v}_{oi}^T \mathbf{x}\) 가 되므로 \(z_{oi}\) 는 상태벡터 \(\mathbf{x}\) 의 \(i\)-번째 고유벡터 \(\mathbf{v}_{oi}\) 방향 성분에 해당한다.

만약 \(W_o\) 의 고유값을 큰값부터 작은값 순으로 정렬한다면, 즉 \(\lambda_{o1} \ge \lambda_{o2} \ge \cdots \ge \lambda_{on} \gt 0\), 작은 관측 에너지를 생성하는 상태는 상태벡터 구성 요소중에서 작은 고유값에 해당하는 고유벡터로 스팬(span)되는 부분공간에 위치한 성분이 된다. 따라서 이러한 성분을 잘라내는 내는 것이 관측 에너지 관점에서의 모델 축소 방법이다.

여기서 문제는 도달하기 어려운 상태가 최대 관측 에너지를 생성하고, 반대로 가장 작은 관측 에너지를 생성하는 상태가 가장 도달하기 쉬운 상태가 될 수도 있다는 데 있다. 이러한 시스템의 경우 어떤 상태를 잘라내야 하는지 알 수 없다.

하지만 만약 두 그래미안이 서로 같고 대각행렬이라면 도달하기 어려운 상태가 가장 작은 관측 에너지를 생성하는 상태가 되고, 반대로 가장 큰 관측에너지를 생성하는 상태가 가장 도달하기 위운 상태가 될 것이다. 이러한 시스템을 밸런스된 시스템(balanced system)이라고 한다. 즉 제어가능성과 관측가능성의 수준이 서로 균형을 이룬 시스템이라는 뜻이다.

 

 

그렇다면 주어진 시스템을 밸런스된 시스템으로 변환할 순 없을까. 여기에는 두가지 전제가 필요하다. 우선 두 시스템의 입출력 관계식은 동일해야 한다. 또한 원래 시스템과 변환된 시스템의 제어가능성과 관측가능성 특성은 동일해야 한다.

원래 시스템과 임의로 변환된 시스템의 전달함수(transfer function)를 계산해 보자. 먼저 변환행렬 \(T\) 를 이용하여 상태변수 \(\mathbf{x}\) 를 상태변수 \(\tilde{\mathbf{x}}\) 로 변환한다.

 

\[ \begin{align} \mathbf{x}=T \tilde{\mathbf{x}} \tag{9} \end{align} \]

 

그러면 식 (1)은 다음과 같이 변환된다.

 

\[ \begin{align} \dot{\tilde{\mathbf{x}}} & ̇= \tilde{A} \tilde{\mathbf{x}}+ \tilde{B} \mathbf{u} \tag{10} \\ \\ \mathbf{y} &= \tilde{C} \tilde{\mathbf{x}} +D \mathbf{u} \end{align} \]

 

여기서

 

\[ \begin{align} \tilde{A}=T^{-1} AT, \ \ \ \tilde{B}=T^{-1} B, \ \ \ \tilde{C}=CT \tag{11} \end{align} \]

 

이다. \(A\) 와 \(\tilde{A}\) 는 상사관계이므로 고유값은 서로 같다. 식 (1)의 전달함수는 다음과 같다.

 

\[ \begin{align} G(s)=C(sI-A)^{-1} B+D \tag{12} \end{align} \]

 

또한 식 (11)을 이용하면 식 (10)의 전달함수는 다음과 같이 된다.

 

\[ \begin{align} \tilde{G}(s) &= \tilde{C}(sI- \tilde{A} )^{-1} \tilde{B}+D \tag{13} \\ \\ &= CT(sI-T^{-1} AT)^{-1} T^{-1} B+D \\ \\ &= CT \left[ T^{-1} (sI-A) T \right]^{-1} T^{-1} B+D \\ \\ &= C(sI-A)^{-1} B+D \\ \\ &= G(s) \end{align} \]

 

식 (13)에 의하면 동일한 입출력 관계식에 대해 무수히 많은 상태공간 방정식이 대응하는 것을 알 수 있다. 여기서 전달함수에 대응하는 상태공간 방정식을 실현(realization)이라고 한다. 즉 동일한 시스템에 대해서 상태공간 방정식의 실현은 무수히 많다. 그 중 밸런스된 시스템으로 실현할 수도 있을 것이다. 또한 다른 상태공간 방정식으로 실현되 있더라도 변환행렬을 통해서 밸런스된 시스템으로 재 실현할 수도 있다.

식 (11)을 이용하면 원래 시스템 (1)과 변환된 시스템 (10)의 그래미안 행렬간의 관계를 계산할 수 있다. 식 (2)에 식 (11)을 대입하면 두 그래미안 행렬식은 다음과 같이 된다.

 

\[ \begin{align} & T \tilde{A} T^{-1} W_c+W_c T^{-T} \tilde{A}^T T^T+T \tilde{B} \tilde{B}^T T^T=0 \tag{14} \\ \\ & T^{-T} \tilde{A}^T T^T W_o+W_o T \tilde{A}T^{-1}+T^{-T} \tilde{C}^T \tilde{C}T^{-1}=0 \end{align} \]

 

위 식의 양변에 각각 \(T^{-1}\) 과 \(T^{-T}\), \(T^T\) 와 \(T\) 를 곱하면 다음과 같이 된다.

 

\[ \begin{align} & \tilde{A} \tilde{W}_c + \tilde{W}_c \tilde{A}^T+\tilde{B}\tilde{B}^T=0 \tag{15} \\ \\ & \tilde{A}^T \tilde{W}_o+\tilde{W}_o \tilde{A}+\tilde{C}^T \tilde{C}=0 \end{align} \]

 

여기서

 

\[ \begin{align} & T^{-1} W_c T^{-T}=\tilde{W}_c \tag{16} \\ \\ & T^T W_o T= \tilde{W}_o \end{align} \]

 

이다. 따라서 \(W_c \gt 0, \ W_o \gt 0\) 이면 \(\tilde{W}_c \gt 0, \ \tilde{W}_o \gt 0\) 이고 그 반대도 성립하므로 두 시스템의 제어가능성과 관측가능성 특성은 동일하다.

이제 문제는 밸런스된 시스템으로 변환하기 위한 변환행렬 \(T\) 를 어떻게 계산할 것인가 하는 것이다.

제어가능하고 관측가능하며 안정한 LTI 시스템이 밸런싱 됐다면 무한 그래미안 행렬은 다음과 같이 된다.

 

\[ \begin{align} \tilde{W}_c= \tilde{W}_o= \Sigma = diag ( \sigma_1, ... , \sigma_n ) \tag{17} \end{align} \]

 

여기서 \(\sigma_i \gt 0\) 은 밸런스된 시스템의 그래미안 행렬의 고유값이다. 식 (16)의 두 식을 곱하면 다음과 같이 된다.

 

\[ \begin{align} \tilde{W}_c \tilde{W}_o= \Sigma^2=T^{-1} W_c W_o T \tag{18} \end{align} \]

 

식 (18)에 의하면

 

\[ \begin{align} W_c W_o T=T \Sigma^2 \tag{19} \end{align} \]

 

이므로 변환행렬 \(T\) 는 \(W_c W_o\) 의 고유벡터를 열로 하는 행렬이며 \(\Sigma^2\) 는 \(W_c W_o\) 의 고유값을 대각성분으로 하는 대각행렬이라는 것을 알 수 있다. 이 변환행렬 \(T\) 를 밸런싱 변환행렬 (balancing transformation matrix)이라고 하며 \(\sigma_i\) 를 시스템의 항켈 특이값 (Hankel singular value)이라고 한다.

앞선 논의에 의하면 제어와 관측이 동시에 어려운 상태는 이 항켈 특이값이 작은 값에 해당하므로 모델의 차원 축소를 위해서는 이 값이 작은 순서대로 제거하면 된다. 밸런스된 절삭(balanced truncation)에서는 큰 항켈 특이값에 해당하는 상태만 유지된다. 이러한 모델 축소 방법을 '밸런싱 방법을 이용한 모델 차원 축소 (balanced model reduction)' 방법이라고 한다.

밸런싱 방법을 이용하기 위해서는 여전히 원래 시스템의 제어가능성 및 관측가능성 그래미안을 계산해야 하는데, 축소하려는 시스템이 매우 고차적인 경우, 예를 들면 \(10^4\) 이상, 문제가 될 수 있기 때문에 이를 효율적으로 계산하는 방안이 필요하다.