상대 궤도운동 방정식 (Relative Orbit Equation of Motion)

우주공학의 미래라고 불리는 분산 우주시스템(distributed space system)은 단일 위성으로는 불가능한 임무를 수행하기 위해서 두 개 이상의 위성을 집단적으로 사용하는 시스템이다.

 

 

분산 우주시스템의 임무 개념의 예로서 궤도상(on-orbit) 서비스, 우주 상황 인식, 분산 군집(swarm) 기반 센싱, 위성 편대비행(formation flying), 랑데부 및 도킹 등을 들 수 있다.

분산 우주 시스템의 장점은 여러 위성 간의 상대 운동을 활용하는 데서 발생한다. 따라서 상대 운동을 표현하기 위한 좌표계와 기준 위성이 필요하다. 보통 분산 우주시스템의 임무가 지구를 중심으로 수행되므로 관성 좌표계로는 지구중심 관성좌표계(ECI, earth-centered inertial frame)를 사용한다 (https://pasus.tistory.com/184). 또한 기준 위성을 중심으로 하는 다른 위성의 상대 운동을 표현하기 위한 좌표계로서는 Hill 좌표계를 사용한다.

Hill 좌표계는 기준 위성의 중심에 원점이 있고 $x$ 축은 지구 중심으로부터 기준 위성의 위치 벡터 방향, $z$ 축은 기준 위성의 각운동량 벡터 방향, $y$ 축은 오른손 법칙에 의해서 정한다. Hill 좌표계는 LVLH(local-vertical-local-horizontal) 좌표계 또는 RTN(radial/along-track/cross-track) 좌표계로도 불린다.

기준 위성은 일반적으로 수동적이고 기동성이 없는 실제 위성이나 또는 가상의 위성으로 정한다. 예를 들면 어떤 우주선이 국제우주정거장과 랑데부한다고 할 때, 국제우주정거장인 타깃을 따라가는 데 필요한 기동을 수행해야 하는 추격 위성(chaser) 보다는 타깃인 국제우주정거장을 기준 위성으로 삼는다. 기준 위성을 chief 위성, 기준 위성에 대해 상대 운동으로 표현되는 위성을 deputy 위성이라고 한다.

아래 그림에서 $\{i\}$ 는 ECI 좌표계, $\{o\}$ 는 Hill 좌표계를 나타내고, $\overrightarrow{h}$ 는 chief 위성의 각운동량 벡터, ${\overrightarrow{r}}_1$ 은 지구 중심을 기준으로 한 chief 위성의 위치벡터, ${\overrightarrow{r}}_2$ 는 deputy 위성의 위치벡터, $\overrightarrow{r}={\overrightarrow{r}}_2-{\overrightarrow{r}}_1$ 은 chief 위성에 대한 deputy 위성의 상대 위치벡터이다.

 

 

Hill 좌표계의 좌표축은 각각 다음과 같이 구할 수 있다.

 

$$\begin{array}{}\text{(1)}& {\hat{o}}_1=\frac{{\overrightarrow{r}}_1}{r_1},{\hat{o}}_3=\frac{\overrightarrow{h}}{h},{\hat{o}}_2={\hat{o}}_3\times {\hat{o}}_1\end{array}$$

 

여기서 $r_1=\Vert {\overrightarrow{r}}_1{\Vert }_2$, $h=\Vert \overrightarrow{h}{\Vert }_2$ 이다.

Chief 위성과 deputy 위성은 각각 다음과 같이 이체문제의 기본 궤도 미분 방정식을 만족한다 (https://pasus.tistory.com/64).

 

$$\begin{array}{r}\text{(2)}& & \frac{{}^i{d}^2{\overrightarrow{r}}_1}{d{t}^2}+\frac{\mu }{r_1^3}{\overrightarrow{r}}_1={\overrightarrow{f}}_1\\ \text{(3)}& & \frac{{}^i{d}^2{\overrightarrow{r}}_2}{d{t}^2}+\frac{\mu }{r_2^3}{\overrightarrow{r}}_2={\overrightarrow{f}}_2\end{array}$$

 

여기서 $\mu =GM$ 은 중력 파라미터, ${\overrightarrow{f}}_1,{\overrightarrow{f}}_2$ 는 섭동력(perturbing specific force, 예를 들면 J2, 대기항력, SRP등)이다. 식 (3)에서 식 (2)를 빼면 다음과 같이 상대 궤도 미분방정식을 얻을 수 있다.

 

$$\begin{array}{}\text{(4)}& \frac{{}^i{d}^2}{d{t}^2}({\overrightarrow{r}}_2-{\overrightarrow{r}}_1)=\frac{\mu }{r_1^3}{\overrightarrow{r}}_1-\frac{\mu }{r_2^3}{\overrightarrow{r}}_2+{\overrightarrow{f}}_2-{\overrightarrow{f}}_1\end{array}$$

 

또는

 

$$\begin{array}{}\text{(5)}& \frac{{}^i{d}^2\overrightarrow{r}}{d{t}^2}=\frac{\mu }{r_1^3}({\overrightarrow{r}}_1-\frac{r_1^3}{r_2^3}({\overrightarrow{r}}_1+\overrightarrow{r}))+\overrightarrow{f}\end{array}$$

 

식 (5)로 주어지는 상대 운동 방정식은 두 위성의 거리의 크기에 제한이 없고 chief 위성의 궤도의 모양에 관계없이 유효하다.

이제 상대 위치벡터 $\overrightarrow{r}$ 의 크기가 chief 위성의 위치벡터 ${\overrightarrow{r}}_1$ 의 크기보다 매우 작은 경우를 가정해 보자. 즉

 

$$\begin{array}{}\text{(6)}& \frac{r}{r_1}\ll 1\end{array}$$

 

이 가정은 랑데부 기동이나 근접 편대 비행의 경우와 같이 두 위성이 서로 근접해 있는 경우에 해당한다. 그러면 식 (5)에서

 

$$\begin{array}{r}\text{(7)}& r_2^{-3}& =\Vert {\overrightarrow{r}}_2{\Vert }_2^{-3}={[({\overrightarrow{r}}_1+\overrightarrow{r})\cdot ({\overrightarrow{r}}_1+\overrightarrow{r})]}^{-3/2}\\ & ={[r_1^2+2\overrightarrow{r}\cdot {\overrightarrow{r}}_1+r^2]}^{-3/2}\end{array}$$

 

이다. 다음 관계식에서

 

$$\begin{array}{}\text{(8)}& (x+y{)}^r=\sum _{k=0}^{\mathrm{\infty }}(\genfrac{}{}{0ex}{}{r}{k})x^{r-k}{y}^k=x^r+r{x}^{r-1}y+\cdots \end{array}$$

 

$x=r_1^2$, $y=2\overrightarrow{r}\cdot {\overrightarrow{r}}_1+r^2$, $r=-3/2$ 로 놓고 식 (7)을 식 (8)에 대입하면 다음과 같다.

 

$$\begin{array}{r}\text{(9)}& r_2^{-3}& =r_1^{-3}-\frac{3}{2}{r}_1^{-5}(2\overrightarrow{r}\cdot {\overrightarrow{r}}_1+r^2)+\cdots \\ & =r_1^{-3}(1-\frac{3}{2}\left(\frac{2\overrightarrow{r}\cdot {\overrightarrow{r}}_1}{r_1^2}\right)-\frac{3}{2}{\left(\frac{r}{r_1}\right)}^2)+\cdots \\ \\ & \approx r_1^{-3}(1-3\left(\frac{\overrightarrow{r}\cdot {\overrightarrow{r}}_1}{r_1^2}\right))\end{array}$$

 

식 (9)를 이용하면 식 (5)에서

 

$$\begin{array}{r}\text{(10)}& & {\overrightarrow{r}}_1-\frac{r_1^3}{r_2^3}({\overrightarrow{r}}_1+\overrightarrow{r})={\overrightarrow{r}}_1-(1-3\left(\frac{\overrightarrow{r}\cdot {\overrightarrow{r}}_1}{r_1^2}\right))({\overrightarrow{r}}_1+\overrightarrow{r})\\ & ={\overrightarrow{r}}_1-({\overrightarrow{r}}_1+\overrightarrow{r}-3\left(\frac{\overrightarrow{r}\cdot {\overrightarrow{r}}_1}{r_1^2}\right)({\overrightarrow{r}}_1+\overrightarrow{r}))\\ \\ & =-\overrightarrow{r}+3\left(\frac{\overrightarrow{r}\cdot {\overrightarrow{r}}_1}{r_1^2}\right){\overrightarrow{r}}_1+\mathcal{O}(\Vert \overrightarrow{r}{\Vert }_2^2)\\ \\ & \approx -\overrightarrow{r}+3\left(\frac{\overrightarrow{r}\cdot {\overrightarrow{r}}_1}{r_1^2}\right){\overrightarrow{r}}_1\end{array}$$

 

이 되므로 식 (5)는 다음과 같이 된다.

 

$$\begin{array}{}\text{(11)}& \frac{{}^i{d}^2\overrightarrow{r}}{d{t}^2}=\frac{\mu }{r_1^3}(-\overrightarrow{r}+3\left(\frac{\overrightarrow{r}\cdot {\overrightarrow{r}}_1}{r_1^2}\right){\overrightarrow{r}}_1)+\overrightarrow{f}\end{array}$$

 

식 (11)은 식 (5)의 선형화된 버전으로서 chief 위성에 대한 deputy 위성의 상대 운동을 나타내는 방정식이다. 식 (11)은 미지수 $\overrightarrow{r}$ 에 대해서 선형식이고 ${\overrightarrow{r}}_1$ 에 관해서는 비선형 식이다. 하지만 ${\overrightarrow{r}}_1$ 은 식 (2)를 풀면 독립적으로 결정되기 때문에 미지수가 아님에 주의해야 한다.

 

 

기본 운동학 방정식(BKE)에 의하면 식 (11)에 있는 관성 좌표계에서의 미분을 Hill 좌표계에서의 미분으로 변환시킬 수 있다 (https://pasus.tistory.com/121).

 

$$\begin{array}{r}\text{(12)}& \frac{{}^{i}d\overrightarrow{r}}{dt}& =\frac{{}^{o}d\overrightarrow{r}}{dt}{+}^i{\overrightarrow{\omega }}^o\times \overrightarrow{r}\\ \frac{{}^i{d}^2\overrightarrow{r}}{d{t}^2}& =\frac{{}^{o}d}{dt}(\frac{{}^{o}d\overrightarrow{r}}{dt}{+}^i{\overrightarrow{\omega }}^o\times \overrightarrow{r}){+}^i{\overrightarrow{\omega }}^o\times (\frac{{}^{o}d\overrightarrow{r}}{dt}{+}^i{\overrightarrow{\omega }}^o\times \overrightarrow{r})\\ \\ & =\frac{{}^o{d}^2\overrightarrow{r}}{d{t}^2}+\frac{{}^o{d}^i{\overrightarrow{\omega }}^o}{dt}\times \overrightarrow{r}+2^i{\overrightarrow{\omega }}^o\times \frac{{}^{o}d\overrightarrow{r}}{dt}{+}^i{\overrightarrow{\omega }}^o\times {(}^i{\overrightarrow{\omega }}^o\times \overrightarrow{r})\end{array}$$

 

여기서 ${}^i{\overrightarrow{\omega }}^o$ 는 관성 좌표계에 대한 Hill 좌표계의 각속도 벡터이다. 식 (12)를 (11)에 대입하면 다음과 같다.

 

$$\begin{array}{r}\text{(13)}& \frac{{}^o{d}^2\overrightarrow{r}}{d{t}^2}& =\frac{\mu }{r_1^3}(-\overrightarrow{r}+3\left(\frac{\overrightarrow{r}\cdot {\overrightarrow{r}}_1}{r_1^2}\right){\overrightarrow{r}}_1)-\frac{{}^o{d}^i{\overrightarrow{\omega }}^o}{dt}\times \overrightarrow{r}\\ & -2^i{\overrightarrow{\omega }}^o\times \frac{{}^{o}d\overrightarrow{r}}{dt}{-}^i{\overrightarrow{\omega }}^o\times {(}^i{\overrightarrow{\omega }}^o\times \overrightarrow{r})+\overrightarrow{f}\end{array}$$

 

한편 chief 위성의 각운동량 벡터 $\overrightarrow{h}$ 는 다음과 같으므로

 

$$\begin{array}{r}\text{(14)}& \overrightarrow{h}& ={\overrightarrow{r}}_1\times {\overrightarrow{v}}_1=r_1{\hat{o}}_1\times (v_r{\hat{o}}_1+v_{\theta }{\hat{o}}_2)\\ & =r_1{v}_{\theta }{\hat{o}}_3=r_1^2{\omega }_{io}{\hat{o}}_3\\ \\ & =h{\hat{o}}_3\end{array}$$

 

각속도 벡터 ${}^i{\overrightarrow{\omega }}^o$ 와 그 미분은 다음과 같이 쓸 수 있다.

 

$$\begin{array}{r}\text{(15)}& {}^i{\overrightarrow{\omega }}^o& ={\omega }_{io}{\hat{o}}_3=\frac{\overrightarrow{h}}{r_1^2}\\ \frac{{}^o{d}^i{\overrightarrow{\omega }}^o}{dt}& =\frac{{}^{o}d}{dt}\left(\frac{\overrightarrow{h}}{r_1^2}\right)=-2\frac{\overrightarrow{h}}{r_1^3}{\dot{r}}_1\\ \\ & =-2\frac{\overrightarrow{h}}{r_1^3}\frac{{\overrightarrow{r}}_1\cdot {\overrightarrow{v}}_1}{r_1}=-2\frac{\overrightarrow{h}}{r_1^4}({\overrightarrow{r}}_1\cdot {\overrightarrow{v}}_1)\end{array}$$

 

여기서 ${\overrightarrow{v}}_1=\frac{{}^{i}d{\overrightarrow{r}}_1}{dt}$ 이다.

 

 

이제 식 (13)을 Hill 좌표계의 각 축 성분으로 풀어 써보자. 먼저 chief 위성의 위치벡터 ${\overrightarrow{r}}_1$ 과 상대 위치벡터 $\overrightarrow{r}$ 을 다음과 같이 Hill 좌표계의 축 성분으로 표시하면 다음과 같다.

 

$$\begin{array}{r}\text{(16)}& & {\overrightarrow{r}}_1=r_1{\hat{o}}_1\\ & \overrightarrow{r}=x{\hat{o}}_1+y{\hat{o}}_2+z{\hat{o}}_3\end{array}$$

 

그러면 식 (15)를 이용하면 식 (12)에서 각 성분은 다음과 같이 표현할 수 있다.

 

$$\begin{array}{r}\text{(17)}& \frac{{}^o{d}^i{\overrightarrow{\omega }}^o}{dt}\times \overrightarrow{r}& =-2\frac{h}{r_1^4}({\overrightarrow{r}}_1\cdot {\overrightarrow{v}}_1){\hat{o}}_3\times (x{\hat{o}}_1+y{\hat{o}}_2+z{\hat{o}}_3)\\ & =\frac{2({\overrightarrow{r}}_1\cdot {\overrightarrow{v}}_1)h}{r_1^4}(y{\hat{o}}_1-x{\hat{o}}_2)\end{array}$$

 

$$\begin{array}{r}\text{(18)}& 2^i{\overrightarrow{\omega }}^o\times \frac{{}^{o}d\overrightarrow{r}}{dt}& =2\frac{h}{r_1^2}{\hat{o}}_3\times (\dot{x}{\hat{o}}_1+\dot{y}{\hat{o}}_2+\dot{z}{\hat{o}}_3)\\ & =2\frac{h}{r_1^2}(-\dot{y}{\hat{o}}_1+\dot{x}{\hat{o}}_2)\end{array}$$

 

$$\begin{array}{r}\text{(19)}& {}^i{\overrightarrow{\omega }}^o\times {(}^i{\overrightarrow{\omega }}^o\times \overrightarrow{r})& =\frac{h}{r_1^2}{\hat{o}}_3\times (\frac{h}{r_1^2}{\hat{o}}_3\times (x{\hat{o}}_1+y{\hat{o}}_2+z{\hat{o}}_3))\\ & =\frac{h}{r_1^2}{\hat{o}}_3\times \frac{h}{r_1^2}(-y{\hat{o}}_1+x{\hat{o}}_2)\\ \\ & =-\frac{h^2}{r_1^4}(x{\hat{o}}_1+y{\hat{o}}_2)\end{array}$$

 

$$\begin{array}{}\text{(20)}& 3\left(\frac{\overrightarrow{r}\cdot {\overrightarrow{r}}_1}{r_1^2}\right){\overrightarrow{r}}_1=3\left(\frac{r_{1}x}{r_1^2}\right)r_1{\hat{o}}_1=3x{\hat{o}}_1\end{array}$$

 

식 (16)~(20)을 식 (13)에 대입하고 각 축 성분으로 나누어 쓰면 다음과 같다.

 

$$\begin{array}{r}\text{(21)}& & \ddot{x}-(\frac{2\mu }{r_1^3}+\frac{h^2}{r_1^4})x+\frac{2({\overrightarrow{r}}_1\cdot {\overrightarrow{v}}_1)h}{r_1^4}y-2\frac{h}{r_1^2}\dot{y}=\overrightarrow{f}\cdot {\hat{o}}_1=f_1\\ & \ddot{y}+(\frac{\mu }{r_1^3}-\frac{h^2}{r_1^4})y-\frac{2({\overrightarrow{r}}_1\cdot {\overrightarrow{v}}_1)h}{r_1^4}x+2\frac{h}{r_1^2}\dot{x}=\overrightarrow{f}\cdot {\hat{o}}_2=f_2\\ \\ & \ddot{z}+\frac{\mu }{r_1^3}z=\overrightarrow{f}\cdot {\hat{o}}_3=f_3\end{array}$$

 

상대 위치 좌표 $x,y,z$ 를 시간 함수로 얻기 위해서는 식 (21)을 풀어야 한다. 식에서 $x$ 와 $y$ 는 결합되어 있지만 $z$ 방향의 상대 운동은 다른 두 방향의 상대 운동과 독립적이다. chief 위성의 궤도가 타원이면 ${\overrightarrow{r}}_1$ 과 ${\overrightarrow{v}}_1$ 은 시간에 따라 변하기 때문에 식 (21)의 계수가 시간의 함수가 되므로 식 (21)은 수치적으로 풀어야 한다. 하지만 원 궤도라면 ${\overrightarrow{r}}_1\cdot {\overrightarrow{v}}_1=0$ 이므로 섭동력을 무시할 경우 해석적인 풀이가 가능하다.