최적 유도법칙 (Optimal Guidance Law): 최종 속도 설정

중력장에서 비행체 또는 미사일의 운동 방정식은 다음과 같이 주어진다.

 

\[ \begin{align} \dot{\mathbf{r}} &= \mathbf{v} \tag{1} \\ \\ \dot{\mathbf{v}} &= -\frac{\mu}{r^3} \mathbf{r}+ \mathbf{a} \\ \\ &= \mathbf{g}( \mathbf{r})+ \mathbf{a} \end{align} \]

 

여기서 \(\mathbf{r}\) 과 \(\mathbf{v}\) 는 각각 관성좌표계에 대한 위치벡터와 속도벡터를 나타낸다. \(\mathbf{a}\) 는 제어 가속도, \(\mu\) 는 중력파라미터, \(\mathbf{g}(\mathbf{r})\) 은 비행체 또는 미사일에 작용하는 중력 가속도로서 위치의 함수이다.

 

 

제어의 목적은 시간 \(t_0\) 에서 출발지를 출발하여 최소의 에너지를 사용하여 시간 \(t_f\) 에 목적지에 비행체를 도착시키는 것이다. 이 제어 목적을 실현시키기 위하여 최소화해야 할 비용함수와 제약조건은 다음과 같다.

 

\[ \begin{align} & J= \frac{1}{2} \int_{t_0}^{t_f} \mathbf{a}^T \mathbf{a} \ dt \tag{2} \\ \\ & \mathbf{r}(t_0 )=\mathbf{r}_0, \ \ \ \mathbf{r}(t_f )=\mathbf{r}_f \tag{3} \\ \\ & \mathbf{v}(t_0 )= \mathbf{v}_0, \ \ \ \mathbf{v}(t_f )= \mathbf{v}_f \end{align} \]

 

식 (1)~(3)으로 주어지는 최적제어 문제는 위성이나 적 비행체의 요격, 달 또는 지구로의 정밀 착륙, 지상 목표물의 타격 임무 등에 적용될 수 있는 유도법칙 설계 문제이기도 하다.

 

 

식 (1)에서 중력 가속도는 일반적으로 위치의 함수이지만 여기서는 해석적인 해를 유도하기 위하여 상수벡터로 가정하겠다. 즉, \(\mathbf{g}(\mathbf{r})= \mathbf{g}\). 이 경우 식 (1), (3)은 다음과 같이 상태변수 방정식으로 표현할 수 있다.

 

\[ \dot{\mathbf{x}} =A \mathbf{x}+B( \mathbf{a}+ \mathbf{g}) \tag{4} \]

여기서

\[ \begin{align} & A= \begin{bmatrix} 0 & I \\ 0 & 0 \end{bmatrix} , \ \ \ B= \begin{bmatrix} 0 \\ I \end{bmatrix} \\ \\ & \mathbf{x}= \begin{bmatrix} \mathbf{r} \\ \mathbf{v} \end{bmatrix} , \ \ \ \mathbf{x}(t_0 )= \begin{bmatrix} \mathbf{r}_0 \\ \mathbf{v}_0 \end{bmatrix}, \ \ \ \mathbf{x}(t_f )= \begin{bmatrix} \mathbf{r}_f \\ \mathbf{v}_f \end{bmatrix} \end{align} \]

 

이다. 식 (2)와 (4)로 주어지는 최적제어 문제는 예제(https://pasus.tistory.com/258, https://pasus.tistory.com/259)에서 풀었던 최종상태가 설정된 LQR 문제와 거의 동일하지만, 여기서 다시 자세히 풀어보도록 한다.

먼저 해밀토니안(Hamiltonian) 함수를 정의한다.

 

\[ \mathcal{H}= \frac{1}{2} (\mathbf{a}^T \mathbf{a})+ \lambda^T \left( A\mathbf{x}+B(\mathbf{a}+\mathbf{g}) \right) \tag{5} \]

 

그러면 다음과 같이 상태변수와 코스테이트(costate) 미분 방정식을 얻을 수 있다.

 

\[ \begin{align} \dot{\mathbf{x}} &= \frac{\partial \mathcal{H}}{\partial \lambda}=A \mathbf{x}+B(\mathbf{a}+\mathbf{g}) \tag{6} \\ \\ - \dot{\lambda} &= \frac{ \partial \mathcal{H}}{\partial \mathbf{x}} =A^T \lambda \end{align} \]

 

코스테이트 미분 방정식은 상태변수와 독립이므로 다음과 같이 쉽게 해를 구할 수 있다.

 

\[ \lambda (t)= e^{A^T (t_f-t) } \lambda (t_f) \tag{7} \]

 

정정조건(stationary condition)은 다음과 같다.

 

\[ 0= \frac{\partial \mathcal{H}}{ \partial \mathbf{a} }=B^T \lambda + \mathbf{a} \tag{8} \]

 

따라서 최적제어는 식 (8)과 (7)로부터 다음과 같이 계산할 수 있다.

 

\[ \begin{align} \mathbf{a}(t) &= -B^T \lambda (t) \tag{9} \\ \\ &= -B^T e^{A^T (t_f-t) } \lambda (t_f) \end{align} \]

 

식 (9)를 식 (4)에 대입하면 상태변수 미분 방정식은 다음과 같이 된다.

 

\[ \dot{\mathbf{x}} =A \mathbf{x}+B \left( \mathbf{g}-B^T e^{A^T (t_f-t) } \lambda (t_f ) \right) \tag{10} \]

 

위 식을 풀면 \(\mathbf{x}(t)\) 를 다음과 같이 계산할 수 있다.

 

\[ \begin{align} \mathbf{x}(t) &= e^{A(t-t_0 )} \mathbf{x}(t_0 )+ \int_{t_0}^t e^{A(t-\tau)} \ d \tau \ B \mathbf{g} \tag{11} \\ \\ & \ \ \ \ \ - \int_{t_0}^t e^{A(t-\tau)} BB^T e^{A^T (t_f-\tau) } \ d\tau \ \lambda (t_f) \end{align} \]

 

위 식에 의하면 \(t=t_f\) 에서 \(\mathbf{x}(t_f )\) 를 다음과 같이 구할 수 있다.

 

\[ \begin{align} \mathbf{x}(t_f) &= e^{A(t_f-t_0 )} \mathbf{x}(t_0 )+ \int_{t_0}^{t_f} e^{A(t_f-\tau)} \ d \tau \ B \mathbf{g} \tag{12} \\ \\ & \ \ \ \ \ - \int_{t_0}^{t_f} e^{A(t_f-\tau)} BB^T e^{A^T (t_f-\tau) } \ d\tau \ \lambda (t_f) \end{align} \]

 

이 값이 \(\mathbf{x}(t_f )= \mathbf{x}_f\) 가 되도록 \(\lambda (t_f)\) 를 계산해야 한다. 여기서 \(G(t)\) 를 다음과 같이 정의하면,

 

\[ G(t) = \int_t^{t_f} e^{A(t_f-\tau)} BB^T e^{A^T (t_f-\tau) } \ d\tau \tag{13} \]

 

\(\lambda (t_f)\) 를 다음과 같이 계산할 수 있다.

 

\[ \lambda (t_f )= -G^{-1} (t_0 ) \left( \mathbf{x}_f- \int_{t_0}^{t_f} e^{A(t_f-\tau)} d\tau B\mathbf{g}-e^{A(t_f-t_0 )} \mathbf{x}(t_0 ) \right) \tag{14} \]

 

마지막으로 식 (14)를 식 (9)에 대입하면 최적제어는

 

\[ \mathbf{a}(t)= B^T e^{A^T (t_f-t) } G^{-1} (t_0 ) \left( \mathbf{x}_f- \int_{t_0}^{t_f} e^{A(t_f-\tau)} d\tau B\mathbf{g}-e^{A(t_f-t_0 )} \mathbf{x}(t_0 ) \right) \tag{15} \]

 

이 된다. 이 최적제어는 현재의 상태변수와 무관하고 상태변수의 초기값 \(\mathbf{x}(t_0 )\) 와 최종값 \(\mathbf{x}_f\) 에만 영향을 받으므로 오픈루프(open-loop) 제어이다. 따라서 만약 어떤 이유로 시스템이 교란(perturbed)된다면 최종값에 도달하지 못할 수도 있다. 따라서 최적제어가 초기 상태변수 \(\mathbf{x}(t_0 )\) 대신에 현재의 상태변수 \(\mathbf{x}(t)\) 의 함수가 되도록 식을 조금 변형해 보자.

 

 

식 (12)는 시간 \(t_0\) 에서 상태변수 \(\mathbf{x}(t_0 )\) 가 주어졌을 때 \(\mathbf{x}(t_f )\) 를 계산하는 식이다. 이 식은 임의의 시간 \(t_0 \le t_f\) 에서 상태변수 \(\mathbf{x}(t_0 )\) 가 주어졌을 때에도 모두 성립하는 식이기 때문에, 시간 \(t \le t_f\) 에서 상태변수 \(\mathbf{x}(t)\) 가 주어졌을 때에도 성립해야 한다.

 

\[ \begin{align} \mathbf{x}(t_f) &= e^{A(t_f-t )} \mathbf{x}(t )+ \int_{t}^{t_f} e^{A(t_f-\tau)} \ d \tau \ B \mathbf{g} \tag{16} \\ \\ & \ \ \ \ \ - \int_{t}^{t_f} e^{A(t_f-\tau)} BB^T e^{A^T (t_f-\tau) } \ d\tau \ \lambda (t_f) \end{align} \]

 

따라서 식 (14)는 다음과 같이 수정된다.

 

\[ \lambda (t_f )= -G^{-1} (t ) \left( \mathbf{x}_f- \int_{t}^{t_f} e^{A(t_f-\tau)} d\tau B\mathbf{g}-e^{A(t_f-t )} \mathbf{x}(t ) \right) \tag{17} \]

 

최적제어도 식 (18)과 같이 수정되므로 최적제어는 이제 최종 상태변수와 함께 현재의 상태변수 \(\mathbf{x}(t)\) 의 피드백 제어가 된다.

 

\[ \mathbf{a}(t)= B^T e^{A^T (t_f-t) } G^{-1} (t ) \left( \mathbf{x}_f- \int_{t}^{t_f} e^{A(t_f-\tau)} d\tau B\mathbf{g}-e^{A(t_f-t )} \mathbf{x}(t ) \right) \tag{18} \]

 

이제 식 (18)을 구성하는 행렬을 차례로 계산해 보자. 먼저 행렬지수함수는 다음과 같이 계산된다.

 

\[ e^{A(t_f-t)} =I+A(t_f-t)= \begin{bmatrix} I & (t_f-t)I \\ 0 & I \end{bmatrix} \tag{19} \]

 

행렬 \(G(t)\) 는 정의와 식 (19)에 의해서 다음과 같이 계산된다.

 

\[ \begin{align} G(t) &= \int_t^{t_f} e^{A(t_f-\tau)} BB^T e^{A^T (t_f-\tau) } \ d \tau \tag{20} \\ \\ &= \int_t^{t_f} \begin{bmatrix} I & (t_f-\tau ) I \\ 0 & I \end{bmatrix} \begin{bmatrix} 0 \\ I \end{bmatrix} \begin{bmatrix} 0 & I \end{bmatrix} \begin{bmatrix} I & 0 \\ (t_f-\tau) I & I \end{bmatrix} \ d \tau \\ \\ &= \int_t^{t_f} \begin{bmatrix} (t_f-\tau)^2 I & (t_f-\tau) I \\ (t_f-\tau) I & I \end{bmatrix} \ d \tau \\ \\ &= \begin{bmatrix} -\frac{(t_f-\tau)^3 }{3} I & -\frac{(t_f-\tau)^2}{2} I \\ - \frac{(t_f-\tau)^2}{2} I & \tau I \end{bmatrix}_t^{t_f} \\ \\ &= \begin{bmatrix} \frac{t_{go}^3}{3} I & \frac{t_{go}^2 }{2} I \\ \frac{t_{go}^2}{2} I & t_{go} I \end{bmatrix} \end{align} \]

 

여기서 \(t_{go}=t_f-t\) 는 잔여시간(time-to-go)을 나타낸다. \(G(t)\) 의 역행렬은 다음과 같다.

 

\[ \begin{align} G^{-1} (t) &= \frac{12}{t_{go}^4 } \begin{bmatrix} t_{go} I & -\frac{t_{go}^2}{2} I \\ - \frac{t_{go}^2}{2} I & \frac{t_{go}^3}{3} I \end{bmatrix} \tag{21} \\ \\ &= \begin{bmatrix} \frac{12}{t_{go}^3} I & - \frac{6}{t_{go}^2} I \\ - \frac{6}{t_{go}^2 } I & \frac{4}{t_{go}} I \end{bmatrix} \end{align} \]

 

식 (19)와 (21)을 이용하면 식 (18)의 구성 부분을 다음과 같이 계산할 수 있다.

 

\[ \begin{align} & B^T e^{A^T (t_f-t) } G^{-1} (t)= \begin{bmatrix} \frac{6}{t_{go}^2} I & -\frac{2}{t_{go}} I \end{bmatrix} \tag{22} \\ \\ & \int_t^{t_f} e^{A(t_f-\tau)} \ d\tau \ B \mathbf{g}= \begin{bmatrix} \frac{t_{go}^2}{2} \mathbf{g} \\ t_{go} \mathbf{g} \end{bmatrix} \tag{23} \\ \\ & \mathbf{x}_f- \int_t^{t_f} e^{A(t_f-\tau)} \ d\tau B \mathbf{g}-e^{A(t_f-t)} \mathbf{x}(t) \tag{24} \\ \\ & \ \ \ \ \ \ \ = \begin{bmatrix} \mathbf{r}_f-\frac{t_{go}^2}{2} \mathbf{g}- \mathbf{r}(t)-t_{go} \mathbf{v}(t) \\ \mathbf{v}_f-t_{go} \mathbf{g}-\mathbf{v}(t) \end{bmatrix} \end{align} \]

 

식 (22), (23), (24)를 식 (18)에 대입하면 최적제어는 다음과 같다.

 

\[ \begin{align} \mathbf{a}(t) & = \frac{6}{t_{go}^2 } \mathbf{r}_f-3 \mathbf{g}- \frac{6}{t_{go}^2 } \mathbf{r}(t)- \frac{6}{t_{go}} \mathbf{v}(t) \tag{25} \\ & \ \ \ \ \ - \frac{2}{t_{go}} \mathbf{v}_f+2 \mathbf{g}+ \frac{2}{t_{go}} \mathbf{v}(t) \\ \\ &= \frac{6}{t_{go}^2 } \left( \mathbf{r}_f- \left( \mathbf{r}(t)+t_{go} \mathbf{v}(t)+ \frac{1}{2} t_{to}^2 \mathbf{g} \right) \right) \\ & \ \ \ \ \ -\frac{2}{t_{go}} \left( \mathbf{v}_f-(\mathbf{v}(t)+t_{go} \mathbf{g}) \right) \\ \\ &= \frac{6}{t_{go}^2} ZEM- \frac{2}{t_{go}} ZEV \end{align} \]

 

여기서 ZEM(Zero-Effort-Miss) 거리와 ZEV(Zero-Effort-Velocity)는 다음과 같이 정의된다.

 

\[ \begin{align} & ZEM = \mathbf{r}_f- \left( \mathbf{r}(t)+t_{go} \mathbf{v}(t)+ \frac{1}{2} t_{to}^2 \mathbf{g} \right) \tag{26} \\ & ZEV = \mathbf{v}_f-(\mathbf{v}(t)+t_{go} \mathbf{g}) \end{align} \]

 

정의에 의하면 ZEM과 ZEV는 현재 시간 \(t\) 에서 임무 종료 시간 \(t_f\) 까지 제어 가속도가 \(0\) 일 때, 즉 zero-effort 일 때 임무 종료 시에 예측되는 거리와 속도 오차를 뜻한다.

식 (25)를 최종 속도가 설정된 경우의 최적 피드백 유도법칙(optimal feedback guidance law)이라고 한다. 최적 피드백 유도법칙은 최종 시간이 설정되어 있으므로 미리 결정된 시간에 요격을 달성하는 데 필요한 가속도 명령을 계산한다.

최종 속도벡터 \(\mathbf{v}_f\) 에 대해서는 랑데부나 요격 등의 임무 조건에 따라 벡터의 방향과 크기가 요구조건으로 제시될 수 있다.

 

 

대부분의 종말 유도 단계에서는 중력 가속도를 상수라고 가정할 수 있지만 중력 가속도를 위치의 함수로 모델링해야 하는 경우도 있을 수 있다. 예를 들면 임무 비행시간 동안 비행체 또는 미사일의 운동에 상당한 고도 변화를 수반하는 경우를 들 수 있겠다. 이러한 경우에 제어 가속도가 중력 가속도를 직접적으로 보상할 수 있다는 가정하에 식 (25)의 최적 유도법칙을 수정할 수 있다. 즉, 식 (1)의 제어 가속도를 다음과 같이 분할한다면,

 

\[ \mathbf{a}= \mathbf{a}_1- \mathbf{g}(\mathbf{r}) \tag{27} \]

 

동일한 유도 방법을 통하여 식 (25)를 다음과 같이 수정할 수 있다.

 

\[ \mathbf{a}(t) = \frac{6}{t_{go}^2 } \left( \mathbf{r}_f- \left( \mathbf{r}(t)+t_{go} \mathbf{v}(t)\right) \right) -\frac{2}{t_{go}} \left( \mathbf{v}_f-\mathbf{v}(t) \right) -\mathbf{g}(\mathbf{r}) \tag{28} \]

 

식 (1)의 운동 방정식은 또한 일정한 속도로 비행하는 목표물과의 상대 운동 방정식과 형태가 동일하므로 식 (25)와 (27)의 유도법칙은 고정된 목적지 뿐만 아니라 일정한 속도로 비행하는 목표물을 요격하는 데도 사용할 수 있다. 그림은 미사일과 적 항공기간의 기하학적인 관계를 보여준다. \(\{i\}\) 는 관성좌표계를 나타내며, \(\mathbf{r}_m, \mathbf{r}_T\) 는 각각 미사일과 적 항공기의 위치벡터이고, \(\mathbf{v}_m, \mathbf{v}_T\) 는 각각 미사일과 적 항공기의 속도벡터이다. \(\mathbf{a}_m\) 은 미사일의 가속도 명령(command)이다. \(\mathbf{r}= \mathbf{r}_T- \mathbf{r}_m, \ \mathbf{v}= \mathbf{v}_T-\mathbf{v}_m\) 은 각각 미사일과 적 항공기간의 상대 거리벡터와 상대 속도벡터이다.

 

 

그러면 다음과 같은 운동 방정식이 성립한다.

 

\[ \begin{align} & \dot{\mathbf{r}}= \mathbf{v} \tag{29} \\ \\ & \dot{\mathbf{v}} = - \mathbf{g}- \mathbf{a}_m \end{align} \]

 

식 (29)와 (1)의 차이는 가속도 벡터의 부호가 서로 반대라는 것 밖에 없다.