# Low-Thrust Nonlinear Orbit Control for Very Low Lunar Orbits
## 2 Orbit Dynamics  
An accurate description of the orbital dynamics requires a proper definition of an inertial frame. As a preliminary step, the Moon Centered Inertial frame (MCI) is introduced, with unit vectors \((\hat{c}_{1}^{M},\hat{c}_{2}^{M},\hat{c}_{3}^{M})\) defined as follows: (i) \(\hat{c}_{1}^{M}\) lies in the Moon equatorial plane and is coplanar with the line that connects the Earth and the Moon at a reference time \(t_{ref}\), directed toward the far side of the Moon; (ii) \(\hat{c}_{3}^{M}\) is aligned with the Moon rotation axis at \(t_{ref}\); (iii) \(\hat{c}_{2}^{M}\) completes the right-handed set. The Local Vertical Local Horizontal (LVLH) frame is instead aligned with unit vectors \((\hat{r},\hat{\theta},\hat{h})\): (i) \(\hat{r}\) is directed from the Moon center to the instantaneous position of the vehicle; (ii) \(\hat{h}\) is aligned with the spacecraft angular momentum; (iii) \(\hat{\theta}\) completes the right-handed set. If \(\mathbf{R}_{i}(\alpha)\) denotes an elementary rotation matrix by angle \(\alpha\) about axis \(i\), then the two reference frames can be linked through three successive rotations,  
\[ \begin{bmatrix}\hat{r}\\ \hat{\theta}\\ \hat{h}\end{bmatrix} =\mathbf{R}_{3}(\omega+\theta)\mathbf{R}_{1}(i)\mathbf{R} _{3}(\Omega) \begin{bmatrix}\hat{c}_{1}^{M}\\ \hat{c}_{2}^{M}\\ \hat{c}_{3}^{M}\end{bmatrix} , \tag{1}\]  
where \(\Omega\) is the Right Ascension of the Ascending Node (RAAN), \(i\) the inclination of the orbital plane, \(\omega\) the argument of periselenium, and \(\theta\) the true anomaly. Furthermore, we introduce the Local Horizontal (LH) frame, associated with unit vectors \((\hat{r},\hat{E},\hat{N})\), where \(\hat{E}\) and \(\hat{N}\) are aligned with the local (lunar) East and North directions, respectively. If \(\zeta\) is the heading angle, then the following relation holds,  
\[ \begin{bmatrix}\hat{r}\\ \hat{\theta}\\ \hat{h}\end{bmatrix} =\mathbf{R}_{1}(\zeta) \begin{bmatrix}\hat{r}\\ \hat{E}\\ \hat{N}\end{bmatrix} . \tag{2}\]  
Figure 1 shows the reference frames and related angles. The spacecraft is modeled as a point mass and its motion is mainly affected by the Moon gravitational field. The presence of perturbing effects is also addressed: (i) the gravitational pull of Earth and Sun; (ii) the relevant harmonics of the selenopotential (down to two orders of magnitude less than the largest harmonic term \(J_{2}\)), that is, those associated with coefficients \(|J_{lm}|>10^{-6}\). Appendix A collects all planetary data employed in this research. The gravitational influence of the Earth and Sun can be modeled as third body perturbations following Battin-Giorgi approach ([25; 26]):  
\[\bar{a}_{3i}=-\frac{\mu_{i}}{||\bar{r}_{i}||^{3}(q_{i}+1)^{3/2}}\left(\hat{r }+\hat{r}_{i}q_{i}\frac{q_{i}^{2}+3q_{i}+3}{(q_{i}+1)^{3/2}+1}\right),\ \ q_{i}:=\frac{\hat{r} \cdot\hat{r}-2\hat{r}\cdot\hat{r}_{i}}{\hat{r}_{i}\cdot\hat{r}_{i}}, \tag{3}\]where \(\hat{\mathbf{r}}_{i}\) denotes the vector position of the \(i\)-body with respect to the Moon and \(\hat{\mathbf{r}}\) the instantaneous position vector of the spacecraft with respect to the central body. Instead, the selenopotential \(U\) can be expanded in terms of Legendre polynomials \(P_{l,m}\),  
\[U=\frac{\mu_{m}}{r}-\frac{\mu_{m}}{r}\sum_{l=2}^{\infty}\Big{(}\frac{R_{m}}{l} \Big{)}^{l}J_{l}P_{l0}(sin\phi)+\frac{\mu_{m}}{r}\sum_{l=2}^{\infty}\sum_{m=1} ^{l}\Big{(}\frac{R_{m}}{r}\Big{)}^{l}J_{l,m}P_{l,m}(\sin\phi)\cos\big{[}m( \lambda_{g}-\lambda_{l,m})\big{]}, \tag{4}\]  
where \(R_{m}\) is the lunar equatorial radius, \(\phi\) the geographical latitude and \(\lambda_{g}\) the geographical longitude, and the coefficients (\(J_{lm}\) and \(\lambda_{lm}\)) of zonal, tesseral, and sectoral harmonics are provided by the Lunar Prospector LP100K model for the Moon gravitational field [27].