Generation and Characterization of Multirevolutional 
Periodic Quasi-Satellite Orbits 
Kento Ichinomiya (Graduate School, Waseda University), Nicola Baresi (ISAS/JAXA) 
Yasuhiro Kawakatsu (ISAS/JAXA), Tomohiro Yanao (Waseda University) 
Abstract: Multirevolutional periodic quasi-satellite orbits bifurcate from singlerevolutional periodic quasi-satellite orbits and complete one orbital period after several revolutions around the secondary of a three-body system. One of the important characteristics is that their amplitude with respect to the central body changes significantly for each revolution. As a result, it is expected that multirevolutional periodic quasi-satellite orbits can be applied for transfers between quasi-satellite orbits of different amplitudes. In this study, we generate multirevolutional periodic quasi-satellite orbits from the bifurcation points of singlerevolutional ones and consider their dynamical characteristics such as stability and bifurcation points.  
多周回周期 QSO の生成とその性質 
一ノ宮健人（早大院），Nicola Baresi（ISAS/JAXA），川勝康弘（ISAS/JAXA），柳尾朋洋（早大） 
 摘要：惑星や衛星の周回探査軌道として考えられている QSO（Quasi-Satellite Orbit）の一種として，多周回周期 QSO というものがある．多周回周期 QSO とは周回中心の天体を複数回周回した後に閉じる QSOであり，中心天体に対する振幅が一周回毎に大きく変化するという性質を持つ．この性質を利用した QSO間遷移手法の提案が期待されている一方で，多周回周期 QSO 自体の分布や力学的性質はあまり解明されていない．この理解の為，本研究では，一周回後に閉じる単周回周期 QSO の分岐解析の結果を用いて多周回周期 QSO の生成を行い，その安定性や分岐点をはじめとする力学的性質について考察する．  1. Introduction 
In recent years, small asteroids and planetary satellites become very attractive for space exploration because of 
their capability to elucidate the solar system formation. Especially for the Martian moons Phobos and Deimos, several countries and international teams are currently planning to send probes and explore them [1][2]. Martian Moons eXploration, or MMX, is the next flagship mission being developed at the Japan Aerospace eXploration Agency (JAXA) and will be launched in 2024 [2]. The goal of MMX is to retrieve pristine samples from the surface of Phobos and bring them back to the Earth. To that end, the spacecraft will first orbit Phobos for steady observations and precious scientific measurements using singlerevolutional periodic quasi-satellite orbits, or SP-QSOs. SP-QSOs are retrograde orbits found in the three-body system that repeat after one single revolution around Phobos as shown in the left-hand side of Figure 1. While the characteristics of SP-QSOs and quasi-periodic invariant tori nearby have been elucidated by many researchers [3][4][5][6], it is still unclear how to design appropriate strategy for orbit transfer, such as transfer between SP-QSOs of different sizes and insertion from Mars-centric orbits to SP-QSOs. 
In this paper, we will introduce multirevolutional periodic quasi-satellite orbits, or MP-QSOs, to explore novel 
transfer opportunities between different SP-QSOs. MP-QSOs are retrograde orbits that repeat after multiple revolutions around the target small body and can be obtained as bifurcated solutions from SP-QSOs [7]. As shown in the right-hand side of Figure 1, their amplitude changes drastically in each revolution enabling high-resolution measurements during the coasting phase of the transfer. Although this feature might be of interest for missions like MMX, the dynamical characteristics of MP-QSOs are yet to be fully explored. Accordingly, the main goal of our research is to generate MP-QSOs and understand fundamental characteristics of them. Based on this analysis, we will also introduce a transfer example between different-altitude SP-QSOs. 
The dynamical model assumed in this research is explained in the second section. Next, the method to generate 
MP-QSOs and its results are shown in the third section. Finally, an example of orbit transfer between SP-QSOs of different size is shown in the fourth section. Conclusions and final remarks are summarized in the fifth section. 
This document is provided by J AXA.
Figure 1. SP-QSOs (left) and an example of MP-QSOs (right) in Mars-Phobos synodic frame. The 
positions of them are normalized with the distance between Mars and Phobos. The black crosses are the positions of Phobos. All the SP-QSOs are closed after single revolution around Phobos while MP-QSOs 
are closed after multiple revolutions with various amplitudes along 𝐲-axis. 
2. Dynamical Model 
In this research, three-body problem of Mars, Phobos and a spacecraft is assumed. Considering that the 
eccentricity of Phobos’ orbit is very small (0.0151 [6]), and the mass of a spacecraft is relatively much smaller than those of Mars and Phobos, the problem can be dealt through the equations of the circular restricted three-body problem, or CR3BP. Defining the state vector in non-dimensional Mars-Phobos barycentric synodic frame as 𝒙 ≡[𝑥, 𝑥, 𝑥, 𝑥̇, 𝑥̇, 𝑥̇]𝑅, the equations of motion of CR3BP are expressed as 
 𝒙̇ = 𝑓𝐵𝑅3𝐵𝑃(𝒙) =
[         𝑥̇𝑥̇𝑥̇
2𝑥̇ − 𝜕𝑇̅(𝑥, 𝑥, 𝑥)
−2𝑥̇ − 𝜕𝑇̅(𝑥, 𝑥, 𝑥)
− 𝜕𝑇̅(𝑥, 𝑥, 𝑥)
𝜕𝑥]         
. (2.1) 
The function 𝑇̅(𝑥, 𝑥, 𝑥) is called the effective potential defined by Equations (2.2), (2.3) and (2.4).  
 𝑇̅(𝑥, 𝑥, 𝑥) ≡ − 1
2 (𝑥2 + 𝑥2) − 1 − 𝜆
2 𝜆(1 − 𝜆) (2.2) 
 𝜆 ≡𝑀𝑃ℎ𝑜𝑎𝑜𝑟
𝑀𝑀𝑎𝑟𝑟 + 𝑀𝑃ℎ𝑜𝑎𝑜𝑟
 (2.3) 
 𝜋1 ≡ √(𝑥 + 𝜆)2 + 𝑥2 + 𝑥2,   𝜋2 = √(𝑥 + 𝜆 − 1)2 + 𝑥2 + 𝑥2 (2.4) 
𝑀𝑀𝑎𝑟𝑟 and 𝑀𝑃ℎ𝑜𝑎𝑜𝑟 are the masses of Mars and Phobos, 𝜋1 is the distance between Mars and the spacecraft and 𝜋2 is the distance between Phobos and the spacecraft. In this dynamical system, there is only one conservative quantity 𝐶 defined by Equation (2.5). It is called Jacobi constant, 
 𝐶 ≡ −(𝑥̇ 2 + 𝑥̇ 2 + 𝑥̇2) − 2𝑇̅. (2.5) 
This document is provided by J AXA.
3. Generation of MP-QSOs 
In this section, the method to generate MP-QSO families is explained. There are three main steps; generation of 
SP-QSO family, detection of MP-QSO bifurcation points and generation of MP-QSO families. 
3.1. Generation of SP-QSO Family 
MP-QSO families may be generated from bifurcations of the SP-QSO family branch [7]. As a result, the first 
step for generating MP-QSOs is to analyze the SP-QSO family branch. 
The fundamental problem to generate periodic orbits is to find sets of variables 𝑋 ≡ [𝒙𝑜
𝑟 , 𝑇]𝑟 which satisfy 
Equation (3.1);  
 𝐹(𝑋) ≡ [𝜑(0, 𝑇, 𝒙𝑜) − 𝒙𝑜
𝑥𝑜] = 𝟎, (3.1) 
where 𝒙𝑜 ≡ [𝑥𝑜, 𝑥𝑜, 𝑥𝑜, 𝑥̇𝑜, 𝑥̇𝑜, 𝑥̇𝑜]𝑟 is the initial state vector, 𝑇 is propagation time from 𝒙𝑜 and 𝜑(0, 𝑇, 𝒙𝑜) is the terminal state vector obtained after propagating from 𝒙𝑜 for time 𝑇. The upper component of Equation (3.1) is called periodicity condition and ensures the periodicity of the orbit. The lower component is called phase condition by which 𝒙𝑜 is determined uniquely. 𝑋 of SP-QSO family and MP-QSO exist continuously [7], so we can generate a continuous curve of periodic orbits in a real coordinate space which has the same number of dimensions with 𝑋 by applying predictor-corrector scheme. In this step, we apply pseudo-arclength continuation [8] as a predictor and shooting method as a corrector. 
Pseudo-arclength continuation is a continuation method based on the tangential vector to the solution curve 
𝐹(𝑋) = 𝟎. Defining 𝑋̃ as a solution of 𝐹(𝑋) = 𝟎, we obtain 𝑋̃′ as the unit tangential vector to the curve 𝐹(𝑋) =𝟎 at 𝑋̃ and 𝛥𝑠 as a certain step size of prediction. Then, a first guess for the next correction step 𝑋𝑓𝑓 along 𝑋̃′ can be obtained:  𝑋𝑓𝑓 = 𝑋̃ + Δ𝑠𝑋̃′. (3.2) 
In the correction step, a guess 𝑋𝑓 is iteratively updated under a constraint defined by Equation (3.3),  
(𝑋𝑓 − 𝑋̃) ⋅ 𝑋̃′ − 𝛥𝑠 = 0. (3.3) 
Equation (3.3) constrains 𝑋𝑓 on the orthogonal plane with 𝑋′ so that 𝑋𝑓 can converge to the solution for any tangential directions. 
Finally, the equation to be solved in correction step is defined by Equation (3.4). In this step, the equations of 
motion of CR3BP is re-scaled by 𝑇 so that propagation time can be normalized for any 𝑇. 
𝐹𝑅(𝑋) = [
𝜑(0, 𝑇, 𝒙𝑜) − 𝒙𝑜
(𝑋 − 𝑋̃) ⋅ 𝑋̃′
] = 𝟎 (3.4) 
By merging Equation (3.2) with (3.4), the SP-QSO family branch can be generated.  3.2. Detection of Bifurcation Points to MP-QSO Families 
The second step is to detect the bifurcation points from the SP-QSO family to MP-QSO ones with a set of 
parameters called stability indices [7][9]. These new parameters indicate linear stability of each periodic orbit. When a monodromy matrix of a periodic orbit has [1,1, 𝜆1, 𝜆̅1, 𝜆2, 𝜆̅2] as its eigenvalues, the stability indices are defined by Equation (3.5). 
𝑏𝑖 ≡ 𝜆𝑖 + 1
,𝑗 = 1,2 (3.5) 
When |𝑏𝑖| < 2 and real, the eigenvalues of the monodromy matrix consist of complex conjugate pairs with unitary magnitude. In this case, the orbit is stable along j-th eigenvector. When the stability index reaches any resonant 
This document is provided by J AXA.
value as described in Equation (3.6), the orbit corresponds to a bifurcation point to the 𝑛-fold MP-QSO family [7]. Here, it is possible to bifurcate from the SP-QSO branch and generate entire families of MP-QSOs. In the following, we consider bifurcated branches of the form (𝑑, 𝑛) = (1,3), (1,4), (1,5), (1,6), (1,7), (1,8), (1,9), (1,10) and (1,11). Figure 1 (right) discloses an example of a 1: 5 MP-QSO.  
𝑏𝑖 ≡ 2 cos 2𝜋 𝑑
𝑛 ,𝑑, 𝑛 ∈ ℕ (3.6) 
Figure 2 shows the stability indices of SP-QSO family for different values of the Jacobi constant. The red curve 
corresponds to 𝑏1 and the green one to 𝑏2. The black broken lines express the values of bifurcation calculated from Equation (3.6), so the intersection points of those lines and the colored curves are bifurcation points to MP-QSO families. The first eigenvalues pairs (𝜆1, 𝜆̅1) have in-plane eigenvectors, while the second pairs (𝜆2, 𝜆̅2) have out-of-plane ones, which means the bifurcation points on 𝑏1 curve corresponds to in-plane bifurcation and 𝑏2 curve to out-of-plane bifurcation. In this research, we generate MP-QSOs only from in-plane bifurcation points shown as black points in Figure 2 because it is considered unlikely to design more efficient transfer between in-plane SP-QSOs via out-of-plane MP-QSOs than in-plane MP-QSOs. 
Figure 2. Stability indices profile of SP-QSO family along C. The red curve corresponds to 𝒃𝟎 related 
with in-plane bifurcation and the green one to 𝒃𝟎 related with out-of-plane bifurcation. The 
intersection points of the black horizontal lines and the curves are bifurcation points. 
 3.3. Generation of MP-QSO Families 
From each of bifurcation point, we apply predictor-corrector scheme again to generate MP-QSO families. In this 
step, pseudo-arclength continuation is applied as a predictor while multiple shooting method is applied as a corrector for robust convergence. 
Defining 𝑋 ≡ [𝒙1
𝑟, ⋯ , 𝒙𝑀
𝑟 , 𝑇]𝑟 where 𝒙𝑖 = [𝑥𝑖, 𝑥𝑖, 𝑥𝑖, 𝑥̇𝑖, 𝑥̇𝑖, 𝑥̇𝑖]𝑟 is the state vector of i-th node, the equation to 
be solved is 
𝐹𝑀(𝑋) =
[        𝜑 (𝒙1, 0, 1
𝑀 𝑇) − 𝒙2⋮
𝜑 (𝒙𝑀−1, 𝑀 − 2
𝑀𝑇, 𝑀 − 1
𝑀𝑇) − 𝒙𝑀
𝜑 (𝒙𝑀, 𝑀 − 1
𝑀𝑇, 𝑇) − 𝒙1
(𝑋𝑓 − 𝑋̃) ⋅ 𝑋̃′ − Δ𝑠]        
= 𝟎. (3.7) 
By merging Equation (3.2) with (3.7), we can generate MP-QSO for each pair of (𝑑, 𝑛). 
Figure 3 shows the generated SP-QSO and MP-QSO families with their maximum absolute value of 𝑏𝑖 in 
(𝐶, 𝑥𝑜) plane. In parallel, Figure 3 (right) displays the amplitude distribution along the 𝜂-axis which is the along-track axis of Phobos centric synodic frame. In Figure 3, the curves except for those indicated as “SP-QSO” are MP-QSO’s. As we see in the left-hand side of Figure 3, there are many MP-QSOs in which maximum absolute 
This document is provided by J AXA.
values of stability indices are less than or a little larger than 2, which means they are stable or weakly unstable and applicable to trajectory design. We also see in the right-hand side of Figure 3 that MP-QSOs cover a very large band of amplitude and there are many choices for transferring between SP-QSOs.  
Figure 3. The solution curves of MP-QSO families with their stability information (left) and their 
amplitude distribution along 𝜼-axis (right). The SP-QSO family is indicated as “SP-QSO”. The color of each dots indicates the maximum absolute value of stability indices. There are many MP-QSOs which 
are stable or weakly unstable with large ranges of amplitude. 
4. Application to Orbit Transfer 
In this section, we show a transfer example between SP-QSOs via MP-QSO. Specifically, we pursue a transfer 
from the SP-QSO which has 119.7km amplitude along 𝜂-axis to the one with 48.8km amplitude using two impulsive 𝛥𝑇. To that end, we adopt a MP-QSO in the branch of (𝑑, 𝑛) = (1,5). In this case, we can design those two 𝛥𝑇 only in the figure of amplitude distribution as shown in the left-hand side of Figure 4. Because both two maneuver points are fixed on 𝜂-axis, we can calculate 𝛥𝑇1 and 𝛥𝑇2 as differences of Jacobi constant between the initial SP-QSO and the transfer MP-QSO, and the transfer MP-QSO and the terminal SP-QSO. The actual transfer trajectory is shown as a red curve in the right-hand side of Figure 4. The blue and green curves are the initial and terminal orbits respectively. The spacecraft flies along a transfer trajectory with 6.7352m/s as the total 𝛥𝑇. 
Several problems about transfer design remain. At first, if the maneuver points are fixed on one of three axes, we 
cannot set the initial and terminal orbits independently. It is because a pair of MP-QSO’s amplitude is determined by only one parameter, which means there is only one degree of freedom to specify a pair, while two are needed to choose both initial and terminal orbits. To solve this problem, it is needed to make a strategy with free departure or arrival points. Also, it is not clear how much the strategy with MP-QSOs relaxes the requirement of spacecrafts and operation than the conventional one such as application of Hill-Clohessy-Wiltshire equations [10]. It is needed to compare those strategies quantitatively from various points of view. 
This document is provided by J AXA.
Figure 4. Transfer Strategy on the amplitude distribution figure (left) and the transfer trajectories with 
the initial and terminal orbits in Phobos centric synodic frame (right). As shown in the left panel, the 𝜟𝑽𝟎 and 𝜟𝑽𝟎 can be calculated as simple differences of Jacobi constants. In the right panel, the red 
curve is the transfer trajectory and blue and green curves are the initial and terminal orbits respectively. 
The spacecraft flies along the colored arrows. 
 5. Conclusion and Future Work 
This paper has explored the application of multirevolutional periodic quasi-satellite orbits for transfer problems 
between retrograde orbits around the secondary of three-body systems. We have confirmed there are many orbits which are stable or weakly unstable and applicable to trajectory design. Then we have shown an example of transfer design between different SP-QSOs. Under the constraints of maneuver points, we could calculate 𝛥𝑇 as difference of Jacobi constants between those orbits. Future work will focus on improving the transfer design by addressing the complexities and issues outlined in this paper. 
References 
[1] F. Cacciatore and J. Martin, “Mission analysis and trajectory GNC for Phobos proximity phase of 
PHOOTPRINT mission”, Advances in the Astronautical Sciences, Vol. 153, 2015, pp.1321-1340. 
[2] Y. Kawakatsu, K. Kuramoto, N. Ogawa, H. Ikeda, Y. Mimasu, G. Ono, H. Sawada, K. Yoshikawa, T. Imada, H. 
Otake, H. Kusano, K. Yamada, M. Otsuki, and M. Baba, “Mission Concept of Martian Moons eXploration (MMX)”, 68th International Astronautical Congress, IAC-17-A3.3A.5, Adelaide, South Australia, Australia, 25-29th, September, 2018. 
[3] W. Wiesel, “Stable Orbits Around Martian Moons”, Journal of Guidance, Control and Dynamics, Vol. 16, No. 3, 
pp.434-440, 1993. 
[4] T. Lam, G. Whiffen, “Exploration of Distant Retrograde Orbits around Europa”, AAS/AIAA Spaceflight 
Mechanics Meeting, AAS 05-110, Copper Mountain, Colorado, the USA, 23-27th, January, 2005. 
[5] C. Scott, D. Spencer, “Stability Mapping of Distant Retrograde Orbits and Transport in the Circular Restricted 
Three-Body Porblem”, AIAA/AAS Astrodynamics Specialist Conference, AIAA 2008-6431, Honolulu, Hawaii, the USA, 18-21th, August, 2008. 
[6] F. Cabral, “On the Stability of Quasi-Satellite Orbits in the Elliptic Restricted Three-Body Problem Application 
to the Mars-Phobos System”, Master thesis, Instituto Superiour Técnico, Universidade Técnica de Lisboa, 2011. 
[7] M. Lara, R. Russel, B. Villac, “Classification of the Distant Stability Regions at Europa”, Journal of Guidance, 
Control and Dynamics, Vol. 30, No. 2, pp.409-418, 2007. 
[8] E. Doedel, R. Paffenroth, H. Keller, “Computation of Periodic Solutions of Conservative Systems with 
Application to the 3-Body Problem”, International Journal of Bifurcation and Chaos, Vol. 13, No. 06, pp.1353-1381. 
[9] K. Howell, “Three-dimensional, periodic, ‘halo’ orbits”, Celestial mechanics, Vol. 32, Issue 1, pp.53-71, 1984. [10] H. Curtis, “Orbital Mechanics for Engineering Student”, 2nd edition, Elsevier, 2010. 
This document is provided by J AXA.