Orbital Mechanics - Ch.3 Orbital Position as a Function of Time

 

2장에서는 이체 문제(Two-body problem)에 대한 위치와 전근점 이각(True anomaly)간의 관계를 다루었다. 시간에 대한 정보가 나타난 유일한 부분은 타원의 주기 부분이었다. 원 궤도를 가지는 경우, 시간에 대한 함수로 위치를 얻는 것은 간단하다. 타원, 포물선, 쌍곡선 경로에 대해서, 시간에 대한 상대적인 위치로 케플러 방정식의 다양한 형태를 유도할 것이다. 해당 초월 방정식은 뉴턴 기법과 같은 과정을 이용하여 반복적으로 풀어야한다. 

 케플러 방정식의 다른 형태는 single universal Kepler's equation과 결합한다.

 

3.2 Time Since Periapsis

궤도 공식은 $r=\left( h^2 /\mu \right) / \left( 1+e\text{cos}\theta \right)$ 강체 $m_1$ 근방의 $m_2$의 전근점 이각(True Anomaly)의 함수로써의 위치를 제공한다. 실질적 이유로 인해 $m_2$의 위치를 시간에 대한 함수로 결정할 수 있을 필요가 있다.

 

$$\frac{d\theta}{dt} = \frac{h}{r^2}$$

$$\frac{\mu^2}{h^3} dt = \frac{d\theta}{\left(1+e \text{cos}\theta \right)^2 }$$

$$\frac{\mu^2}{h^3} \left( t -t_p \right) = \int_{0}^{\theta} \frac{d\theta}{\left(1+e \text{cos}\theta \right)^2 }$$

 

3.3 Circular Orbits 원형 궤도

원에 대해서, $e=0$, 적분은 $\int d\theta$으로 바뀐다.

$$t=\frac{h^3}{\mu^2}\theta$$

$$t=\frac{r^{3/2}}{\sqrt{\mu}} \theta$$

원 궤도의 주기 $T=2\pi r^{3/2} /\sqrt{\mu}$에 대해서

$$t=\frac{\theta}{2\pi} T \hspace{10mm} \theta=\frac{2\pi}{T} t$$

 

원은 어떤 직경이든 대칭이므로, 장축 선 apse line(그러므로 장축단)은 임의로 선택 가능하다.

3.4 Elliptical Orbits 타원 궤도

$$\int_{0}^{\theta} \frac{d\theta}{\left(1+e\text{cos}\theta \right)^2 } = \frac{1}{\left(1-e^2 \right)^{3/2}}   \left[ 2 \text{tan}^{-1} \left( \sqrt{\frac{1-e}{1+e} \text{tan}\frac{\theta}{2}} \right) - \frac{e \sqrt{1-e^2} \text{sin} \theta}{1+e \text{cos} \theta} \right]$$

$$\frac{\mu^2}{h^3} t= \frac{1}{\left(1-e^2 \right)^{3/2}}   \left[ 2 \text{tan}^{-1} \left( \sqrt{\frac{1-e}{1+e} \text{tan}\frac{\theta}{2}} \right) - \frac{e \sqrt{1-e^2} \text{sin} \theta}{1+e \text{cos} \theta} \right]$$

$$M_e = 2 \text{tan}^{-1} \left( \sqrt{\frac{1-e}{1+e} \text{tan}\frac{\theta}{2}} \right) - \frac{e \sqrt{1-e^2} \text{sin} \theta}{1+e \text{cos} \theta}$$

$M_e$는 평균 근점 (Mean anomaly)라 부른다. 이심률 $e$에 대해 관찰해보면, $M_e$는 전근점 이각 true anomaly에 대해 단조 증가 함수(기울기가 항상 0보다 큰)이다.(3.5)

$$M_e = \frac{\mu^2}{h^3} \left( 1-e^2 \right)^{3/2} t$$

타원 궤도의 위치 벡터의 각속도 벡터는 일정하지 않다. 그러나 주기 $T$ 마다 $2\pi$ 라디안만큼 지나기 때문에, 비율 $2\pi/T$는 평균 각속도가 될 수 있으며 이를 $n$이라 표기하고 평균 운동(Mean motion)이라 부르자. (3.6)

$$n=\frac{2\pi}{T}$$

평균 운동의 관점에서 식 (3.5)는 다음처럼 간단해질 수 있다.

$$M_e = nt$$

평균 근점(Mean anomaly)는 일정한 각속력 $n$으로 타원을 도는 움직이는 가상의 강체의 방위각(azimuth)이다. 원 궤도에 대해서는 평균 근점 mean anomaly $M_e$과 전근점 이각 true anomaly $\theta$는 동일하다.

식 (3.3)을 이심 근점 eccentric anomaly라 부르는 축 각 auxiliary angle $E$를 이용해서 간단하게 만들면 편리하다. 이는 타원의 장반경 $a$와 동일한 반경을 가지는 동심 보조원이 타원에 외접하도록 하면 된다.

$E$를 $theta$에 대한 함수로 표현하기 위해 기하 관계로부터 유도하면,

$$a \text{cos}E = ae + r \text{cos} \theta$$

반경 $r$의 식 (2.62)으로부터, (3.7a/b)

$a \text{cos}E = ae + \frac{a \left( 1-e^2 \text{cos} \theta \right) }{1 + e \text{cos}\theta }$$

$$\text[cos}E = \frac{e+\text{cos}\theta}{1+ e \text{cos}\theta}$$

$$\text{cos}\theta = \frac{e-\text{cos}E}{e \text{cos}E - 1}$$

$$\text{cos}\E = \frac{\sqrt{1-e^2} \text{sin}\theta}{1 + e\text{cos}\theta}$$

 

$$E = 2 \text{tan}^-1 \left( \sqrt{\frac{1-e}{1+e} \text{tan} \frac{\theta}{2}} \right)$$

위 식으로부터 식 (3.3)은 다음의 케플러 공식을 유도한다. (3.11)

$$M_e = E-e \text{sin} E$$

어떤 진근점 원각 true anomaly $\theta$ 이심 근점 $E$ 을 식 (3.10)으로부터 계산할 수 있다.

케플러의 공식인 식 (3.11)에 $E$를 대입하면 바로 평균 근점 mean anomaly를 얻을 수 있다. 이로부터, 주기 $T$는 식 (3.5)으로부터 얻을 수 있다. (3.12)

$$t = \frac{M_e}{2\pi}T$$

달리 표현하면 어떤 시간이 주어졌을 때, 식 (3.12)로부터 평균 근점 mean anomaly $M_e$를 얻을 수 있다는 것이다. $M_e$를 케플러 공식에 대입함으로써 다음과 같은 표현의 이심 근점 eccentric anomaly 을 얻을 수 있다.

$$E- e\text{sin}E = M_e$$

$E$에 대한 초월함수를 바로 풀 순 없으나 Newton's method와 같은 반복법으로 근사 해는 얻을 수 있다.

 

달을 대상으로 예시를 들어보면.. (km)

$\mu_E = G\left( m_{Earth}) = 0.3986 \cdot \text{E+6km}$

Moon's Perigee $r_p = 0.3633 \cdot \text{E+6km}$

Moon's Apogee $r_a = 0.4055\cdot \text{E+6km}$

$e=\frac{r_a - r_p}{r_a + r_p} =0.054891$

$$r_p = \frac{h^2}{\mu} \frac{1}{1+e\text{cos}0} \rightarrow h=0.039085 \cdot \text{E+6km}^2/s$$

$$T=\frac{2 \pi}{\mu^2}\left(\frac{h}{\sqrt{1-e^2}}\right)^3 = 2.3718 \cdot \text{E+6sec}$$

 

 

Q. 달이 perigee에서 전근점 이각 true anomaly 120도가 되는데 걸리는 시간은..?

 

% Moon
rp = 0.3633E6;
ra = 0.4055E6;
mu = 0.3986E6;
 
sec2day = 1/(24*60*60);
R2D = 180/pi;
D2R = pi/180;
 
% True anomaly
theta = 120;
 
e = (ra-rp)/(ra+rp);
h = sqrt(rp*mu*(1+e*cos(0)));
T = (2*pi/(mu*mu))*(h/sqrt(1-e^2))^3;
fprintf("e = %8.6f \n",e);
fprintf("h = %10.1f km2/s\n",h);
fprintf("T = %10.1f sec, %10.6f days\n",T, T*sec2day);
E = 2*atan(sqrt((1-e)/(1+e))*tand(theta/2))*R2D;
Me = (E - e * sind(E)*R2D);
fprintf("theta %7.3f deg / %7.5f rad, E %7.3f deg / %7.5f rad\n",theta, theta*D2R, E, E*D2R);
fprintf("Me %7.3f deg / %7.6f rad\n",Me ,Me*D2R);
t = Me*D2R*T/(2*pi);
fprintf("Time since perigee %10.3f sec, %6.4f hr, %6.4f day\n",t, t/60/60, t*sec2day);

 

결과는...

e = 0.054891 
h =   390845.5 km2/s
T =  2371844.9 sec,  27.451909 days
theta 120.000 deg / 2.09440 rad, E 117.237 deg / 2.04617 rad
Me 114.441 deg / 1.997365 rad
Time since perigee 753986.936 sec, 209.4408 hr, 8.7267 day