# Collider search of light dark matter model with dark sector decay

Yu Cheng<EMAIL_ADDRESS>School of Physics, East China University of
Science and Technology, Shanghai 200237, China Tsung-Dao Lee Institute, and
School of Physics and Astronomy, Shanghai Jiao-Tong University, Shanghai,
China  Wei Liao<EMAIL_ADDRESS>School of Physics, East China University of
Science and Technology, Shanghai 200237, China Qi-Shu Yan  School of Physical
Sciences, University of Chinese Academy of Sciences, Beijing 100049, China
Center for Future High Energy Physics, Institute of High Energy Physics,
Chinese Academy of Sciences, Beijing 100049, China

###### Abstract

We explore the possibility that the dark matter relic density is not produced
by thermal mechanism directly, but by the decay of other heavier dark sector
particles which on the other hand can be produced by the thermal freeze-out
mechanism. Using a concrete model with a light dark matter from dark sector
decay, we study the collider signature of the dark sector particles in
association with Higgs production processes. We find that the future lepton
colliders can be a better place to probe the signature of this kind of light
dark matter model than the hadron collider such as LHC. Meanwhile, it is found
that a Higgs factory with center of mass energy 250 GeV has a better potential
to resolve the signature of this kind of light dark matter model than the
Higgs factory with center of mass energy 350 GeV.

##  1 Introduction

The present direct and indirect search of Dark Matter(DM) have so far given a
lot of constraints on the DM interaction [1, 2, 3, 4, 5, 6, 7, 8] . In
particular, the direct detection experiments [1, 2, 3] have pushed the DM-
nucleus interaction cross section to be less than around $10^{-46}$ cm2 for
GeV-TeV scale DM particles. These direct detection measurements have raised
serious concerns that the GeV-TeV scale DM should have interactions much
weaker than we thought, such that for the DM candidate of weakly interacting
massive particle, an over-production of the relic density in thermal
production mechanism seems hard to avoid.

One way out of this difficulty is to note that the dark sector particle
directly involved in the thermal production mechanism may not be the DM
candidate, but rather the mother particle of a decay process with the DM
particle as the children particle. In such a scenario, the DM relic density is
derived from the relic density of another dark sector particle obtained from
the thermal production mechanism. Since the mass of the DM particle could be
much smaller than the mass of other dark sector particle, one can obtain a
relic density of DM much smaller than the relic density directly obtained in
the thermal production mechanism. Therefore, in such a scenario, one can still
have a weak scale interaction between the dark sector and the visible sector,
and meanwhile obtain a DM relic density consistent with cosmological
observations.

This possibility of light DM from dark sector decay has been explored in a
previous work [9] . In this article, we further study some collider
phenomenology of this scenario of DM. We will focus on the potential collider
signature of a DM model which can have light DM from dark sector decay. We
study the feasible signatures at future lepton colliders, such as CEPC, FCC-ee
and ILC. For the sake of comparison, we also study the signature and
constraints from the hadron collider such as the LHC.

The article is organized as follows. In the next section, we quickly review
the DM model we are considering and the various constraints on the model. Then
we study the signature of the DM model at future lepton colliders. After a
detailed analysis of signature on the lepton collider, we give a brief
analysis on the signature on the hadron collider. We summarize in the last
section.

##  2 The model of light DM and dark particle search at colliders

The model we consider is an extension of the Type-I seesaw model where the
right-handed neutrinos $N_{R}$ are coupled to a dark sector consisting of a
real scalar $\phi$ and a Dirac fermion $\chi$ [9]. The model is similar to DM
models in which the $N_{R}$ works as a portal to the DM[10, 11, 13, 14, 12].
Previous works on this type of models have focused on the possibility that the
DM relic density is achieved by the thermal freeze-out mechanism[13, 14, 15]
or by the freeze-in mechanism[16, 19, 17, 18]. Here, we consider the
possibility that the DM relic density is given mainly by the decay of other
heavier particles in the dark sector.

We impose a $Z_{2}$ symmetry in the Lagrangian with the dark sector particles
$\chi$ and $\phi$ odd and SM fermions even under the $Z_{2}$ operation. This
guarantees the lightest of the dark sector particles to be stable, and makes
it a DM candidate. We can write the full Lagrangian as the sum of four parts

$\mathcal{L}=\mathcal{L}_{\mathrm{SM}}+\mathcal{L}_{\text{Seesaw
}}+\mathcal{L}_{\mathrm{DS}}+\mathcal{L}_{\text{int}}$ (2.1)

where the first term is the Standard Model(SM) Lagrangian and the other three
terms are as follows.

$\displaystyle\mathcal{L}_{\text{Seesaw }}$
$\displaystyle=-Y_{\alpha\beta}\overline{L^{\alpha}}\tilde{H}N_{R\beta}-\frac{1}{2}M_{N}\overline{N_{R}^{c}}N_{R}+\text{
h.c. }$ (2.2) $\displaystyle\mathcal{L}_{\text{DS }}$
$\displaystyle=\bar{\chi}\left(i\not\partial-
m_{\chi}\right)\chi+\frac{1}{2}\partial_{\mu}\phi\partial^{\mu}\phi-\frac{1}{2}\mu_{\phi}^{2}\phi^{2}-\frac{1}{4}\lambda_{\phi}{\phi}^{4}-\lambda_{H\phi}\left(H^{\dagger}H\right)\phi^{2}.$
(2.3) $\displaystyle\mathcal{L}_{\text{int }}$
$\displaystyle=-\left(y_{\text{DS }}\phi\bar{\chi}N_{R}+h.c\right).$ (2.4)

In these equations, $L^{\alpha}$ is the left-handed lepton doublets, $H$ the
SM Higgs doublet with $\tilde{H}=i\tau_{2}H^{*}$, $m_{\chi}$ the mass of
$\chi$, $M_{N}$ the mass of the right-handed neutrinos. There are three
generations of right-handed neutrinos and we consider the simplest case where
three generations of right-handed neutrinos are degenerate. After the
electroweak symmetry breaking, the scalar $\phi$ gets a further contribution
to its mass and we have

$m^{2}_{\phi}={\mu}^{2}_{\phi}+\lambda_{H\phi}v_{ew}^{2},$ (2.5)

where $v_{ew}/\sqrt{2}$ is the vacuum expectation of Higgs doublet $H$ after
the electroweak symmetry breaking.

As shown in [9], the right relic abundance of DM can be obtained in this model
by one scenario that the dark fermion $\chi$ serves as the DM candidate and is
produced mainly through the decay of heavier dark scalar $\phi$. In this
scenario, the coupling $y_{DS}$ is assumed to be very small, typically at
order $10^{-12}$. For such a small value of $y_{DS}$, the fermion $\chi$ is
never in thermal equilibrium. We also assume that the coupling
$\lambda_{H\phi}$ is not so small so that the scalar $\phi$ is in the thermal
equilibrium with the SM particles at high temperature. As temperature drops
down to the mass scale of $\phi$, the heavier $\phi$ freezes out and decouples
from the thermal bath. Then, $\phi$ after freeze-out decays into fermion
$\chi$ and the right-handed neutrinos. The right-handed neutrinos would later
decay to leptons and other particles, and we are left with a relic density of
$\chi$ which is to be compared with the DM relic density found in cosmological
observations.

In the previous work [9], only the case for $m_{\phi}>200$ GeV was considered
for obtaining the right relic density. The right relic density of DM can also
be obtained for smaller $\phi$ mass. In this article, we consider here the
parameter space with $m_{\phi}<100$ GeV. The allowed parameter space of
$m_{\phi}$ and $\lambda_{H\phi}$ for $m_{\phi}$ below 100 GeV is shown in Fig.
1(a). We can see that we can obtain the right relic abundance, and
$\lambda_{H\phi}$ is demanded to be in the range $10^{-4}-10^{-2}$ for
$m_{\phi}$ in the range 10-100 GeV. Here we choose $y_{DS}=10^{-12}$,
$m_{N}=10$ GeV, and $m_{\chi}$ in the range 5-1000 MeV. Note that for these
parameters, the $\phi$ can decay sufficiently fast so that the DM relic
density can be produced at sufficient high temperature, as can be shown
similar to the analysis given in [9]. Typically for $m_{\phi}=50$ GeV and
$m_{\chi}=1$ GeV with chosen parameters, the lifetime of $\phi$ is around 0.7s
which would not give much impact on the BBN processes. We note that the peak
at around 75 GeV is due to the opening of the $\phi\phi\rightarrow W^{+}W^{-}$
annihilation channel and the numerical results shown in the figure is obtained
by using the public code MicrOmegas[20, 21] and CALCHEP[22] to solve the
Boltzmann equation and calculate relic abundance. In numerical calculation,
all the contributions from $2\rightarrow 2$ process have been taken into
account. The model implementation for MicrOmegas was done using FeynRules[23]
package.

For $m_{\phi}<m_{H}/2$, Higgs would decay to $\phi$ and would contribute to
the Higgs invisible decay. The experimental data will give an additional
constraint on model parameters $\lambda_{H\phi}$ and $m_{\phi}$. The branching
ratio of Higgs invisible decay is demanded by latest experiment[24] to be less
than 11%. In our model

$\Gamma_{h\rightarrow\phi\phi}=\frac{\lambda_{H\phi}^{2}v_{ew}^{2}}{8\pi
m_{H}}\sqrt{1-\frac{4m_{\phi}^{2}}{m_{H}^{2}}}$ (2.6)

Only the parameters space for
$\Gamma_{h\rightarrow\phi\phi}/(\Gamma_{SM}+\Gamma_{h\rightarrow\phi\phi})\leq
11\%$ is allowed[24]. The gray region in Fig. 1(a) is excluded by the
constraint from the branching ratio of Higgs invisible decay which is demanded
to be less than 11%.

The dark particle $\chi$ does not couple to SM sector directly, so it’s hard
to find relevant signals of $\chi$ on the collider. But the scalar $\phi$ in
the dark sector couples to the SM sector directly with a not so small coupling
constant $\lambda_{H\phi}$, which gives us an opportunity to search this dark
scalar $\phi$ at colliders. In general, the production of this dark sector
particle $\phi$ will lead to a large $E_{\mathrm{T}}^{\mathrm{miss}}$ which
can be a good signal for dark $\phi$ at colliders. Note here that the dark
sector particle $\phi$ is a long-lived particle due to the very small coupling
$y_{DS}$, which can lead to a displaced vertex signal in the collider detector
or long baseline experiment such as MATHUSLA[25], we will leave this topic in
the future work.

Dark sector particle $\phi$ is mainly produced by the gluon gluon fusion
process through Higgs at hadron colliders, or by the ZH($\rightarrow\phi\phi$)
process at future Higgs factory which are shown in Fig. 2(b). The cross
section for gluon gluon fusion production at hadron colliders can be obtained
numerically by MadGraph5_aMC@NLO[26, 27] package and the result is shown in
Fig. 1(b) for $\lambda_{H\phi}=0.01$ and center-of-mass energy $\sqrt{s}=13$
TeV. The total cross section has been rescaled to 48.61 pb which is calculated
at next-to-next-to-next-to-leading order (N3LO) in QCD[29]. We can see in Fig.
1(b) that the cross section drops sharply around $m_{\phi}=m_{H}/2\simeq 62.5$
GeV. This is because for $m_{\phi}\lower
2.71246pt\hbox{$\;\stackrel{{\scriptstyle\textstyle<}}{{\sim}}\;$}62$ GeV, the
production process is basically a $2\to 1$ process with the Higgs produced on-
shell and then decays to $\phi\phi$, as can be seen in the left panel in Fig.
2(b), while for $m_{\phi}\lower
2.71246pt\hbox{$\;\stackrel{{\scriptstyle\textstyle>}}{{\sim}}\;$}63$ GeV on-
shell Higgs can not decay to $\phi\phi$ and the production of $\phi\phi$ shown
in Fig. 2(b) is mediated by off-shell Higgs. So there is a transition from a
$2\rightarrow 1$ gluon gluon fusion production process to a $2\rightarrow 2$
production process at around $m_{\phi}\simeq 62.5$ GeV. We can see in Fig.
1(b) that the cross section for $m_{\phi}>m_{H}/2$ is around five orders of
magnitude smaller than that for $m_{\phi}<m_{H}/2$. For the integrated
luminosity at $300fb^{-1}$, it is estimated that there are only few pairs of
$\phi$ generated for $m_{\phi}>m_{H}/2$ at the LHC. For designed future
$e^{+}e^{-}$ colliders, a rapid decrease of Higgs production cross section
will also happen when $m_{\phi}\simeq 62.5$ GeV which is associated with a
transition from a $2\rightarrow 2$ production process to a $2\rightarrow 3$
production process. Also for $e^{+}e^{-}$ collider running at 250-350 GeV, the
energy conservation will suppress the production of the interested events with
a $\phi$ mass much heavier than $m_{H}$/2. Therefore, in this work, we mainly
focus on the parameter region with $m_{\phi}<m_{H}/2$.

(a)  (b)

Figure 1:  (a) The allowed parameter space of $m_{\phi}$ and $\lambda_{H\phi}$
for $m_{\chi}$ in the range from 5 MeV to 1000 MeV. $M_{N}=10$ GeV and
$y_{DS}=10^{-12}$. The gray region is excluded by Higgs invisible decay. (b)
Cross section for gluon gluon fusion production of $\phi$ pairs at pp
collider, here we choose $\lambda=0.01$ and center-of-mass energy 13 TeV.

(a)  (b)

Figure 2: Major production processes of $\phi$ at pp and $e^{+}e^{-}$
colliders.

##  3 Dark sector search via the Higgstrahlung process at $e^{+}e^{-}$
colliders

For the planned future Higgs factories like CEPC,FCC-ee and ILC with an
integrated luminosity at the order of $ab^{-1}$, around $10^{6}$ Higgs bosons
can be produced, which offers us an opportunity to search for the dark
particle signal. In our model, the major signal process is $e^{+}e^{-}\to ZH$
with the Higgs boson decaying into a pair of $\phi$ which can fly off the
detectors as missing energy. Below we provide a feasibility study for
$2\ell+{E\\!\\!\\!\\!/}$ final states at $e^{+}e^{-}$ colliders with a
detailed MC analysis.

Before presenting our MC analysis, one remarkable aspect of $e^{+}e^{-}$
colliders is worthy of highlighting, i.e. the full four momentum of the
missing energy can be precisely reconstructed at $e^{+}e^{-}$ colliders. In
contrast, only transverse momentum of the missing energy can be reconstructed
at pp colliders.

Meanwhile, another remarkable feature is that the charged leptons $\ell$ both
in the final states and in the initial states can be measured accurately,
which is crucial to reconstruct the mass parameters of mother particles in the
signal and background events. In principle, our analysis can be extended to 2
jets plus missing energy case. Due to the large errors in the jet energy
reconstruction, in this work we restrict to show the 2 charged lepton plus
missing energy case.

In this work, both signal and background events at parton level are generated
by using the MadGraph5_aMC@NLO[26, 27]. The main Feynman diagrams for signal
and background processes are displayed in Fig. 3(b). The signal process and
the main background process with decays are listed below

$\displaystyle e^{+}e^{-}\rightarrow ZH\rightarrow
l^{+}l^{-}+E^{\mathrm{miss}}(\phi\phi)\quad(\textbf{signal process})\,,$ (3.1)
$\displaystyle e^{+}e^{-}\rightarrow ZZ\rightarrow
l^{+}l^{-}+E^{\mathrm{miss}}(\nu\bar{\nu})\ \quad(\textbf{background
process})\,,$

where the missing energy in the signal process originates from Higgs boson’s
invisible decay, while the missing energy in the background process comes from
Z boson decay.

(a)  (b)

Figure 3: (a) Signal process at Higgs factory; (b)Major background process at
Higgs factory.

To fully utilize the experimental precision and to generate the most relevant
events, we introduce the following geometric acceptance cuts at the parton
level. The leptons in the final state must be energetic, be near the central
region of the detectors, and be isolated, i.e. $p_{T}(l)>20$ GeV,
$|\eta(\ell)|<2.5$, and $R(\ell)>0.4$. The leptons’ momentum and energy in the
non-central region has worse accuracy and it is better to neglect those events
with leptons at the non-central region.

All acceptance cuts for the final state $\ell\ell+E^{\text{miss }}$ are
provided in Table 1. Such acceptance cuts are introduced so as to measure the
full momentum of leptons in the final state precisely and to guarantee the
successful reconstruction of the masses of mother particles of missing
children particles.

Geometric Acceptance Cuts
---
Energetic Charged lepton |  $P_{T}(l)>20$ GeV
Central Region | $|\eta(\ell)|<2.5$
Isolation condition | $\Delta R_{\ell\ell}>0.4$
Table 1: Geometric acceptance cuts for the $2\ell+E^{\text{miss }}$ final
states search at Higgs factory are shown.

For signal events, we choose two benchmark points with theoretical parameters
(BP1, $m_{\phi}=30$ GeV, $\lambda_{H\phi}=0.005$) and (BP2, $m_{\phi}=50$ GeV,
$\lambda_{H\phi}=8\times 10^{-4}$) as two show cases to demonstrate how the
recoil mass method can work to separate signal and background events. It
should be mentioned that both of these two benchmark points can lead to the
correct relic density for a cosmological DM candidate.

(a)  (b)

Figure 4:  Missing invariant mass for different center-of-mass frame energy of
the collider for BP1 are displayed: (a) $\sqrt{s}=250$ GeV, around 7100 signal
events and 350 background events between $125\pm 2$ GeV (b) $\sqrt{s}=350$
GeV, around 3850 signal events and 210 background events between $125\pm 2$
GeV. Here we set the bin width to be 2 GeV, which is reasonable when the
errors of the collision energy and momentum of charged leptons are considered.

As it is known that the total energy at the center-of-mass frame for
$e^{+}e^{-}$ colliders can be determined with a remarkable precision (say 0.5
GeV due to the beam spreading) [28], therefore four momentum of the mother
particles (say $H^{0}$ or $Z$ here) can be measured precisely in the central
region of detectors by accurate measurements of the four momentum of the final
state lepton pair $p(l_{1})$ and $p(l_{2})$ from Z decay (say $\frac{\delta
p}{p}=10^{-5}$). By using these measurements, the four momentum of the
invisible particle can be fully reconstructed via the so-called recoil mass
method.

The four momentum of the invisible mother particle denoted as
$P^{\mu}_{inv}=q^{\mu}-P^{\mu}(Z)$ can be used to find its invariant mass,
where $q^{\mu}=(\sqrt{s},0,0,0)$ and
$P^{\mu}(Z)=p(\ell_{1})^{\mu}+p(\ell_{2})^{\mu}$. The invariant mass can be
defined as $m_{\mathrm{inv}}^{\mathrm{miss}}=\sqrt{P_{inv}^{2}}$, which can be
expressed as a function of collision energy and measured lepton four momenta

$m_{\mathrm{inv}}^{\mathrm{miss}}=\sqrt{(\sqrt{s}-E_{l_{1}}-E_{l_{2}})^{2}-(\vec{p}(l_{1})+\vec{p}(l_{2}))^{2}}\,.$
(3.2)

It should be stressed that the error of the observable
$m_{\mathrm{inv}}^{\mathrm{miss}}$ only depends upon the error of $\sqrt{s}$
and the error of ($E_{\ell}$ and $p(\ell)$). For the signal events, the
invariant mass $P^{\mu}_{inv}$ will yield a sharp peak at the Higgs mass
region. In contrast, the events of the major background process
$e^{+}e^{-}\rightarrow Z(\rightarrow l^{+}l^{-})Z(\rightarrow\nu\nu)$ will
lead to a peak at the Z mass region.

Meanwhile, as it is well-known that the decay width of a Higgs boson is around
4 MeV, which is much narrower than that of a Z boson (2.49 GeV). Therefore,
the signal events will enrich at a narrower Higgs boson mass window, i.e.
$|m_{\mathrm{inv}}^{\mathrm{miss}}-m_{H}|\sim 2$ GeV, as demonstrated in Fig.
4 for BP1. Therefore, it clearly shows that it is possible to distinguish a
signal/background event by setting a cut around the Higgs mass window
$|m_{\mathrm{inv}}^{\mathrm{miss}}-125|=2$ GeV.

In Fig. 4 for BP1, even though the events via off-shell Z and $\gamma$ can
contribute significantly even out of the $m_{Z}$ mass window, the signal
events condense in the Higgs boson mass window, and a significance
$\sigma=\frac{N_{S}}{\sqrt{N_{B}}}=400$ for the BP1 can be established
($N_{S}$ and $N_{B}$ denote the number of signal and background events in the
Higgs mass window ). For the collision energy $\sqrt{s}=350$ GeV, a
significance $\sigma=265$ for the BP1 can be established. Compared with
$\sqrt{s}=250$ GeV case, the loss of significance for $\sqrt{s}=350$ GeV can
mainly be attributed to the fact that the cross section of the signal process
is suppressed by collision energy since the signal process is s-channel
dominant, while the cross section of background process is less suppressed
since it is t/u-channel dominant.

(a)  (b)

Figure 5: The recoiled mass for different center-of-mass frame energy of the
collider for BP2 are displayed: (a) $\sqrt{s}=250$ GeV, around 140 signal
events and 350 background events between $125\pm 2$ GeV. Significance of the
signal is around 7.5 $\sigma$ (b) $\sqrt{s}=350$ GeV, around 75 signal events
and 210 background events between $125\pm 2$ GeV. Significance of the signal
is around 4.7 $\sigma$. Here we set the bin width to be 2 GeV, which is
reasonable when the errors of the collision energy and momentum of charged
leptons are considered.

In Fig. 5 for BP2, in the case with $\sqrt{s}=250$ GeV and an assumed 1/ab
dataset, a significance $\sigma=7.48$ can be achieved. While for the case with
$\sqrt{s}=350$ GeV and the 1/ab dataset, a significance $\sigma=4.70$ is
feasible. For both BP1 and BP2, we can arrive at the conclusion that an
$e^{+}e^{-}$ collider with $\sqrt{s}=250$ GeV has a better sensitivity for
signal events. Note here that we can also choose the Higgs boson mass window
to be $|m_{\mathrm{inv}}^{\mathrm{miss}}-125|=1$ GeV. Comparing with the 2 GeV
window, this scheme cut more background events and lead to a higher background
rejection rate, so a higher significance can be obtained. While for the 2 GeV
window, a higher signal efficiency will be achieved.

We can extend our analysis for BP1 and BP2 to other points in the theoretical
parameter space formed by $\lambda_{H\phi}$ and $m_{\phi}$ as demonstrated in
Table. 2. We can obtain the lower bounds of the values of coupling
$\lambda_{H\phi}$ which can be probed at $3\sigma$ and $5\sigma$ significance
level for integrated luminosity of both $\sqrt{s}=250$ GeV and $\sqrt{s}=350$
GeV and for some different values of mass parameter $m_{\phi}$.

$\sqrt{s}=250$ GeV
---
|  $m_{\phi}$ = 30 GeV |  $m_{\phi}$ = 50 GeV |  $m_{\phi}$ = 60 GeV
|  3 $\sigma$ |  5 $\sigma$ |  3 $\sigma$ |  5 $\sigma$ |  3 $\sigma$ |  5 $\sigma$
5 $ab^{-1}$ | $4.30\times 10^{-4}$ | $5.60\times 10^{-4}$ | $5.20\times 10^{-4}$ | $6.80\times 10^{-4}$ | $7.60\times 10^{-4}$ | $9.80\times 10^{-4}$
1 $ab^{-1}$ | $6.60\times 10^{-4}$ | $8.60\times 10^{-4}$ | $7.90\times 10^{-4}$ | $1.03\times 10^{-3}$ | $1.16\times 10^{-3}$ | $1.51\times 10^{-3}$
500 $fb^{-1}$ | $7.90\times 10^{-4}$ | $1.10\times 10^{-3}$ | $9.50\times 10^{-4}$ | $1.26\times 10^{-3}$ | $1.40\times 10^{-3}$ | $1.84\times 10^{-3}$
$\sqrt{s}=350$ GeV
5 $ab^{-1}$ | $5.20\times 10^{-4}$ | $6.70\times 10^{-4}$ | $6.30\times 10^{-4}$ | $8.10\times 10^{-4}$ | $9.20\times 10^{-4}$ | $1.19\times 10^{-3}$
1 $ab^{-1}$ | $7.90\times 10^{-4}$ | $1.05\times 10^{-3}$ | $9.60\times 10^{-4}$ | $1.26\times 10^{-3}$ | $1.40\times 10^{-3}$ | $1.85\times 10^{-3}$
500 $fb^{-1}$ | $9.60\times 10^{-4}$ | $1.27\times 10^{-3}$ | $1.15\times 10^{-3}$ | $1.53\times 10^{-3}$ | $1.71\times 10^{-3}$ | $2.23\times 10^{-3}$
Table 2: The minimum values of coupling $\lambda_{H\phi}$ which can be probed
at 3 $\sigma$ and 5 $\sigma$ confidence levels with different masses of
$m_{\phi}$ with collision energy $\sqrt{s}=250$ GeV and $\sqrt{s}=350$ GeV are
shown.

When collision energies are fixed, the minimum integrated luminosity to obtain
$3\sigma$ and $5\sigma$ significance for different values of the coupling
$\lambda_{H\phi}$ and mass $m_{\phi}$ are shown in Table 3.

From the results given in Table 3, we can observe that the future $e^{+}e^{-}$
colliders can probe the signal process with the coupling $\lambda_{H\phi}$
down to $6\times 10^{-4}$ or so. With the results of more points in
theoretical parameter space, we can conclude that $e^{+}e^{-}$ colliders with
collision energy $\sqrt{s}=250$ can have a better sensitivity to our model.

$\sqrt{s}=250$ GeV
---
|  $m_{\phi}=30$ GeV |  $m_{\phi}=50$ GeV |  $m_{\phi}=60$ GeV
$\lambda_{H\phi}$ | $6\times 10^{-4}$ | $8\times 10^{-4}$ | $1\times 10^{-3}$ | $6\times 10^{-4}$ | $8\times 10^{-4}$ | $1\times 10^{-3}$ | $8\times 10^{-4}$ | $1\times 10^{-3}$ | $1.2\times 10^{-3}$
3 $\sigma$ |  1.41 $ab^{-1}$ |  480 $fb^{-1}$ |  230 $fb^{-1}$ |  2.90 $ab^{-1}$ |  940 $fb^{-1}$ |  430 $fb^{-1}$ |  4.10 $ab^{-1}$ |  1.70 $ab^{-1}$ |  890 $fb^{-1}$
5 $\sigma$ |  3.80 $ab^{-1}$ |  1.27 $ab^{-1}$ |  560 $fb^{-1}$ |  7.80 $ab^{-1}$ |  2.60 $ab^{-1}$ |  1.13 $ab^{-1}$ |  11.20 $ab^{-1}$ |  4.70 $ab^{-1}$ |  2.38 $ab^{-1}$
$\sqrt{s}=350$ GeV
$\lambda_{H\phi}$ | $6\times 10^{-4}$ | $8\times 10^{-4}$ | $1\times 10^{-3}$ | $6\times 10^{-4}$ | $8\times 10^{-4}$ | $1\times 10^{-3}$ | $8\times 10^{-4}$ | $1\times 10^{-3}$ | $1.2\times 10^{-3}$
3 $\sigma$ |  2.9 $ab^{-1}$ |  990 $fb^{-1}$ |  420 $fb^{-1}$ |  5.80 $ab^{-1}$ |  1.97 $ab^{-1}$ |  860 $fb^{-1}$ |  8.40 $ab^{-1}$ |  3.60 $ab^{-1}$ |  1.80 $ab^{-1}$
5 $\sigma$ |  7.85 $ab^{-1}$ |  2.70 $ab^{-1}$ |  1.16 $ab^{-1}$ |  16.10 $ab^{-1}$ |  5.30 $ab^{-1}$ |  2.31 $ab^{-1}$ |  23.00 $ab^{-1}$ |  9.80 $ab^{-1}$ |  4.80 $ab^{-1}$
Table 3: The minimum integrated luminosity for 3 $\sigma$ and 5 $\sigma$
confidence level for different masses $m_{\phi}$ and values of coupling
$\lambda_{H\phi}$ with two collision energy $\sqrt{s}=250$ GeV and
$\sqrt{s}=350$ GeV are provided.

##  4 Dark sector particle $\phi$ search at LHC

For the sake of comparison, we study the feasibility to detect dark particle
$\phi$ at proton-proton colliders. We focus on both mono-jet and
$2\ell+E_{\mathrm{T}}^{\mathrm{miss}}$ final states at pp colliders.

It is well-known that the mono-jet process serves as a unique probe to search
for the DM production processes which lead to a large
$E_{\mathrm{T}}^{\mathrm{miss}}$ and one energetic jet in the final state. In
our model, the signal process is $pp\to H+j$ with Higgs boson decaying into a
$\phi$ pair. The dominant background process is $pp\to Z+j$ with the Z boson
decaying into a neutrino pair. The detailed signal/background processes and
decays are given below

$\displaystyle pp\rightarrow H+\text{jet}\rightarrow
E_{\mathrm{T}}^{\text{miss }}(\phi\phi)+\text{jet}\quad(\textbf{signal})$
(4.1) $\displaystyle pp\rightarrow Z+\text{jet}\rightarrow
E_{\mathrm{T}}^{\text{miss }}(\nu\bar{\nu})+\text{jet}\
\quad(\textbf{background})$

We generate the relevant the mono-jet events at parton level by using
MadGraph5_aMC@NLO[26, 27] generator for signal and background processes. For
simplicity, only one energetic jet in the final state is considered. The PDF
set used for the generation is NNPDF23LO, and the cross section for the
process with the hard jet in the final state is obtained by using a dataset at
parton level (no hadronization or no detector effects has been considered) .
We set the following selection cuts to the leading jet: $p_{T}(j)>150$ GeV and
$|\eta(j)|<2.4$.

In Fig. 6, the distributions of the transverse momentum of the leading jet
$p_{T}(j)$ of signal/background events are shown. For the signal events, we
have chosen the theoretical parameters $m_{\phi}=30$ GeV in the left panel and
$m_{\phi}=50$ GeV in the right panel with $\lambda_{H\phi}=5\times 10^{-3}$.
The blue(red) bins represent the distribution of signal(background) events.

To examine the reliability of our MC simulation, the experimental data are
denoted by the black points, which is normalized to the experimental results
for mono-jets search at LHC with integrated luminosity $139fb^{-1}$[30, 31]
and simply rescaled to $300fb^{-1}$. It is noticed that there exists a small
gap between the experimental data and the results of MC simulation which can
be attributed to the W+jets processes and the NLO contribution, as has been
taken into account by [30].

We can find in both Fig. 6(a) and Fig. 6(b) that the background events are
approximately 3 orders of magnitude larger than signal events for each a bin,
even when the NLO correction for the $gg\to H^{0}j$ is taken into account.
This indicates that it is really challenging for the LHC to find the signal
process from the mono-jets events with a total integrated luminosity
$300fb^{-1}$ after run3.

Selection cuts for mono-jet final states
---
Leading jet |  $p_{T}(j)>150$ GeV and $|\eta|<2.4$
Selection cuts for $ZH\rightarrow\ell\ell+E_{\mathrm{T}}^{\text{miss}}$ final
states
Two leptons |  leading(subleading) $p_{T}(\ell)$ ¿ 30(20) GeV, $|\eta(\ell)|<2.5$
$m_{\ell\ell}$ |  $76$ GeV $<m_{\ell\ell}<106$ GeV
$E_{\mathrm{T}}^{\text{miss }}$ |  $E_{\mathrm{T}}^{\text{miss }}>90$GeV
$\Delta R_{\ell}$ | $\Delta R_{\ell}\geq 0.4$
Table 4: Selection cuts for mono-jet and $2\ell+E_{\mathrm{T}}^{\text{miss }}$
final states at pp colliders are presented.

(a)  (b)

Figure 6: (a) $p_{T}(j)$ distribution of signal events and background events
in the monojet processes. The black points represent the experimental data.
(a)$\lambda_{H\phi}=5\times 10^{-3}$, $m_{\phi}=30$ GeV.
(b)$\lambda_{H\phi}=5\times 10^{-3}$, $m_{\phi}=50$ GeV.

The process for production of Higgs boson in association with a Z boson with
the Higgs boson invisibly decay is also a good channel to search for new
physics at pp colliders. The signal process and major background process is
similar with the $e^{+}e^{-}$ collider case with proton instead of electron in
the initial state.

$\displaystyle pp\rightarrow ZH\rightarrow
l^{+}l^{-}+E_{\mathrm{T}}^{\mathrm{miss}}(\phi\phi)\quad(\textbf{signal})$
(4.2) $\displaystyle pp\rightarrow ZZ\rightarrow
l^{+}l^{-}+E_{\mathrm{T}}^{\mathrm{miss}}(\nu\bar{\nu})\
\quad(\textbf{background})$

Since we can not reconstruct the full momentum of the missing particles in the
final state at pp colliders, instead we can only use the transverse missing
momentum $E_{\mathrm{T}}^{\mathrm{miss}}$ as an observable.

In Fig. 7, we show the distribution of $E_{\mathrm{T}}^{\mathrm{miss}}$ for
signal/background in blue/red at LO and compare the experimental data with the
Monte Carlo results. Here the ATLAS results for invisibly decaying Higgs boson
in association with a Z boson data at 36.1 $fb^{-1}$[32] are quoted and scaled
to 300 $fb^{-1}$.

In our analysis, we demand that the leptons $p_{T}(\ell)>20$ GeV and the
leading lepton $p_{T}(\ell)>30$ GeV. Also
$E_{\mathrm{T}}^{\mathrm{miss}}>90\mathrm{GeV}$ is demanded and the invariant
mass of two charged leptons is required to be in the range $76$ GeV
$<m_{\ell\ell}<106~{}\mathrm{GeV}$ to reject background processes with two
leptons that do not originate from the decay of a Z boson. All the cuts we use
are shown in Table 4.

Note here that the result has a deviation of about 10% from the experimental
data where the NLO corrections have been taken into account. We can find in
the figure that the signal is 2 orders of magnitude smaller than the
background. So the signal events for dark sector particle $\phi$ is also very
hard to be found from ZH process at pp colliders.

Comparing these results with those at $e^{+}e^{-}$ colliders given in Sec. 3,
we can conclude that $e^{+}e^{-}$ colliders have more advantages in searching
for such dark sector particle $\phi$ from Higgs invisible decay.

(a)  (b)

Figure 7: Distribution of $E_{\mathrm{T}}^{\mathrm{miss}}$ for signal process
$(H\rightarrow\phi\phi)$ \+ $(Z\rightarrow l^{+}l^{-})$ and background process
$(Z\rightarrow\nu\bar{\nu})$ \+ $(Z\rightarrow l^{+}l^{-})$. The black points
represent the experimental data. (a)$\lambda_{H\phi}=5\times 10^{-3}$,
$m_{\phi}=30$ GeV. (b)$\lambda_{H\phi}=5\times 10^{-3}$, $m_{\phi}=50$ GeV.

##  5 Conclusions

In summary, we continued the research on the scenario that the DM can be
produced by the decay of other heavier particles in the dark sector. For the
specific model which include one singlet scalar $\phi$ and one singlet fermion
$\chi$ in the dark sector, we have studied the possibility for a lighter dark
sector particle $\phi$ with $m_{\phi}<100$ GeV to obtain the right relic
density. We find this scenario with lower $m_{\phi}$ is also appropriate in a
wide range of parameter space.

We have studied the possibility to find the signal of the dark sector particle
$\phi$ on colliders. We mainly focus on the region that $m_{\phi}<m_{H}/2$ and
have studied the production rate and possible signals on the colliders. For
future Higgs factory running at 250-350 GeV, the dark particle $\phi$ is
mainly produced from ZH process and will result in the dilepton final state
with large $E_{\mathrm{T}}^{\mathrm{miss}}$. This kind of signature can also
be generated by SM process from diboson(ZZ) decay, which is the domant
background events. From reconstructing the full missing three momentum and
energy we can obtain the invariant mass of the missing part which give us a
powerful tool to distinguish the signal and background events. Using the
missing invariant mass as an observable, we can probe this kind of dark sector
signal with the coupling of $\phi$ and Higgs down to $5\times 10^{-4}$ at
future Higgs factories. We have studied the signal and background processes
for some benchmark points in our model and have found that the signal has
strong significance for a wide range of the parameter space. We have also
studied the sensitivity for searching this dark sector particle signal at
Higgs factories with different center-of-mass frame energy. We find that a
center-of-mass frame energy 250 GeV is slightly better to find the dark sector
signal.

For comparison, we have also studied the possibility of looking for the same
dark sector signals on the pp colliders. Higgs production and then decay to
$\phi$ pair associated with an energetic jet in our model will lead to the
Mono-jet signals, for which $Z\rightarrow\nu\bar{\nu}+$ jets process is the
major SM background. We have studied $p_{T}(j)$ distribution for signal and
background and have found that the signal is 3 orders of magnitude smaller the
background. We have further studied the dark sector particle produced by ZH
process at pp colliders, where $(ZZ\rightarrow
l^{+}l^{-})+E_{\mathrm{T}}^{\mathrm{miss}}$ is the major background in the SM.
Since we can not reconstruct the full missing momentum, we could only choose
transverse missing energy as an observable. We find that for some benchmark
points of the model, the signal is 2 orders of magnitude smaller than the
background in the $E_{\mathrm{T}}^{\mathrm{miss}}$ plot. We find that for the
missing energy signals considered in this article the pp collider is not a
good place to detect the signal of the dark sector $\phi$ particle.

Note that the dark sector particle $\phi$ may lead to signals of displaced
vertex in the detector, which depends upon the dark sector and right-handed
neutrinos coupling $y_{DS}$. We will leave this topic in the future work.

## Acknowledgements

W. Liao is supported by National Natural Science Foundation of China under the
grant No. 11875130. Q.S.Yan is supported by the National Natural Science
Foundation of China under the grant No. 11475180 and No. 11875260.

## References

  * [1] XENON Collaboration, E. Aprile et al., “Dark Matter Search Results from a One Ton-Year Exposure of XENON1T,” Phys. Rev. Lett. 121 no. 11, (2018) 111302, arXiv:1805.12562 [astro-ph.CO].
  * [2] X. Cui, A. Abdukerim, W. Chen, X. Chen, Y. Chen, B. Dong, D. Fang, C. Fu, K. Giboni, F. Giuliani, et al., “Dark matter results from 54-ton-day exposure of PandaX-II experiment,” Phys. Rev. Lett. 119 no. 18, (2017) 181302.
  * [3] PandaX-II Collaboration, Q. Wang et al., “Results of dark matter search using the full PandaX-II exposure,” Chin. Phys. C 44 no. 12, (2020) 125001, arXiv:2007.15469 [astro-ph.CO].
  * [4] LUX Collaboration, D. Akerib et al., “Results from a search for dark matter in the complete LUX exposure,” Phys. Rev. Lett. 118 no. 2, (2017) 021303, arXiv:1608.07648 [astro-ph.CO].
  * [5] M. D. Campos, F. S. Queiroz, C. E. Yaguna, and C. Weniger, “Search for right-handed neutrinos from dark matter annihilation with gamma-rays,” JCAP 07 (2017) 016, arXiv:1702.06145 [hep-ph].
  * [6] IceCube Collaboration, M. Aartsen et al., “Multipole analysis of IceCube data to search for dark matter accumulated in the Galactic halo,” Eur. Phys. J. C 75 no. 99, (2015) 20, arXiv:1406.6868 [astro-ph.HE].
  * [7] CMS Collaboration, V. Khachatryan et al., “Search for dark matter, extra dimensions, and unparticles in monojet events in proton–proton collisions at $\sqrt{s}=8$ TeV,” Eur. Phys. J. C 75 no. 5, (2015) 235, arXiv:1408.3583 [hep-ex].
  * [8] ATLAS Collaboration, G. Aad et al., “Search for dark matter candidates and large extra dimensions in events with a jet and missing transverse momentum with the ATLAS detector,” JHEP 04 (2013) 075, arXiv:1210.4491 [hep-ex].
  * [9] Y. Cheng and W. Liao, “Light dark matter from dark sector decay,” Phys. Lett. B 815 (2021) 136118, arXiv:2012.01875 [hep-ph].
  * [10] V. Gonzalez Macias and J. Wudka, “Effective theories for Dark Matter interactions and the neutrino portal paradigm,” JHEP 07 (2015) 161, arXiv:1506.03825 [hep-ph].
  * [11] M. Blennow, E. Fernandez-Martinez, A. Olivares-Del Campo, S. Pascoli, S. Rosauro-Alcaraz, and A. V. Titov, “Neutrino Portals to Dark Matter,” Eur. Phys. J. C 79 no. 7, (2019) 555, arXiv:1903.00006 [hep-ph].
  * [12] S. Baek, P. Ko, and W.-I. Park, “Singlet Portal Extensions of the Standard Seesaw Models to a Dark Sector with Local Dark Symmetry,” JHEP 07 (2013) 013, arXiv:1303.4280 [hep-ph].
  * [13] M. Escudero, N. Rius, and V. Sanz, “Sterile Neutrino portal to Dark Matter II: Exact Dark symmetry,” Eur. Phys. J. C 77 no. 6, (2017) 397, arXiv:1607.02373 [hep-ph].
  * [14] V. González-Macías, J. I. Illana, and J. Wudka, “A realistic model for Dark Matter interactions in the neutrino portal paradigm,” JHEP 05 (2016) 171, arXiv:1601.05051 [hep-ph].
  * [15] P. Bandyopadhyay, E. J. Chun, R. Mandal, and F. S. Queiroz, “Scrutinizing Right-Handed Neutrino Portal Dark Matter With Yukawa Effect,” Phys. Lett. B 788 (2019) 530–534, arXiv:1807.05122 [hep-ph].
  * [16] M. Chianese and S. F. King, “The Dark Side of the Littlest Seesaw: freeze-in, the two right-handed neutrino portal and leptogenesis-friendly fimpzillas,” JCAP 09 (2018) 027, arXiv:1806.10606 [hep-ph].
  * [17] M. Chianese, B. Fu, and S. F. King, “Minimal Seesaw extension for Neutrino Mass and Mixing, Leptogenesis and Dark Matter: FIMPzillas through the Right-Handed Neutrino Portal,” JCAP 03 (2020) 030, arXiv:1910.12916 [hep-ph].
  * [18] M. Chianese, B. Fu, and S. F. King, “Interplay between neutrino and gravity portals for FIMP dark matter,” arXiv:2009.01847 [hep-ph].
  * [19] P. Bandyopadhyay, E. J. Chun, and R. Mandal, “Feeble neutrino portal dark matter at neutrino detectors,” JCAP 08 (2020) 019, arXiv:2005.13933 [hep-ph].
  * [20] G. Belanger, F. Boudjema, A. Pukhov, and A. Semenov, “micrOMEGAs$\\_$3: A program for calculating dark matter observables,” Comput. Phys. Commun. 185 (2014) 960–985, arXiv:1305.0237 [hep-ph].
  * [21] G. Bélanger, F. Boudjema, A. Goudelis, A. Pukhov, and B. Zaldivar, “micrOMEGAs5.0 : Freeze-in,” Comput. Phys. Commun. 231 (2018) 173–186, arXiv:1801.03509 [hep-ph].
  * [22] A. Belyaev, N. D. Christensen, and A. Pukhov, “CalcHEP 3.4 for collider physics within and beyond the Standard Model,” Comput. Phys. Commun. 184 (2013) 1729–1769, arXiv:1207.6082 [hep-ph].
  * [23] A. Alloul, N. D. Christensen, C. Degrande, C. Duhr, and B. Fuks, “FeynRules 2.0 - A complete toolbox for tree-level phenomenology,” Comput. Phys. Commun. 185 (2014) 2250–2300, arXiv:1310.1921 [hep-ph].
  * [24] ATLAS Collaboration, “Combination of searches for invisible Higgs boson decays with the ATLAS experiment”, ATLAS-CONF-2020-052.
  * [25] MATHUSLA Collaboration, H. Lubatti et al., “Explore the lifetime frontier with MATHUSLA,” JINST 15 no. 06, (2020) C06026, arXiv:1901.04040 [hep-ex].
  * [26] J. Alwall, M. Herquet, F. Maltoni, O. Mattelaer, and T. Stelzer, “MadGraph 5 : Going Beyond,” JHEP 06 (2011) 128, arXiv:1106.0522 [hep-ph].
  * [27] J. Alwall, R. Frederix, S. Frixione, V. Hirschi, F. Maltoni, O. Mattelaer, H. S. Shao, T. Stelzer, P. Torrielli, and M. Zaro, “The automated computation of tree-level and next-to-leading order differential cross sections, and their matching to parton shower simulations,” JHEP 07 (2014) 079, arXiv:1405.0301 [hep-ph].
  * [28] T. Behnke, J. E. Brau, P. N. Burrows, J. Fuster, M. Peskin, M. Stanitzki, Y. Sugimoto, S. Yamada, H. Yamamoto and H. Abramowicz, et al. arXiv:1306.6329 [physics.ins-det].
  * [29] M. Cepeda et al., “Report from Working Group 2: Higgs Physics at the HL-LHC and HE-LHC,” CERN Yellow Rep. Monogr. 7 (2019) 221–584, arXiv:1902.00134 [hep-ph].
  * [30] ATLAS Collaboration, “Constraining the Dark Sector with the monojet signature in the ATLAS experiment”, ATL-PHYS-PUB-2021-020.
  * [31] ATLAS Collaboration, G. Aad et al., “Search for new phenomena in events with an energetic jet and missing transverse momentum in $pp$ collisions at $\sqrt{s}$ =13 TeV with the ATLAS detector,” Phys. Rev. D 103 no. 11, (2021) 112006, arXiv:2102.10874 [hep-ex].
  * [32] ATLAS Collaboration, M. Aaboud et al., “Search for an invisibly decaying Higgs boson or dark matter candidates produced in association with a $Z$ boson in $pp$ collisions at $\sqrt{s}=$ 13 TeV with the ATLAS detector,” Phys. Lett. B 776 (2018) 318–337, arXiv:1708.09624 [hep-ex].


These authors contributed equally to this work.

These authors contributed equally to this work.

[1]Carmine Giordano

1]Department of Aerospace Science and Technology, Politecnico di Milano, via
La Masa, 34, Milan, 20156, Italy

# Initial Trajectory Assessment of the RAMSES Mission to (99942) Apophis

Andrea C. Morelli<EMAIL_ADDRESS>Alessandra Mannocchi
<EMAIL_ADDRESS><EMAIL_ADDRESS>Fabio Ferrari
<EMAIL_ADDRESS>Francesco Topputo<EMAIL_ADDRESS>[

###### Abstract

(99942) Apophis is a potentially hazardous asteroid that will closely approach
the Earth on April 13, 2029. Although the likelihood of an impact has been
ruled out, this close encounter represents a unique opportunity for planetary
science and defense. By investigating the physical and dynamical changes
induced by this interaction, valuable insights into asteroid cohesion,
strength, and internal structure can be obtained. In light of these
circumstances, a fast mission to Apophis holds great scientific importance and
potential for understanding potentially hazardous asteroids. To this aim, ESA
proposed the mission RAMSES (Rapid Apophis Mission for SEcurity and Safety) to
reach Apophis before its close encounter. In this context, the paper focuses
on the reachability analysis of (99942) Apophis, examining thousands of
trajectories departing from Earth and reaching the asteroid before the fly-by,
using a low-thrust spacecraft. A two-layer approach combining direct
sequential convex programming and an indirect method is employed for fast and
reliable trajectory optimization. The results reveal multiple feasible launch
windows and provide essential information for mission planning and system
design.

###### keywords:

Apophis, Trajectory Design, Convex Optimization, Indirect Methods, RAMSES
mission

##  1 Introduction

(99942) Apophis is a potentially hazardous asteroid with a diameter of about
370 metres that caused a brief period of concern in December 2004 when initial
observations indicated a probability up to 2.7% that it would hit the Earth on
April 13, 2029 [4]. Subsequent observations, however, improved predictions and
ruled out the possibility of impact [9]. Nevertheless, Apophis will pass
within 31,000 km of Earth’s surface, closer than geosynchronous satellites
[17].
While it is certain that Apophis will miss the Earth, the exact consequences
on the asteroid itself remain uncertain [26], with the main question lying in
our limited knowledge of Apophis’ internal structure. As a matter of fact, the
close encounter with the Earth will subject Apophis to significant tidal
torques. Earth’s gravitational interactions with Apophis will likely modify
its physical and dynamical properties, providing valuable insights into
asteroid cohesion and strength [18, 20]. Depending on those, the asteroid
could experience different outcomes, from measurable seismic waves and real-
time surface disturbance [7], to local surface effects [16], up to the
complete surface reshaping [24, 12]. For this reason, this encounter presents
an unprecedented planetary defense and science opportunity, and a mission to
Apophis has the potential to revolutionize our understanding of the internal
structure of asteroids, which is crucial for effective planetary defense
missions [1].
As direct consequence of this celestial event, NASA has recently announced
that Apophis will be the selected target of Osiris-REx extended mission [6].
On the other hand, ESA proposed the mission RAMSES (Rapid Apophis Mission for
SEcurity and Safety) having the aim to reach Apophis before the Earth close
encounter111https://esastar-publication-
ext.sso.esa.int/ESATenderActions/details/56945 (last accessed: July 7, 2023).
Its objective is to accompany the asteroid during the fly-by in order to
determine the physical and dynamical changes induced by the interaction with
the gravity of the Earth.
In the context of the RAMSES mission, the present work aims to scrutinize the
reachability of the Apophis asteroid. Reachability analysis is an essential
task in the preliminary assessment of asteroid missions [21, 13, 19]. In this
work, all the possible trajectories departing from the Earth and reaching
Apophis before the Earth fly-by are computed and their overall costs are
estimated. The best suitable options are identified, and relevant figures that
support the system design are evaluated. The variation of the required
propellant mass with the departure date and Time of Flight (ToF) have to be
produced, with thousands of trajectories to be computed. For this reason, a
fast and reliable approach is necessary, in order to have the full envelope of
the transfers to Apophis in a relatively small time. To this aim, a two-layer
approach has been employed. It first exploits the flexibility of a sequential
convex programming algorithm and then an indirect method, the latter
guaranteeing the optimality of the solutions.
The paper is structured as follows. Section 2 presents an overview of the
assumptions made on the mission and on the platform. Section 3 shows the
approach used to fast solve the Apophis reachability problem, with details on
the employed methodology. The results are presented in Section 4. Conclusions
are drawn in Section 5.

##  2 Mission overview and assumptions

RAMSES is a small satellite mission that aims at the characterization of the
dynamical and physical properties of asteroid Apophis before its close
encounter with the Earth, on April 13, 2029, and their changes during the fly-
by. Considering the Earth fly-by’s close temporal proximity, the whole mission
design and implementation is strictly driven by time constraints. As a matter
of fact, in order to properly map the pre-fly-by characteristics of Apophis,
the spacecraft shall reach the asteroid at least 2 months before the
encounter. On the other hand, a time frame of at least 3 years is required to
design and manufacture the spacecraft, placing the earliest launch date on
November 1, 2026. These stringent time constraints open the possibility to use
the electric thruster as main propulsion system. Indeed, even though electric
propulsion systems require more complex operations, they allow wider departure
and arrival windows. In addition, the amount of propellant they require is
significantly lower and single-point failure effects are not as catastrophic
as they would be when chemical propulsion is considered, since maneuvers are
not strongly time-driven. Overall, electric propulsion is a suitable solution
for a fast mission to Apophis.
The probe wet initial mass $m_{0}$ has been assumed to be 500, with 73 of
available propellant $m_{p}$. The dry mass represents more than the 85% of the
total mass, corresponding to a typical value for this class of spacecraft. The
selected values for the maximum-thrust-to-initial-mass ratio
$T_{\text{max}}/m_{0}$ and for the specific impulse $I_{\mathrm{sp}}$ for the
electric propulsion subsystem are 1.2 e-4 and 1500, respectively, which are
compatible with state-of-the-art Hall-effect thrusters [5]. Accordingly, the
maximum value for the thrust $T_{\mathrm{max}}$ is 60.
The spacecraft is inserted directly into an interplanetary transfer to Apophis
by a launcher. Considering the characteristics of Ariane 6.2 and 6.4
launchers222Ariane 6 User Manual I2R0: https://www.arianespace.com/wp-
content/uploads/2021/03/Mua-6_Issue-2_Revision-0_March-2021.pdf (last
accessed: July 7, 2023), the infinity velocity $v_{\infty}$ at the Earth
sphere-of-influence interface has been constrained to be lower than 4, with
free declination $\delta$ and right ascension $\alpha$ in the J2000 reference
frame. The analysis considers a two-body problem dynamics, with the Sun as the
central body. The main assumptions for the Apophis reachability analysis are
listed in Table 1.

Time Constraints | Earliest departure date | November 1, 2026
---|---|---
| Max. Time of Flight $\text{ToF}_{\text{max}}$ | 800 days
| Latest arrival | February 13, 2029
Spacecraft | $m_{0}$ |  500
| Propulsion | Continuous
| $T_{\text{max}}/m_{0}$ |  1.2 e-4
| $I_{\mathrm{sp}}$ |  1500
| $T_{\mathrm{max}}$ |  60
| $m_{p}$ |  73
Launcher | $v_{\infty}$ |  $\leq$ 4
| $\alpha$ | $\in[-90,+90]\deg$
| $\delta$ | $\in[-180,+180]\deg$
Table 1: Search space assumptions to perform the reacheability analysis.

##  3 Approach

The procedure followed to analyze the reachability of Apophis has to be fast
and reliable, since thousands of fuel-optimal (FO) trajectories must be
computed. The output of this analysis is a porkchop plot that shows the
propellant mass required to reach Apophis and the reachability of the target
for different departure dates and times of flight, considering the assumptions
in Table 1. This assessment has the aim to evaluate the launching conditions
and compute the trajectories spanning the whole temporal search space. The
approach for this analysis exploits two layers:

  1. 1.

A direct sequential convex programming (SCP) algorithm [15], used for each
departure date and time of flight. It gives a first assessment of the optimal
launching conditions and trajectory, as it allows to easily include the free
infinite velocity given by the launcher;

  2. 2.

The outputs of the first step are used to feed an indirect method [25], able
to guarantee optimality of the solution in terms of propellant mass.

This two-layer approach is employed to exploit the easy handling of the free
launching conditions of the first method and the guaranteed optimality of the
second one. The minimum-fuel problem is solved with fixed initial and final
boundary conditions, and time of flight. The initial boundary condition is
variable in the first step of the methodology and fixed in the second.

###  3.1 Convex Optimization with Free Infinite Velocity

In the first step, the nonconvex low-thrust trajectory optimization problem is
solved by considering a sequence of convex subproblems whose solutions
eventually converge, under certain hypotheses, to the solution of the original
one [14]. At each iteration, the following problem is solved [23, 10]

$\displaystyle\underset{\mathbf{u}(t)}{\text{minimize}}\hphantom{abc}$
$\displaystyle-w(t_{f})+\lambda\max(0,\eta(t))+\lambda\|\bm{\nu}(t)\|_{1}$
(1a) subject to:
$\displaystyle\dot{\mathbf{x}}(t)=\mathbf{f}(\mathbf{\bar{x}}(t),\mathbf{\bar{u}}(t))+\mathbf{A}(\mathbf{\bar{x}}(t))\,(\mathbf{x}(t)-\mathbf{\bar{x}}(t))+\mathbf{B}(\mathbf{u}(t)-\mathbf{\bar{u}}(t))+\bm{\nu}(t)$
(1b) $\displaystyle\Gamma(t)\leq
T_{\text{max}}\textrm{e}^{-\bar{w}(t)}\left(1-w(t)+\bar{w}(t)\right)+\eta(t)$
(1c) $\displaystyle\left\lVert\bm{\tau}(t)\right\rVert_{2}\leq\Gamma(t)$ (1d)
$\displaystyle\left\lVert\mathbf{x}(t)-\mathbf{\bar{x}}(t)\right\rVert_{1}\leq
R$ (1e)
$\displaystyle\mathbf{r}(t_{i})=\mathbf{r}_{i},\;\mathbf{v}(t_{i})=\mathbf{v}_{i},\;w(t_{i})=w_{i}$
(1f)
$\displaystyle\mathbf{r}(t_{f})=\mathbf{r}_{f},\;\mathbf{v}(t_{f})=\mathbf{v}_{f}$
(1g)
$\displaystyle\mathbf{x}_{l}\leq\mathbf{x}\leq\mathbf{x}_{u},\;\mathbf{u}_{l}\leq\mathbf{u}\leq\mathbf{u}_{u}$
(1h)

where $\mathbf{x}=[\mathbf{r},\mathbf{v},w]$ and
$\mathbf{u}=[\tau_{x},\tau_{y},\tau_{z},\Gamma]=[\bm{\tau},\Gamma]$ are the
state and control variables, respectively. The quantities $\eta(t)$ and
$\bm{\nu}(t)$ are slack variables to avoid the so-called artificial
infeasibility, and the constant parameter $\lambda$ in the objective function
is a user-defined weight. The times $t_{i}$ and $t_{f}$ are the initial and
final transfer times. $R$ is the radius of the trust region assuring that the
convexification of the problem is valid. The reader can refer to previous
works [2] for a more comprehensive explanation of the SCP technique. In Eq.
1b, the matrices $\mathbf{A}$ and $\mathbf{B}$ are defined as

$\mathbf{A}(\mathbf{\bar{x}}(t))\mathrel{\vbox{\hbox{\scriptsize.}\hbox{\scriptsize.}}}=\frac{\partial\mathbf{f}}{\partial\mathbf{x}}\Bigg{\rvert}_{\mathbf{\bar{x}}(t)},\;\mathbf{B}\mathrel{\vbox{\hbox{\scriptsize.}\hbox{\scriptsize.}}}=\frac{\partial\mathbf{f}}{\partial\mathbf{u}}\Bigg{\rvert}_{\mathbf{\bar{u}}(t)}\;$
(2)

where

$\mathbf{f}(\mathbf{x},\mathbf{u})=\begin{bmatrix}\mathbf{v}(t)\\\
-\frac{\mu}{r^{3}}\mathbf{r}(t)+\bm{\tau}(t)\\\
-\frac{\Gamma(t)}{I_{\text{sp}}g_{0}}\end{bmatrix}$ (3)

with $r=\left\lVert\mathbf{r}(t)\right\rVert_{2}$. The problem in Eqs. 1 is a
Second-Order Cone Program (SOCP) and can be handled by efficient convex
solvers. In this work, it has been solved using an Hermite–Simpson
discretization scheme [15] and the Embedded COnic Solver (ECOS) [8].
In order to account for the additional degree of freedom provided by the
launcher, the problem is enhanced as follows:

  1. 1.

Three additional variables $v_{\infty}^{x},v_{\infty}^{y}$, and
$v_{\infty}^{z}$ are introduced and collected in the vector
$\mathbf{v}_{\infty}\in\mathbb{R}^{3}$.

  2. 2.

The initial boundary conditions on the velocity of the spacecraft in Eq. 1f
are expressed as

$\mathbf{v}_{i}=\mathbf{v}_{\text{E}}(t_{i})+\mathbf{v}_{\infty}$ (4)

where $\mathbf{v}_{\text{E}}(t_{i})$ indicates the velocity of the Earth at
the initial time $t_{i}$.

  3. 3.

In J2000, it is possible to express the constraints on the maximum magnitude
of the infinite velocity and on its maximum declination in Table 1 as

$\displaystyle\|\mathbf{v}_{\infty}\|_{2}\leq v_{\infty}^{\text{max}}$ (5a)
$\displaystyle-v_{\infty}^{z,\text{max}}\leq v_{\infty}^{z}\leq
v_{\infty}^{z,\text{max}}$ (5b) $\displaystyle
v_{\infty}^{z,\text{max}}=v_{\infty}^{\text{max}}\sin{\tilde{\delta}}$ (5c)

where $\tilde{\delta}$ is the maximum allowable declination. Note that the
constraint in Eq. 5a is a Second-Order Cone Constraint (SOCC). Therefore,
considering the constraints in Eqs. 5 still makes the problem in Eqs. 1
convex. It is worth to explicitly mention that the right ascension of the
infinite velocity is instead free to vary.

Although the SCP algorithm is rather robust [10], the solution of a given
problem depends on the initial guess. Therefore, to obtain an homogeneous
porkchop plot, a continuation method has been used, which considers the
solutions of already-solved problems as initial guess for their neighbours.
This also helps reducing significantly the computational time [19]. The
continuation procedure can be summarized as follows.

  1. 1.

First, the trajectory associated with the earliest departure date and longest
duration is computed with a simple shape-based initial guess. This trajectory
corresponds to the green box of the matrix in Fig. 1.

  2. 2.

The search space is discretized using time grids. In particular, a number $M$
of times of flight and a number $N$ of departure dates are considered.

  3. 3.

Then, for all $i=M-1,\dots,1$ and earliest departure date, the solution of the
(already solved) problem associated with higher time of flight is interpolated
to obtain the initial guess. Since an Hermite–Simpson discretization scheme is
used, the state in each trajectory segment is a third-order polynomial. If the
step in the time of flight direction is fairly small, the state of Apophis
correspondent to the cases $i-1$ and $i$ is similar, and therefore the
considered trajectory represents a consistently good initial guess. In the
case when the optimal control problem $(i-1,\ 1)$ did not converge, a cubic-
based initial guess is instead used for the problem $(i,\ 1)$.

  4. 4.

Once the trajectories associated with all the time of flights and first
departure date have been computed, an additional continuation is performed. In
particular, for all departure dates $j=2,\dots,N$, the solution of the problem
with time of flight $i$ and departure date $j-1$ is used from the new
departure date $j$ to the arrival date associated with the trajectory $(i,\
j-1)$ to obtain the initial guess for the problem $(i,\ j)$. In case the case
$(i,\ j-1)$ did not converge, a shape-based initial guess is used instead.

  5. 5.

For each of the departure dates, the FO trajectory with lower time of flight
that reached convergence is finally used as initial guess to solve the time-
optimal (TO) problem as well [15]. This continuation corresponds to the last
separated row of the matrix in Fig. 1.

In order to further speed-up the convergence process, the number of nodes used
by the algorithm to solve the problems is considered to be directly
proportional to the time of flight. In particular, a number of nodes
$P_{\textrm{max}}=150$ for the largest time of flight
$\text{ToF}_{\textrm{max}}=800$ days is chosen, and therefore

$P=\left\lceil
P_{\textrm{max}}\frac{\text{ToF}}{\text{ToF}_{\textrm{max}}}\right\rceil$ (6)

Figure 1: Continuation scheme for the feasibility plot adopted in the two
layers.

To properly detect the variations in propellant mass $m_{p}$ and having a more
efficient continuation, a variable step grid in the time-of-flight dimension
is exploited. From ToF of 244 days up to 460 days, the time step has been
considered to be 1 day, while a 20-day interval is considered from 460 on.
This step is necessary to refine the porkchop in the zone closer to the TO
solution, where both the convex and the indirect solvers benefit from denser
discretization grids. Contrarily, the time step for the departure dates is
fixed and set to 5 days. This grid requires about 16,000 trajectories to be
computed, confirming that a fast methodology is needed to assess the
feasibility of a low-thrust mission to Apophis.

###  3.2 Indirect Optimization Refinement

The results obtained with the convex step are used as inputs for the indirect
formulation and the generation of the final porkchop. Differently from direct
methods, indirect ones aim at satisfying the necessary optimality conditions
of the low-thrust trajectory optimization problem, which are derived through
the calculus of variations [11]. The methodology developed in previous works
[25] is used. It solves the shooting problem that stems from the imposition of
the necessary conditions as specified here after. In particular, the augmented
dynamical equations of states $[\mathbf{r},\mathbf{v},m]$ and costates
$[\bm{\lambda}_{r},\bm{\lambda}_{v},\lambda_{m}]$ for the two-body problem are
considered [22]

$\dot{\mathbf{y}}=\mathbf{F}(\mathbf{y})\Rightarrow\begin{cases}\dot{\mathbf{r}}=\mathbf{v}\\\
\dot{\mathbf{v}}=-\frac{\mu}{r^{3}}\mathbf{r}+u\frac{T_{max}}{m}\bm{\alpha}\\\
\dot{m}=-u\frac{T_{\text{max}}}{I_{\mathrm{sp}}g_{0}}\\\
\dot{\bm{\lambda}}_{r}=-\frac{3\mu}{r^{5}}(\mathbf{r}\cdot\bm{\lambda}_{v})\mathbf{r}+\frac{\mu}{r^{3}}\bm{\lambda}_{v}\\\
\dot{\bm{\lambda}}_{v}=-\bm{\lambda}_{r}\\\
\dot{\lambda}_{m}=-u\frac{T_{\mathrm{max}}}{m^{2}}\lambda_{v}\end{cases}$ (7)

where $u(t)\in[0,1]$ is the throttle factor and $\bm{\alpha}(t)$ is the thrust
direction vector, and have to be imposed to their optimal value
$\mathbf{u}^{*}(t)=[u^{*}(t),\bm{\alpha}^{*}(t)]$ such that the Pontryagin
Maximum Principle [3] is respected. To integrate these equations, the
algorithm enforces the initial conditions at $t_{i}$ for both the TO and FO
problems

$\mathbf{r}(t_{i})-\mathbf{r}_{i}=\mathbf{0},\quad\mathbf{v}(t_{i})-\mathbf{v}_{i}=\mathbf{0},\quad
m(t_{i})-m_{0}=0$ (8)

The initial mass $m_{0}$ is constant and considered equal to 500 kg as per the
assumptions in Table 1, whereas $\mathbf{r}_{i}$ is the initial position of
the Earth, which varies with the departure date considered in each point of
the search space mesh. The initial velocity is defined from the results of the
convex optimization step. In particular, according to Eq. (4), it is

$\mathbf{v}_{i}=\mathbf{v}_{\text{E}}(t_{i})+\mathbf{v}_{\infty}=\begin{bmatrix}v_{\text{E},x}(t_{i})+v_{\infty}\cos{{\alpha}}\cos{\delta}\\\
v_{\text{E},y}(t_{i})+v_{\infty}\sin{{\alpha}}\cos{\delta}\\\
v_{\text{E},z}(t_{i})+v_{\infty}\sin{\delta}\end{bmatrix}$ (9)

where $v_{\infty}$, $\alpha$, and $\delta$, come from the results of the first
convex layer.
To obtain the solutions in the search space mesh, the same continuation scheme
described in Fig. 1 is exploited for the indirect formulation. We first solve
the FO problem: $\mathbf{y}(t)=\bm{\varphi}(\mathbf{y}_{i},t_{i},t)$ being the
solution flow integrated from $t_{i}$ to a generic time instant $t$, the FO
shooting problems aim to find $\bm{\lambda}_{i}^{*}=\bm{\lambda}^{*}(t_{i})$
such that the solution at the final time $t_{f}$, namely
$\mathbf{y}(t_{f})=\bm{\varphi}([\mathbf{x}_{i},\,\bm{\lambda}_{i}^{*}],t_{i},t_{f})$,
satisfies the boundary conditions

$\begin{bmatrix}\mathbf{r}(t_{f})-\mathbf{r}_{f}\\\
\mathbf{v}(t_{f})-\mathbf{v}_{f}\\\
\lambda_{m}(t_{f})\end{bmatrix}=\mathbf{0}$ (10)

where $\mathbf{r}_{f}$ and $\mathbf{v}_{f}$ are the known final position and
velocity of Apophis, respectively, which depend on on the arrival date
considered for the specific point of the search space mesh. In the integration
of the flow, it can be proved [25] that the optimal thrust vector can be
expressed as $\bm{\alpha}^{*}=-\frac{\bm{\lambda_{v}}}{\lambda_{v}}$, whereas
the optimal throttle factor $u^{*}$ has the following structure:

$u^{*}=\begin{cases}0\quad&\mbox{if}\quad S_{FO}>0\\\
\in[0,1]\quad&\mbox{if}\quad S_{FO}=0\\\ 1\quad&\mbox{if}\quad
S_{FO}<0\end{cases}$ (11)

where $S_{FO}$ is the FO switching function,
$S_{FO}=1-\lambda_{m}-\lambda_{v}\frac{I_{sp}g_{0}}{m}$. When solving the TO
problem, the objective is to find the pair $[\lambda_{i}^{*},\,t_{f}^{*}]$
such that
$\mathbf{y}(t_{f}^{*})=\bm{\varphi}([\mathbf{x}_{i},\,\bm{\lambda}_{i}^{*}],t_{i},t_{f}^{*})$
obtained by the integration of Eq. (7) satisfies the boundary conditions

$\begin{bmatrix}\mathbf{r}(t_{f})-\mathbf{r}_{T}(t_{f})\\\
\mathbf{v}(t_{f})-\mathbf{v}_{T}(t_{f})\\\ \lambda_{m}(t_{f})\\\
H_{t}(t_{f})-\bm{\lambda}_{r}(t_{f})\cdot\mathbf{v}_{f}-\bm{\lambda}_{v}(t_{f})\cdot\mathbf{a}_{f}\end{bmatrix}=\mathbf{0}$
(12)

where $H_{t}$ is the Hamiltonian of the TO problem, i.e.,
$H_{t}=1+\bm{\lambda}^{T}\cdot\mathbf{f}$ with
$\mathbf{f}=\mathbf{f}(\mathbf{x}(t),\mathbf{u}(t))$ being the dynamics of the
spacecraft, and $\mathbf{a}_{f}=\dot{\mathbf{v}}_{f}$ is the acceleration of
Apophis at the final time. In the integration for the TO problem, the optimal
$\bm{\alpha}^{*}$ and $u^{*}$ follow the same logic of the FO problem, apart
from the switching function which is replaced with the TO switching function
$S_{TO}=-\lambda_{m}-\lambda_{v}\frac{I_{sp}g_{0}}{m}$. To solve the TO
problem the values of $v_{\infty}$, $\alpha$, and $\delta$ obtained from the
convex layer are exploited, as well as the value of the ToF, which is used as
initial guess for the optimal $t_{f}^{*}$.

##  4 Results

Considering the mission constraints in Table 1 and the approach in Section 3,
three different analyses are carried out. First, a reachability assessment
over the full departure and arrival windows is performed with a nominal
thrust-to-mass ratio. Three convenient launch windows are identified and
therefore three correspondent representative trajectories are analysed.
Finally, it is shown that a fourth launch window can open up if engines with
higher thrust-to-mass ratios are considered. Details about these analyses are
provided in the remainder.

###  4.1 Reachability assessment

Figure 3 shows the porkchop plot related to the reachability of the Apophis
asteroid under nominal conditions. This and the following plots are to be read
as follows:

  * •

The $x$-axis represents the departure date;

  * •

The $y$-axis represents the time of flight;

  * •

The color code represents the quantity of interest, indicated on the colorbar
on the right of the plot;

  * •

The patched red area indicates rendezvous happening after the latest arrival
accordingly to Table 1, and that are therefore unfeasible;

  * •

The oblique dashed segments represent rendezvous dates coinciding with, from
up to down respectively, -3, -4, -5, and -6 months from Apophis’ closest
encounter with the Earth;

  * •

The bold blue line indicates time-optimal solutions, i.e., the interplanetary
transfers for which the thruster is kept on for the whole time of flight. Note
that for a given departure date, solutions with lower times of flight are
unfeasible (corresponding to the blue-dashed area);

  * •

The dashed black line represents the available propellant of 73 kg.

The irregular behaviour at the center of the porkchop is related to the fly-by
at the Earth. In that point, a quasi-impulsive velocity change is experienced
by the asteroid in the heliocentric reference frame, making it difficult for
the optimizer finding a solution for the rendezvous. Figure 3 shows the
required $\Delta V$, computed using the Tsiolkovsky equation as $\Delta
V=I_{sp}g_{0}\ln{\frac{m_{0}}{m_{f}}}$.
Results show that, under the considered hypotheses, there are three feasible
launch windows, namely:

  * •

A first one, with launch dates between November and December 2026 and times of
flight that span from 540 and 800 days approximately (upper-left part of Fig.
3);

  * •

A second one, with launch dates between April and May 2027 and times of flight
of 620-660 days;

  * •

A third one, consistently larger, with launch dates from September to November
2027, and times of flight spanning from 250 to 500 days.

Figure 2: Propellant mass porkchop plot.

Figure 3: $\Delta V$ porkchop plot.

Figures 4 present the optimal magnitude, declination, and right ascension,
respectively, of the excess velocity of the spacecraft with respect to the
Earth in the equatorial, heliocentric reference frame (J2000) as a function of
the departure date and time of flight. In all the aforementioned windows, the
norm of the excess velocity reaches the maximum allowable value. Moreover, the
optimal declination of the excess velocity in the regions of interest never
assumes large values and it is included between +-40deg, values which are
compatible with the Ariane 6.2 and 6.4 launchers.

(a) Excess velocity magnitude $v_{\infty}$.

(b) Excess velocity declination $\delta$. (c) Excess velocity right ascension
$\alpha$.

Figure 4: Characteristics of the optimal excess velocity provided by the
launcher for each departure date and time of flight.

Figures 5 show the computational time required for the two layers of the
methodology. Results refer to the first departure date. A single optimal
solution requires at most 15 to run both the layers, with a mean runtime of
less than 9. Additionally, it is worth noting that the continuation scheme
proves beneficial in expediting the identification of solutions with lower
ToF. Typically, these solutions require longer computation times due to their
proximity to the time-optimal solution. As a matter of fact, the convex layer
requires less than 1 to compute the optimal solution. The indirect layer
requires slightly longer times due to the internal continuation from the
energy- to the fuel-optimal problem [25]. In conclusion, only few hours are
required to compute the porkchop plot and complete the reachability analysis.

(a) CPU time for the convex layer. (b) CPU time for the indirect layer.

Figure 5: CPU time required for all converged times of flight at the earliest
departure date.

###  4.2 Candidate Trajectories

From the reachability analysis presented in Fig. 3, three sample trajectories
have been selected, indicated with three dots:

  1. A.

(dark orange dot): Departure date on May 10, 2027;

  2. B.

(yellow dot): Departure date on November 11, 2027;

  3. C.

(white dot): Departure date on December 11, 2027.

The three trajectories have been selected because they are respectively
relevant due to  A. an early launch date that would allow backup solutions in
case of failure,  B. a convenient required propellant mass significantly lower
than the others, and  C. a later launch date that would still allow an early
arrival at Apophis with required propellant mass similar to the solution A.
No trajectory has been selected in the first available window (i.e., between
November and December 2026) as a launch in this period could be too early for
the completion of the mission design phases. The characteristics of the
selected trajectories have been reported in Table 2. Note that each of them
requires a propellant mass value lower than 73, and arrives two-to-three
months before the close encounter to ensure sufficient time for pre-fly-by
characterization of Apophis. The thrusting profiles of the selected
trajectories are shown in Fig. 6. In Figs. 7–9 some preliminary geometrical
considerations relative to the selected transfers are reported, as well as the
interplanetary trajectories in the heliocentric equatorial frame.

ID | Departure date | Arrival date |  ToF (d) |  $\mathbf{m}_{p}$ (kg) |  $\bm{\Delta}\mathbf{V}$ (km/s)
---|---|---|---|---|---
A | May 10, 2027 | Feb 8, 2029 | 640 | 70.43 | 2.21
B | Nov 11, 2027 | Jan 23, 2029 | 439 | 54.95 | 1.69
C | Dec 11, 2027 | Jan 14, 2029 | 400 | 72.23 | 2.29
Table 2: Characteristics of the selected representative trajectories.

As expected, solution A exhibits a rather empty thrust profile, while
trajectory C requires longer and more frequent thrusting arcs. This behaviour
correlates with the distance from the time-optimal solution. On the other
hand, there is no significant differences in the distance and in the geometric
quantities. This finding can simplify the system design, since communications
with the ground, power generation, and relative navigation with the target are
not severely affected by the trajectory of choice.

(a) Thrust profile of solution A. (b) Thrust profile of solution B. (c) Thrust
profile of solution C.

Figure 6: Thrust profiles of trajectories A, B, and C.

(a) Spacecraft-Apophis distance. (b) Spacecraft-Earth distance.

(c) Spacecraft-Sun distance. (d) Target phase angle.

(e) Trajectory representation.

Figure 7: Geometrical analysis and representation in the J2000 reference frame
of the solution A.

(a) Spacecraft-Apophis distance. (b) Spacecraft-Earth distance.

(c) Spacecraft-Sun distance. (d) Target phase angle.

(e) Trajectory representation.

Figure 8: Geometrical analysis and representation in the J2000 reference frame
of the solution B.

(a) Spacecraft-Apophis distance. (b) Spacecraft-Earth distance.

(c) Spacecraft-Sun distance. (d) Target phase angle.

(e) Trajectory representation.

Figure 9: Geometrical analysis and representation in the J2000 reference frame
of the solution C.

###  4.3 Analyses with Higher Thrust-to-Mass Ratios

The nominal analysis has been conducted considering a $T_{\text{max}}/m_{0}$
of 1.2 e-4, corresponding to a thrust of 60. In order to explore if additional
opportunities are present when different propulsive systems are considered, a
sensitivity analysis is performed, considering higher thrust-to-mass ratios,
on the bottom right portion of Fig. 3. A possible window in this region will
allow a late launch to the asteroid, giving more time for the design and
integration of the spacecraft.
Figs. 10 show that an additional window that satisfies the time and propellant
constraints opens up in the timespan March–May 2028 with 270-360 days of time
of flight, when thrusters with slightly higher thrust-to-mass ratios are
considered (namely, 1.8 e-4 and 2 e-4 corresponding to, respectively, thrust
levels of 90 and 100 and the nominal initial mass of 500). For these plots,
the patched red area indicates rendezvous happening later than one month
before the close encounter of the asteroid with the Earth. In conclusion,
higher thrust levels may open up new feasible launch windows and thus improve
the asteroid reachability.

(a) $T_{\text{max}}/m_{0}=1.8\times 10^{-4}m/s^{2}$. (b)
$T_{\text{max}}/m_{0}=2\times 10^{-4}m/s^{2}$.

Figure 10: Porkchop plots for higher thrust-to-mass ratios.

##  5 Conclusion

In the context of the RAMSES mission, a fast and robust methodology has been
introduced to evaluate the reachability of Apophis with a low-thrust satellite
before its close encounter with the Earth on April 13, 2029. To this aim,
porckchop plots from the Earth to Apophis have been computed using a two-layer
approach, exploiting a direct sequential convex programming algorithm followed
by an indirect method. This method allows an easy handling of the free
launching condition, while guaranteeing optimality of the solution, and it is
able to compute more than 16,000 trajectories in few hours. Under the
considered hypotheses, three feasible launching windows are identified. Each
window lasts about 3 months and is separated from the other by about half a
year, with the last one opening in September 2027. A fourth late short window,
between April and May 2029, can be found if higher thrust-to-mass ratios are
considered.

Acknowledgments

This work has been performed in response to ESA call AO/1-11659/23/NL/GLC:
Pre-phase A/Phase A System Study of an Apophis Mission (RAMSES Definition
Phase).

Conflict of interest The authors declare that they have no known competing
financial interests or personal relationships that could have appeared to
influence the work reported in this paper.

## References

  *   * Binzel et al [2021] Binzel R, Barbee B, Barnouin O, et al (2021) Apophis 2029: Decadal opportunity for the science of planetary defense. Bulletin of the AAS 53(4). 10.3847/25c2cfeb.f87e0599
  * Bonalli et al [2019] Bonalli R, Cauligi A, Bylard A, et al (2019) Gusto: Guaranteed sequential trajectory optimization via sequential convex programming. In: 2019 International Conference on Robotics and Automation (ICRA), pp 6741–6747, 10.1109/ICRA.2019.8794205
  * Bryson and Ho [1975] Bryson A, Ho YC (1975) Applied optimal control. Taylor & Francis, London, 10.1201/9781315137667
  * Chesley [2005] Chesley SR (2005) Potential impact detection for near-Earth asteroids: The case of 99942 Apophis (2004 MN4). Proceedings of the International Astronomical Union 1(S229):215–228. 10.1017/S1743921305006769
  * Dannenmayer and Mazouffre [2009] Dannenmayer K, Mazouffre S (2009) Elementary scaling laws for the design of low and high power Hall effect thrusters. In: Progress in Propulsion Physics, EDP Sciences, pp 601–616
  * DellaGiustina et al [2022] DellaGiustina D, Golish D, Guzewich S, et al (2022) Osiris-APEX: A proposed OSIRIS-REx extended mission to Apophis. LPI Contributions 2681:2011
  * DeMartini et al [2019] DeMartini JV, Richardson DC, Barnouin OS, et al (2019) Using a discrete element method to investigate seismic response and spin change of 99942 Apophis during its 2029 tidal encounter with Earth. Icarus 328:93–103. 10.1016/j.icarus.2019.03.015
  * Domahidi et al [2013] Domahidi A, Chu E, Boyd S (2013) ECOS: An SOCP Solver for Embedded Systems. In: European Control Conference, Zurich, Switzerland, pp 3071–3076, doi: 10.23919/ECC.2013.6669541
  * Giorgini et al [2008] Giorgini JD, Benner LA, Ostro SJ, et al (2008) Predicting the Earth encounters of (99942) Apophis. Icarus 193(1):1–19. 10.1016/j.icarus.2007.09.012
  * Hofmann et al [2023] Hofmann C, Morelli AC, Topputo F (2023) Performance assessment of convex low-thrust trajectory optimization methods. Journal of Spacecraft and Rockets 60(1). 10.2514/1.A35461
  * Kechichian [1997] Kechichian JA (1997) Optimal low-Earth-orbit-geostationary-Earth-orbit intermediate acceleration orbit transfer. Journal of Guidance, Control, and Dynamics 20(4):803–811. 10.2514/2.4116
  * Kim et al [2023] Kim Y, DeMartini JV, Richardson DC, et al (2023) Tidal resurfacing model for (99942) Apophis during the 2029 close approach with Earth. Monthly Notices of the Royal Astronomical Society 520(3):3405–3415. 10.1093/mnras/stad351
  * Machuca et al [2020] Machuca P, Sanchez JP, Masdemont JJ, et al (2020) High-fidelity trajectory design to flyby near-Earth asteroids using CubeSats. Acta Astronautica 167:146–163. 10.1016/j.actaastro.2019.09.041
  * Malyuta et al [2022] Malyuta D, Reynolds TP, Szmuk M, et al (2022) Convex optimization for trajectory generation: A tutorial on generating dynamically feasible trajectories reliably and efficiently. IEEE Control Systems Magazine 42(5):40–113. 10.1109/MCS.2022.3187542
  * Morelli et al [2021] Morelli AC, Hofmann C, Topputo F (2021) Robust low-thrust trajectory optimization using convex programming and a homotopic approach. IEEE Transactions on Aerospace and Electronic Systems 58(3):2103–2116. 10.1109/TAES.2021.3128869
  * Scheeres et al [2005] Scheeres D, Benner L, Ostro S, et al (2005) Abrupt alteration of Asteroid 2004 MN4’s spin state during its 2029 Earth flyby. Icarus 178(1):281–283. 10.1016/j.icarus.2005.06.002
  * Sokolov et al [2012] Sokolov L, Bashakov A, Borisova T, et al (2012) Impact trajectories of the asteroid apophis in the 21st century. Solar System Research 46:291–300. 10.1134/S0038094612040077
  * Souchay et al [2014] Souchay J, Souami D, Lhotka C, et al (2014) Rotational changes of the asteroid 99942 Apophis during the 2029 close encounter with Earth. Astronomy & Astrophysics 563:A24. 10.1051/0004-6361/201322364
  * Topputo et al [2021] Topputo F, Wang Y, Giordano C, et al (2021) Envelop of reachable asteroids by M-ARGO CubeSat. Advances in Space Research 67(12):4193–4221. 10.1016/j.asr.2021.02.031
  * Valvano et al [2022] Valvano G, Winter OC, Sfair R, et al (2022) Apophis–effects of the 2029 Earth’s encounter on the surface and nearby dynamics. Monthly Notices of the Royal Astronomical Society 510(1):95–109. 10.1093/mnras/stab3299
  * Wagner et al [2015] Wagner S, Wie B, Barbee BW (2015) Target selection for a hypervelocity asteroid intercept vehicle flight validation mission. Acta Astronautica 107:247–261. 10.1016/j.actaastro.2014.11.037
  * Wang and Topputo [2022] Wang Y, Topputo F (2022) Indirect optimization of power-limited asteroid rendezvous trajectories. Journal of Guidance, Control, and Dynamics 45(5):962–971. 10.2514/1.G006179
  * Wang and Grant [2018] Wang Z, Grant MJ (2018) Minimum-Fuel Low-Thrust Transfers for Spacecraft: A Convex Approach. IEEE Transactions on Aerospace and Electronic Systems 54(5):2274–2290. 10.1109/TAES.2018.2812558
  * Yu et al [2014] Yu Y, Richardson DC, Michel P, et al (2014) Numerical predictions of surface effects during the 2029 close approach of asteroid 99942 Apophis. Icarus 242:82–96. 10.1016/j.icarus.2014.07.027
  * Zhang et al [2015] Zhang C, Topputo F, Bernelli-Zazzera F, et al (2015) Low-thrust minimum-fuel optimization in the circular restricted three-body problem. Journal of Guidance, Control, and Dynamics 38(8):1501–1510. 10.2514/1.G001080
  * Zhang and Michel [2020] Zhang Y, Michel P (2020) Tidal distortion and disruption of rubble-pile bodies revisited-soft-sphere discrete element analyses. Astronomy & Astrophysics 640:A102. 10.1051/0004-6361/202037856