\section{Introduction}
Asymmetric dark matter (ADM) scenario sheds light on 
the coincidence problem between the observed baryon and dark matter (DM) abundances in the universe
~\cite{Nussinov:1985xr,Barr:1990ca,Barr:1991qn,Dodelson:1991iv,Kaplan:1991ah,Kuzmin:1996he,Foot:2003jt,Foot:2004pq,Kitano:2004sv,Gudnason:2006ug,Kaplan:2009ag} 
(see also~\cite{Davoudiasl:2012uw,Petraki:2013wwa,Zurek:2013wia} for reviews).
If the DM abundance is provided by a mechanism which is unrelated to the baryogenesis, it is quite puzzling why those 
abundances are close with each other despite the fact that the baryon abundance is dominated by 
the contribution from the matter-antimatter asymmetry. 
In the ADM scenario the coincidence problem can be explained when the DM mass is of ${\cal O}(1)$\,GeV,
where the matter-antimatter asymmetry is thermally distributed between the dark and the Standard Model (visible) sectors.

Among various ADM scenarios, composite baryonic DM in QCD-like dynamics is particularly motivated
since it can naturally provide a large annihilation cross section and  the DM mass in the GeV range 
simultaneously~\cite{Foot:2003jt,Foot:2004pq,Berezhiani:2005ek,Alves:2009nf,An:2009vq,Alves:2010dd,Gu:2012fg,
Buckley:2012ky,Detmold:2014qqa,Gu:2014nga,Lonsdale:2018xwd}.
Recently, a minimal composite ADM model and its ultraviolet (UV) completion~\cite{Ibe:2018juk,Ibe:2018tex,Ibe:2019ena} have been proposed where the asymmetry 
generated by the thermal leptogenesis~\cite{Fukugita:1986hr} 
(see also \cite{Giudice:2003jh,Buchmuller:2005eh,Davidson:2008bu} for review) 
is thermally distributed between the two sectors through 
a portal operator associated with the seesaw 
mechanism~\cite{Minkowski:1977sc, Yanagida:1979as, GellMann:1980vs, Glashow:1979nm, Mohapatra:1979ia}.
The dark sector of the model consists of QCD-like dynamics and QED-like interaction, 
which are called as dark QCD and dark QED, respectively.
The lightest baryons of dark QCD play the role of ADM.
The dark QED photon (dark photon) obtains a mass of $\order{10\mbox{-}100}$ MeV,
which plays a crucial role to transfer the excessive entropy of the dark sector
into the visible sector before neutrino decoupling~\cite{Ibe:2018juk,Blennow:2012de}.

In this paper, we discuss the indirect detection of the composite ADM model in \cite{Ibe:2018juk,Ibe:2018tex,Ibe:2019ena}.
The portal operator in this model is generated in association with the seesaw mechanism.
In this model, the dark-neutron, one of the lightest dark baryons, inevitably obtains a tiny Majorana mass. 
Such a tiny Majorana mass induces the oscillation between DM particle and the antiparticle,
which induces a pair-annihilation of ADM at late times~\cite{Cai:2009ia,Buckley:2011ye,Cirelli:2011ac,Tulin:2012re,Okada:2012rm,Hardy:2014dea,Chen:2015yuz}.
A pair of DM particle and the antiparticle annihilates into multiple dark pions, 
and the (neutral) dark pion subsequently decays into a pair of the dark photons.
The dark photon eventually decays into an electron-positron pair. 
Thus, the late time annihilation of ADM results in multiple 
soft electrons/positrons.
In addition, soft photons are also emitted as final state radiation.
As we will see, the model can be efficiently tested 
by the searches for the $\gamma$-ray from the dwarf spheroidal galaxies (dSphs) by the Fermi-LAT.
We also discuss the constraints from 
the observations of the interstellar electron/positron flux by the Voyager-1.

The organization of the paper is as follows.
In section~\ref{sec:model}, we review the composite ADM model in \cite{Ibe:2018juk,Ibe:2018tex,Ibe:2019ena} and show 
how the tiny Majorana mass of the dark neutron appears associated with the seesaw mechanism.
In section~\ref{sec:signal}, we derive the expected $\gamma$-ray flux from the dSphs and discuss the constraints 
on the model by comparing the flux with the Fermi-LAT results. 
We also estimate the interstellar electron/positron flux in cosmic ray from the late time annihilation and compare it with the Voyger-1 result.
The final section is devoted to the conclusions. 

\section{DM anti-DM oscillation in the composite ADM model}\label{sec:model}
\subsection{A Model of Composite ADM}
In this subsection, we briefly review the composite ADM model in~\cite{Ibe:2018juk,Ibe:2018tex,Ibe:2019ena}.
The model is based on $N_g$-generation dark quarks with $SU(3)_D \times U(1)_D$ gauge symmetry.
$SU(3)_D$ provides the dark QCD dynamics and $U(1)_D$ the dark QED interaction.
The dark quarks are the fundamental representations of $SU(3)_D$.
They are charged under the dark QED and the $B-L$ in analogy to
the up-type and the down-type quarks in the visible sector (see Tab.~\ref{tab:Dmatter}).
They have tiny masses, 
\begin{align}
\label{eq:mass}
{\cal L}_\mathrm{mass} = m_{U}' \overline{U}' U' + m_{D}'\overline{D}' D' + \mathrm{h.c.}\ ,
\end{align}
with $m_{U}'$ and $m_{D}'$ being the mass parameters.
Hereafter, we put primes on the parameters and the fields in the dark sector
when there  are counterparts in the visible sector.

The dark QCD exhibits confinement below the dynamical scale of $SU(3)_D$, 
$\Lambda_{\rm QCD}'$, which leads to the emergence of the dark baryons and the dark mesons.
Throughout this paper, we assume that only one generation of the dark quarks have masses 
smaller than $\Lambda_{\rm QCD}'$.%
\footnote{For $N_g >1$, we assume the heavier dark quarks decay into the lighter 
ones by emitting the dark Higgs boson which has the dark QED charge of $1$.
It should be noted that the dark quark masses are not generated by the vacuum expectation 
value of the dark Higgs boson~\cite{Ibe:2018juk,Ibe:2018tex,Ibe:2019ena},
and hence, the dark Higgs couplings generically violate the flavor symmetry in the dark sector.}
The lightest dark baryons, i.e. the dark nucleons, 
\begin{align}
p' \propto U'U'D'\ , \quad
\bar{p}' \propto \bar{U}'\bar{U}'\bar{D}'\ , \quad
n' \propto U'D'D'\ , \quad
\bar{n}' \propto \bar{U}'\bar{D}'\bar{D}'\ , 
\end{align}
are stable in the decoupling limit from the visible sector due to their $B-L$ charges.
Once the $B-L$ asymmetry is shared between the visible and the dark sector, 
the dark nucleon abundance is dominated by the asymmetric component due to their large annihilation cross section.
Therefore, the dark nucleon with a mass in the GeV range is a good candidate for ADM.



\begin{table}
	\centering
	\caption{The charge assignment of dark quarks. We assume $N_g$ generations of the dark quarks,
although only one generation has a mass smaller than $\Lambda_{\rm QCD}'$.
The $U(1)_{B-L}$ symmetry is the global symmetry which is shared with the visible sector.
	}
	\label{tab:Dmatter}
	\begin{tabular}{|c|c|c||c|}
		\hline
		& $SU(3)_D$ & $U(1)_D$ & $U(1)_{B-L}$ \\ \hline
		$U'$ & $\mathbf{3}$ & $2/3$ & $1/3$ \\
		$\overline U'$ & $\overline{\mathbf{3}}$ & $-2/3$ & $- 1/3$ \\
		$D'$ & $\mathbf{3}$ & $-1/3$ & $1/3$ \\
		$\overline D'$ & $\overline{\mathbf{3}}$ & $1/3$ & $- 1/3$ \\ \hline
	\end{tabular}
\end{table}

When the $B-L$ asymmetry is thermally distributed between the visible and the dark sectors,
the ratio of the $B-L$ asymmetry stored in each sector is given by $A_{\rm DM}/A_{\rm SM} = 44 N_g/237$
for the $B-L$ charges given in Tab.~\ref{tab:Dmatter}~\cite{Fukuda:2014xqa}.%
\footnote{In the presence of additional $B-L$ charged fields in the dark sector, such as dark leptons,
the ratio can be modified.
Besides, the neutrality condition of $U(1)_D$ and the contributions from the dark Higgs sector
also change the ratio by some tens percent for a given $N_g$.}
Thus, the observed ratio of the DM and the baryon abundance can be reproduced when the dark nucleon mass is
\begin{align}
m_{N}' \simeq \frac{\Omega_{\rm DM}}{\Omega_{\rm B}} \frac{A_{\rm B}}{A_{\rm SM}}\frac{A_{\rm SM}}{A_{\rm DM}}
\times m_N
\simeq \frac{8.5\,\rm GeV}{N_g}\ .
\end{align}
Here, we have used the ratio of the baryon asymmetry to the $B-L$ asymmetry in the visible sector, $A_{\rm B}/A_{\rm SM} = 30/97$~\cite{Harvey:1990qw}.
The dark nucleon mass in this range can be naturally realized when $\Lambda_{\rm QCD}'$ is in the GeV range.





The lightest dark mesons,
\begin{eqnarray}
\pi'^0\propto U' \bar U' - D' \bar D' \ , \quad \pi'^+ \propto U'\bar D'
\ , \quad \pi'^- \propto D'\bar U'\ ,
\end{eqnarray}
annihilate or decay into the dark photons.
As a result, they do not contribute to the effective number of neutrino degrees of freedom
nor to the dark matter abundance significantly even if they are stable.
In the following analysis, we assume that the dark charged pions are stable for simplicity.%
\footnote{If $U(1)_D$ is broken by the vacuum expectation value of  a dark Higgs with the dark QED charge of $2$, 
a ${\mathbb Z}_2$ symmetry remains unbroken which makes the dark charged pion stable.
If $U(1)_D$ is broken by the dark Higgs with the charge $1$, the neutral and the charged pions can mix each other, and hence, the charged pions decay.}
The decay of the dark neutral pion into a pair of dark photons, on the other hand, is inevitable due to the chiral anomaly.
As we will see, the decay of the neutral pion plays a central role for the indirect detection of ADM.

The dark photon obtains its mass by the dark Higgs mechanism, 
and it decays into the visible fermions thought the kinetic mixing with the visible QED photon, 
\begin{align}
{\cal L}_{\gamma^\prime} = \frac\epsilon2 F_{\mu\nu} F^{\prime\mu\nu} + \frac{m_{\gamma^\prime}^2}{2} A_\mu^\prime A^{\prime \mu} \,.
\end{align}
Here, $F_{\mu\nu}$ and $F^\prime_{\mu\nu}$ denote the field strengths of 
the visible and the dark QED with $A^\prime_\mu$ being the dark photon gauge field. 
In the following, we assume the kinetic mixing parameters of $\epsilon = 10^{-10}$--$10^{-8}$
and the dark photon mass in $\order{10\mbox{-}100}$ MeV range which 
satisfies all the constraints~\cite{Ibe:2018juk} (see also \cite{Bauer:2018onh,Chang:2016ntp, Chang:2018rso}).%
\footnote{See~\cite{Ibe:2018juk} for discussion on the origin of the tiny kinetic mixing parameters.}
In this parameter range, the dark photon decays when the cosmic temperature is above ${\cal O}(1)$\,MeV.

Finally, let us comment on the ratio between the abundances of the dark protons and the dark neutrons.
In the present model, there is no dark leptons nor dark weak gauge bosons.
Besides, it is expected that 
the mass difference 
between the dark neutron and the dark proton is 
smaller than the mass of the dark pion
when the dark quark masses are smaller than the dynamical
scale of $SU(3)_D$.
Thus, the dark neutron is stable in the limit of the 
vanishing $B-L$ portal interactions (see below).
The ratio between the dark proton abundance 
and the dark neutron abundance is given by~\cite{Ibe:2018tex},
\begin{align}
n_{n}'/n_{p}' \sim e^{-(m_{n}'-m_p')/T_F}\ .
\end{align}
Here, $n_{n,p}'$ and $m_{n,p}'$ are the number densities 
and the masses of the dark neutron and the dark proton, respectively.
$T_F$ denotes the freeze-out temperature of the dark pion annihilation,
$T_F \simeq m_\pi'/{\cal O}(10)$.
Thus, for $m_n'-m_p' \ll m_\pi'$, the dark neutron abundance is comparable 
to that of the dark proton.
In the following, we take $n_n' = n_p'$.


\subsection{The \texorpdfstring{$B-L$}{} portal operator}
The $B-L$ asymmetry generated by thermal leptogenesis is thermally distributed between the visible
and the dark sectors.
For this purpose, there need to be portal interactions which connect the $B-L$ symmetry in the two sectors.
In the model in \cite{Ibe:2018juk} (see also \cite{Ibe:2011hq,Fukuda:2014xqa}), the following operators are assumed 
as the portal operators,
\begin{align}
{\cal L}_{\text{portal}} \sim
\frac{1}{M_*^3} (\overline U' \overline D' \overline D') (L H)
+ \frac{1}{M_*^3} (U'^\dag D'^\dag \overline D') (L H) + \mathrm{h.c.} \,,
\label{eq:portal}
 \end{align}
where $L$ and $H$ are the lepton and the Higgs doublets in the visible sector, and $M_*$ is a dimensional parameter.%
\footnote{The portal operators require the gauge invariant operators which are charged under the $B-L$ symmetry. 
This is the reason why we need both the up-type and the down-type quarks in the dark sector.}
Here, we omit the $\order{1}$ coefficients.
The effects of the above operators decouple at the cosmic temperature below 
$T_{*} \sim M_{*} (M_{*} / M_{\rm PL})^{1/5}$.
Here, $M_{\rm PL} = 2.4 \times 10^{18}~\mathrm{GeV}$ denotes the reduced Planck mass.
For successful ADM with thermal leptogenesis, the decoupling 
temperature, $T_*$, is required to be lower than the temperature, $T_{B-L}$, at which leptogenesis completes.
In the following, we consider the so-called strong washout regime of thermal leptogenesis, 
where the leptogenesis completes at the temperature about $T_{B-L} \simeq M_{R}/z_{B-L}$ with
$z_{B-L} \simeq 10$~\cite{Buchmuller:2005eh}.

In \cite{Ibe:2018juk,Ibe:2018tex}, the UV model has been proposed 
in which the portal operators in Eq.\,\eqref{eq:portal} 
are generated by integrating out the right-handed neutrinos, $\bar{N}$, and the dark colored Higgs boson, $H_C'$.
The gauge charges of $H_C'$ are identical to those of $D'$, while $H_C'$
has the $B-L$ charge $-2/3$.
The right-handed neutrinos couple to both sectors via,
\begin{align}
{\cal L} = \frac{M_R}{2} 
\bar N \bar N + y_N L H \bar N +\frac{1}{2}M_C^2 |H_C'|^2-Y_N H_C'\overline D'\bar{N} 
-Y_C H_C' U' D' - Y_{\bar{C}} H_C'^\dagger \overline U' \overline D'
+ \mathrm{h.c.} 
\label{eq:UVportal}
\end{align}
Here, $M_C$ denotes the dark colored Higgs mass, $M_R$ the mass of the right-handed neutrinos, 
and $y_N$ and $Y$'s are the Yukawa coupling constants.
The flavor and the gauge indices are suppressed.
It should be noted that the mass terms of the right-handed neutrino 
break $B-L$ symmetry explicitly.
The first two terms are relevant for the seesaw mechanism.


By integrating out $\bar{N}$ and $H_C'$ from Eq.\,\eqref{eq:UVportal}, the portal operators in Eq.\,\eqref{eq:portal} are obtained 
where $M_*$ corresponds to 
\begin{align}
    \frac{1}{M_*^3} = \frac{y_NY_NY_{\bar C}}{2M_C^2 M_R}\ , \quad 
    \frac{1}{M_*^3} = \frac{y_NY_NY_{C}^*}{2M_C^2 M_R}\ ,
    \label{eq:mstar}
\end{align}
for each term of Eq.\,\eqref{eq:portal}, respectively.
From the condition of $T_* < T_{B-L}$, the mass of the dark colored Higgs should satisfy,%
\footnote{Hereafter, we take $M_R >0$ and neglect the complex phases of $Y_N$, $Y_C$ and $Y_{\bar C}$.
We also assume $Y_C = Y_{\bar C}$ for simplicity.}
\begin{align}
\label{eq: upper bound on M_C}
    \frac{M_R}{z_{B-L}}\lesssim M_C\lesssim \frac{10}{z_{B-L}^{5/4}}\left(\frac{\hat{m}_\nu}{0.1\rm eV}\right)^{1/4}\sqrt{Y_N Y_C} M_R\ .
\end{align}
The first inequality comes from a consistency condition of the decoupling limit of 
the dark colored Higgs
at the temperature $T_{B-L}$.
In the right hand side, we have reparameterized the neutrino Yukawa coupling by using a tiny neutrino mass parameter, 
$\hat{m}_\nu$,
\begin{align}
|y_{N}^{2}| \sim 10^{-5} \left( \frac{\hat{m}_{\nu}}{0.1 \, \mathrm{eV}} \right) \left( \frac{M_{R}}{10^{9} \, \mathrm{GeV}} \right) \ .
\label{eq:seesaw}
\end{align}


Incidentally, 
the dark nucleon can decay into the dark pion and the anti-neutrino in the visible sector through the $B-L$ portal operator in Eq.\,\eqref{eq:portal}~\cite{Fukuda:2014xqa}.
The lifetime is roughly given by,
\begin{align}
\label{eq:lifetime}
    \tau_N' \simeq 10^{33}\,{\rm sec} 
     \left(\frac{2\,{\rm GeV}}{\Lambda_{\rm QCD}'}\right)^4
      \left(\frac{0.1\,{\rm eV}}{\hat{m}_\nu}\right)
      \left(\frac{M_R}{10^9\,{\rm GeV}}\right)
    \left(\frac{\tilde{M}_C}{3\times 10^9\,{\rm GeV}}\right)^4
        \left(\frac{10\,{\rm GeV}}{m_{\rm DM}}\right)\ ,   
          \end{align}
where $\tilde{M}_C = M_C/\sqrt{Y_NY_C}$.
Thus, the lifetime of the dark nucleons is 
much longer than the age of the universe 
for $M_R \sim M_C \sim 10^9$\,GeV.

\subsection{The Majorana mass of the dark neutron}
The portal operators in Eq.\,\eqref{eq:portal} are generated in association with the seesaw mechanism.
As a notable feature of the UV completion model in Eq.\,\eqref{eq:UVportal}, it 
also leads to the Majorana mass term of the dark neutron.
This can be observed by integrating out $H_C'$ and $\bar{N}$ one by one.
In the case of $M_C>M_R$, we first integrate out $H_C'$ from Eq.\,\eqref{eq:UVportal}, which reads
	\begin{align}
		\label{eq: integrate out Hc'}
		{\cal L} =\,& \bar{N}^\dagger i\sigma^\mu\partial_\mu \bar{N} + 
		\frac{M_R}{2} \bar{N} \bar{N} + y_N L H \bar{N}\notag \\
		&-\frac{\abs{Y_N}^2}{2M_C^2} \bar{D}'\bar{N}\qty(\bar{D}'\bar{N})^\dagger -\frac{Y_N}{2M_C^2} \bar{D}'\bar{N}\qty[ \qty(Y_CU'D')^\dagger + Y_{\bar{C}}\bar{U}'\bar{D}']
		+ \mathrm{h.c.}\notag\\
		&+(\text{quartic in dark quark fields}).
	\end{align}
Here, we show the kinetic term of $\bar{N}$ explicitly which were implicit in Eq.\,\eqref{eq:UVportal}.
This formula is of the form
	\begin{align}
	\label{eq: simplyfied N integrate out}
		{\cal L}=\qty(A\bar{N}\bar{N} + B\bar{N}+\mathrm{h.c.})-C\bar{N}\bar{N}^\dagger,
	\end{align}
where%
\footnote{Here, $\chi^\dagger\sigma^\mu\overset{\leftrightarrow}{\partial}_\mu\eta 
=\chi^\dagger\sigma^\mu\partial_\mu\eta
-\partial_\mu\chi^\dagger\sigma^\mu\eta$
for the Weyl fermions, $\chi$  and $\eta$.
}
	\begin{align}
	\label{eq: formula for ABC}
		A=\frac{M_R}{2},\quad B=y_NLH-\frac{Y_N}{2M_C^2} \bar{D}'\qty[ Y_{\bar{C}}\bar{U}'\bar{D}' + \qty(Y_CU'D')^\dagger],\quad C=-\frac{i}{2}\sigma^\mu\overset{\leftrightarrow}{\partial}_\mu+\frac{\abs{Y_N}^2}{2M_C^2} \bar{D}'\bar{D}'^\dagger.
	\end{align}
To make $\bar N$ integrated out, it is convenient to complete the square of Eq.\,(\ref{eq: formula for ABC}) with respect to $\bar{N}$.
For this purpose, we shift $\bar{N}$ by $\bar{N}\to\bar{N}+\bar{\psi}$, with which we can eliminate the linear term 
in Eq.\,\eqref{eq: simplyfied N integrate out}.
The condition $\bar{\psi}$ must satisfy is $2A\bar{\psi}+B-C\bar{\psi}^\dagger=0$, which reads
	\begin{align}
		\bar{\psi}=-\frac{2A^*B+CB^\dagger}{4\abs{A}^2-C^2}\simeq -\frac{1}{2\abs{M_R}^2}\qty(1+\frac{C^2}{\abs{M_R}^2})\qty(M_R^*B+CB^\dagger )
	\end{align}
After the shift, we integrate out $\bar{N}$ to obtain
	\begin{align}
	\label{eq: fully integrated out portal}
		{\cal L}&=\qty(A\bar{\psi}\bar{\psi} + B\bar{\psi}+\mathrm{h.c.})-C\bar{\psi}\bar{\psi}^\dagger\notag \\
					%
					&=\frac{1}{2}B\bar{\psi}+\mathrm{h.c.}\notag \\
					%
					&=-\qty(1+\frac{C^2}{\abs{M_R}^2})\qty(\frac{1}{2M_R}BB+\frac{1}{2\abs{M_R}^2}CB^\dagger B)+\mathrm{h.c.}
	\end{align}
From Eq.\,\eqref{eq: formula for ABC}, we find that $BB$ term includes the Mojorana mass term of the dark neutron
	\begin{align}
	\label{eq: Majorana mass term}
		\frac{1}{2M_R} BB \supset \frac{Y_N^2Y_{\bar{C}}^2}{8M_R M_C^4}\qty(\bar{U}'\bar{D}'\bar{D}')^2\sim
		\frac{Y_N^2Y_{\bar{C}}^2\Lambda_{\rm QCD}'^6}{8M_R M_C^4}\bar{n}'\bar{n}'.
	\end{align}
In this way, Eq.\,\eqref{eq:UVportal} leads to the Majorana mass,
\begin{align}
\label{eq:MajoranaMass}
    m_M = \frac{Y_N^2Y_{\bar{C}}^2\Lambda_{\rm QCD}'^6}{4M_R M_C^4}= \frac{\Lambda_{\rm QCD}'^6}{4M_R \tilde{M}_C^4}\ ,
\end{align}
in addition to the $B-L$ portal  operators in Eq.\,\eqref{eq:portal}.

Once the dark neutron obtains the Majorana mass, the dark neutron and the anti-dark neutron
oscillate with a time scale of $t_{\rm osc} = m_M^{-1}$~\cite{Cai:2009ia,Buckley:2011ye,Cirelli:2011ac,Tulin:2012re,Okada:2012rm,Hardy:2014dea}.
The probability to find an anti-dark neutron at a time $t$
is given by,
\begin{align}
    P(n'\leftrightarrow \bar{n}') = \sin^2(m_Mt) \ .
\end{align}
Here, we assume that the initial state at $t = 0$ is a pure dark neutron state.
As we will see in the next section, the oscillation 
induces a pair-annihilation of ADM which ends up with multiple 
soft electrons/positrons/photons.

\subsection{Washout Interactions and On-Shell Portal}
Before closing this section, let us discuss the $B-L$ washout interactions which 
are also induced from Eq.\,\eqref{eq:UVportal}.
In fact, the term $CB^\dagger B$  in 
Eq.\,\eqref{eq: fully integrated out portal}
includes
	\begin{align}
	\label{eq:washout}
		\frac{1}{2\abs{M_R}^2} CB^\dagger B\supset \frac{y_NY_N^*}{2M_C^2 \abs{M_R}^2}
																				\qty[Y_{\bar C}^*(\overline U'^\dagger \overline D'^\dagger \overline D'^\dagger)
																			(i\sigma^\mu\partial_\mu)(L H)
																				%
																				+Y_C(U' D' \overline D'^\dagger)(i\sigma^\mu\partial_\mu)(L H)]\ .
	\end{align}
	In these interaction terms, and those in Eq.\,\eqref{eq:portal}, $L$ couples to the dark sector operators which have the opposite $B-L$ charges with each other.
	Thus, if these operators are also in equilibrium
	at $T_{B-L}$, the $B-L$ asymmetry generated by leptogenesis is washed out.
	To avoid such problems, it is required that
	\begin{align}
	\label{eq: lower bound on M_C}
	 \frac{10}{z_{B-L}^{7/4}}\left(\frac{\hat{m}_\nu}{0.1\rm eV}\right)^{1/4} M_R	\lesssim \tilde{M}_C\ .
	\end{align}
	By comparing Eqs.\,\eqref{eq: upper bound on M_C} and \eqref{eq: lower bound on M_C}, we find that the allowed parameter region for the ADM scenario is highly restricted due to the washout interaction when the portal operators are generated from the UV model in Eq.\,\eqref{eq:UVportal}.
	
	This constraint can be easily relaxed by introducing additional $B-L$ portals.
	For example, we may introduce a pair of gauge singlet fermions, $(X,\bar{X})$ 
	with new scalar fields, $H_p$, and $H_{Cp}'$, whose gauge and $B-L$ charges are 
    the same with those of the Higgs doublet of the SM and the dark colored Higgs, respectively.
	In this case, there can be additional operators,
	\begin{align}
	\label{eq:onshellportal}
	    {\cal L}= 
	    M_X X \bar X +\frac{1}{2}M_{H}^2 |H_{p}|^2+\frac{1}{2}M_{Cp}^2 |H_{Cp}'|^2 
	    +y_X L H_p \bar X -Y_X H_{Cp}'\overline D'\bar X \ .
	\end{align}
	Here, $M_X$, $M_{H}$ and $M_{Cp}$ are the mass parameters of $(X,\bar{X})$, $H_p$ and $H_{Cp}'$, respectively, and $y_X$ and $Y_X$ are Yukawa coupling constants.%
	\footnote{$\bar{N}$ and $\bar{X}$ can be distinguished by an approximate discrete symmetry under which $(X,\bar{X})$, $H_p$ and $H_{Cp}'$ are charged. With the discrete symmetry, we can avoid unnecessarily mixing between $\bar{N}$ and $\bar{X}$. 
	}
	As the mass of $\bar{X}$ is the Dirac type, the interaction terms in Eq.\,\eqref{eq:onshellportal} do not violate the $B-L$ symmetry.
	Thus, these interactions do not washout the asymmetry generated by leptogenesis 
	but thermally distribute the asymmetry between the visible and the dark sector for $M_{X,H,Cp} < T_{B-L}$.
	
	
	
	In the following analysis, we divide the parameter region into two. 
	\begin{itemize}
	    \item Off-shell $B-L$ portal scenario:
	    \begin{align}
	 \frac{10}{z_{B-L}^{7/4}}\left(\frac{\hat{m}_\nu}{0.1\rm eV}\right)^{1/4} M_R	\lesssim \tilde M_C\lesssim \frac{10}{z_{B-L}^{5/4}}\left(\frac{\hat{m}_\nu}{0.1\rm eV}\right)^{1/4}M_R\ .
	 \end{align}
	 \item On-Shell $B-L$ portal scenario:
	 	    \begin{align}
	 	    \label{eq:scenario2}
    \frac{10}{z_{B-L}^{5/4}}\left(\frac{\hat{m}_\nu}{0.1\rm eV}\right)^{1/4} M_R\lesssim \tilde M_C \ .
	 \end{align}
	\end{itemize}
In the on-shell portal scenario, we assume that there are lighter particles than $T_{B-L}$ which mediate the $B-L$ asymmetry between two sectors as in Eq.\,\eqref{eq:onshellportal}.%
\footnote{In the on-shell scenario, we may take $Y_N = 0$, and hence, the Majorana dark neutron mass is not inevitable.}
It should be emphasized that the $B-L$ asymmetries in the two sectors are thermally distributed in both the scenarios.%
\footnote{In the absence of the on-shell portal, 
the region with $M_C < T_{B-L}$ results 
in a dark sector asymmetry which depends on the branching ratio
of $\bar{N}$ for small $Y_N$'s~\cite{Falkowski:2011xh}.
If $Y_N$'s are large for $M_C < T_{B-L}$, on the other hand, the $B-L$ asymmetry 
is washed out very strongly and results in too small asymmetry. }
	
	
	
\section{gamma-ray and electron/positron fluxes }\label{sec:signal}


As we have seen in the previous section, the dark neutron obtains a Majorana mass 
when the portal operator is generated in association with the seesaw mechanism.  
Due to the Majorana mass of the dark neutron, the dark neutron can oscillate into the anti-dark neutron. 
The typical time scale of the oscillation, $t_{\rm{osc}}= m_M^{-1}$, is estimated as
\begin{align}
t_{\rm{osc}}
&\simeq
{3.3\times 10^{21}}\,\mbox{sec}\,
\left(\frac{\Lambda'_{\rm{QCD}}}{2\,\mbox{GeV}}\right)^{-6}
\left(\frac{\tilde{M}_C}{3\times 10^9\,\mbox{GeV}}\right)^{4}
\left(\frac{M_R}{10^9\,\mbox{GeV}}\right)\ .
\end{align}
We now see that some fraction of $n'$ can convert into $\bar{n}'$ at  late time,
and then $n'/p'$ and $\bar{n}'$ annihilate into the dark pions. The neutral dark pions decay into the dark photons, and the dark photons finally decay into $e^+e^-$ pairs. $\gamma$ can be also emitted by the final state radiation (FSR) process as depicted in figure \ref{fig:cascade}. 
In this section, we discuss the constraints 
on the late-time annihilation from the 
observations of the $\gamma$-ray from the dSphs and the interstellar $e^++e^-$ flux.

\subsection{Gamma-ray flux from the Dwarf Spheroidal Galaxies}
The $\gamma$-ray signal is one of the most promising channels
to search for dark matter annihilation (e.g., \cite{Gunn:1978gr,Bergstrom:2012fi} for review). 
In particular, dSphs in our galaxy
are the ideal targets to search for the $\gamma$-ray signal,
since they have high dynamical mass-to-light ratios, ($M/L \sim 10-1000$), 
while they lack contaminating astrophysical $\gamma$-ray sources~\cite{Gilmore:2007fy,McConnachie:2012vd}.
In this subsection, we estimate the $\gamma$-ray 
fluxes from the dSphs and compare them with the upper limits 
on the fluxes put by the Fermi-LAT.

First, we calculate the $\gamma$-ray spectrum at production by the $n'\bar{n}'$ annihilation processes:
\begin{align}
n'\bar{n'}
\to
m\pi'^0
+
l\pi'^+ 
+ 
l\pi'^-\ 
,~~~
(m,l = 0,1,2,\cdots)\ .
\end{align}
The cascade spectrum can be calculated by using the technique developed in~\cite{Mardon:2009rc,Elor:2015tva,Elor:2015bho}.

We start to calculate the $\gamma$-ray spectrum at the rest frame of $\gamma'$. 
For $m_{\gamma'} \gg m_e$, the spectrum is given by the Altarelli-Parisi approximation 
formula~\cite{Mardon:2009rc},%
\footnote{In the appendix \ref{sec:FSR}, we compare the direct calculation of
the FSR with the Altarelli-Parisi approximation formula, and confirm the validity of 
the approximation in the parameter region we are interested in.}
\begin{align}
\frac{d\tilde{N}_\gamma}{dx_0}
=
\frac{\alpha_{\rm{EM}}}{\pi}
\frac{1+(1-x_0)^2}{x_0}
\left[-1+\ln\left(\frac{4(1-x_0)}{\epsilon^2_0}\right)\right]\ ,
\label{eq:dNdx0}
\end{align}
where $\epsilon_0={2m_e}/{m_{\gamma'}}$ and
$x_0={2E_0}/{m_{\gamma'}}$ with $E_0$ being the energy of $\gamma$ at the rest frame of $\gamma'$. $\alpha_{\rm{EM}}$ denotes the fine structure constant of SM QED.


\begin{figure}
	\centering
	\includegraphics[width=10cm,clip]{cascade.pdf}
	\caption{
	ADM annihilation which happens at  late time:
	$\bar{n}'$ can be generated from the ADM oscillation.
	Once the $\bar{n}'$ is generated, dark nucleons ($n'/p'$)
	and $\bar{n}'$ annihilate into dark pions ($\pi'^{\pm}$ and $\pi'^0$). 
	$\pi'^0$ subsequently decays into a pair of dark photons ($\gamma'$). 
	$\gamma'$ eventually decays into $e^++e^-$, and emits $\gamma$ through the FSR process.
	}
	\label{fig:cascade}
\end{figure}


The next step is to translate the spectrum in the rest frame of $\gamma'$ to that in the rest frame of $\pi'^0$. 
For the case where $m_{\pi'} \gg m_{\gamma'}$, the spectrum is calculated as
\begin{align}
\frac{d\tilde{N}_{\gamma}}{dx_1}
=
2\int^1_{x_1}
\frac{dx_0}{x_0}
\frac{d\tilde{N}_{\gamma}}{dx_0}
f\left(\frac{2x_1}{x_0}-1\right)
+
\mathcal{O}\left(\frac{m^2_{\gamma'}}{m^2_{\pi'}}\right)\ ,
\label{eq:dNdx1}
\end{align}
where $x_1={2E_1}/{m_{\pi'}}$ with $E_1$ being the energy of $\gamma$ at the rest frame of $\pi'$. 
The function $f$ represents the effect of the anisotropy of the $\gamma'$ decay. 
According to~\cite{Gao:2010qx,Elor:2015tva}, we take
\begin{align}
f(\cos\theta)
=
\frac{3}{8}(1+\cos^2\theta)\ ,
\end{align} 
with $\theta$ being the angle between the $\gamma$ emission line and the boost axis of $\gamma'$.
Note that the angle $\theta$ is kinematically constrained as
\begin{align}
\cos\theta
=
\frac{2x_1}{x_0}-1
+
\mathcal{O}\left(\frac{m^2_{\gamma'}}{m^2_{\pi'}}\right)\ .
\end{align}
This is the reason why we put $f(2x_1/x_0-1)$ in Eq.\,\eqref{eq:dNdx1}.


We next translate the spectrum Eq.\,\eqref{eq:dNdx1} to that in the center of mass (CM) frame for the ADM annihilation.
In order to do that, we need to know how much $\pi'^{0}$ is boosted.
If the total number of the dark pions is two $(m+2l=2)$, 
we can exactly know the energy/boost of the dark pions since they should be emitted back to back in the CM frame.
In this case, the $\gamma$ spectrum is calculated as
\begin{align}
\frac{d\tilde{N}^{(m,l)}_{\gamma}}{dx_2}
=
2\int^1_{x_2}
\frac{dx_1}{x_1}
\frac{d\tilde{N}_{\gamma}}{dx_1}
+
\mathcal{O}\left(\frac{m^2_{\pi'}}{m^2_{\rm{DM}}}\right)\ ,
~~~
\mbox{for}~~m+2l=2\ ,
\label{eq:dNdx2n2}
\end{align}
where $x_2={E_2}/{m_{\rm{DM}}}$ with $E_2$ being the energy of $\gamma$ at the CM frame.
 
On the other hand, in the case of $m+2l\geq3$, 
it becomes highly non-trivial to know how much the $\pi'^0$ can be boosted 
even when we assume that the matrix element of the annihilation is constant as a function 
of the final state momenta.
This is because, in this case, the energy spectrum of the dark pion is given as
\begin{align}
\frac{dN_{\pi'}}{d\xi}
=
\frac{1}{R_n}\frac{dR_n}{d\xi}
 \ ,
\end{align}
where $\xi=E_{\pi'}/m_{\rm{DM}}$ and $R_n$ is the $n=m+2l$ body phase space integration~\cite{Liu:2014cma}.
$E_{\pi'}$ denotes the energy of the dark pion in the CM frame.
In general, it is difficult to perform the phase space integration for $n\geq 3$. 
However, as discussed in \cite{Liu:2014cma,Elor:2015bho}, under the assumption that $m_{\pi'^0}=m_{\pi'^+}\equiv m_{\pi'}\ll m_{\rm{DM}}$, we can perform the phase space integrations analytically as
\begin{align}
\frac{dN_{\pi'}}{d\xi}
=
(n-1)(n-2)(1-\xi)^{n-3}\xi
+
\mathcal{O}\left(\frac{m^2_{\pi'}}{m^2_{\rm{DM}}}\right),
\end{align}
for $n=m+2l\geq3$. Using the results, we finally obtain
\begin{align}
\frac{d\tilde{N}^{(m,l)}_{\gamma}}{dx_2}
=
2(n-1)(n-2)
\int^1_{x_2}d\xi
(1-\xi)^{n-3}
\int^1_{x_2/\xi}
\frac{dx_1}{x_1}
\frac{d\tilde{N}_{\gamma}}{dx_1}
+
\mathcal{O}\left(\frac{m^2_{\pi'}}{m^2_{\rm{DM}}}\right)
 ,
\label{eq:dNdx2}
\end{align}
for $n=m+2l \geq 3$ where we assume $m_{\pi'^0}=m_{\pi'^+}\equiv m_{\pi'}$.

Finally, we sum over the possible intermediate states
and take into account the number of the final states.
It turns out that the total $\gamma$ spectrum from the $n'\bar{n}'$ annihilation is expressed as
\begin{align}
\frac{dN^{(n'\bar{n}')}_\gamma}{dx_2}
=
\sum_{m,l}2m
\left(
\mbox{Br}^{(n'\bar{n}')}(m,l)
\frac{d\tilde{N}^{(m,l)}_{\gamma}}{dx_2}
\right)
\ ,
\label{eq:dNnn}
\end{align}
where $\mbox{Br}^{(n'\bar{n}')}(m,l)$ denotes the branching ratio for the $n'\bar{n'}\to m\pi'^0+l\pi'^+ + l\pi'^-$ annihilation process. The factor $2m$ corresponds to the number of $e^+e^-$ pairs in the annihilation process.

In the same way, we can estimate the $\gamma$ spectrum from the $p'\bar{n}'$ annihilation processes:
\begin{align}
p'\bar{n'}
\to
m\pi'^0
+
l\pi'^+ 
+ 
(l-1)\pi'^-
\ ,~~~
(m = 0,1,2,\cdots, l=1,2,\cdots)\ .
\end{align}
The $\gamma$ spectrum is calculated as
\begin{align}
\frac{dN^{(p'\bar{n}')}_\gamma}{dx_2}
=
\sum_{m,l}2m
\left(
\mbox{Br}^{(p'\bar{n}')}(m,l)
\frac{d\tilde{N}^{(m,l)}_{\gamma}}{dx_2}
\right)\ ,
\label{eq:dNpn}
\end{align}
with replacing $n=m+2l$ by $n=m+2l-1$ in the calculation of ${d\tilde{N}^{(m,l)}_{\gamma}}/{dx_2}$.

In the following analysis, we simply assume that the branching ratio of the dark nucleon annihilation can be estimated as that of nucleon-antinucleon annihilation. According to~\cite{Orfanidis:1973ix}, 
we approximate the branching ratios by the fireball model,%
\footnote{In this approximation, the Parity violating mode, $(m,l)=(2,0)$, is allowed, 
although it is not significant numerically.} 
\begin{align}
\mbox{Br}^{(n'\bar{n}')}(m,l)
&=
\frac{2\alpha^{2l}{} }{(1+\alpha)^n+(1-\alpha)^n}\,
_n C _{2l}\,
P_n\ ,~~\mbox{with}~~
n={m+2l}\ ,\\
\mbox{Br}^{(p'\bar{n}')}(m,l)
&=
\frac{2\alpha^{2l-1}{} }{(1+\alpha)^n+(1-\alpha)^n}\,
_n C _{2l-1}\,
P_n\ ,~~\mbox{with}~~
n={m+2l-1}\ ,
\end{align}
where 
\begin{align}
P_n
=
\frac{1}{\sqrt{2\pi}\sigma}\exp
\left(-\frac{(n-\langle{n}\rangle)^2}{2\sigma^2}\right)\ ,
\end{align}
with $a=1/4$, $\langle{n}\rangle=5.05$, $\sigma^2=a\langle{n}\rangle$ and
\begin{align}
\alpha
 =
   \left\{
     \begin{array}{l}
       \sqrt{2}~~~\mbox{for}~~n=2\ ,\\
       1.5~~~\mbox{for}~~n\neq2 \ .
     \end{array}
   \right.
\end{align}


We are now ready to estimate the $\gamma$-ray spectrum emitted from the ADM annihilation. 
Figure \ref{fig:Spectra} shows the value of the $\gamma$-ray spectrum. 
Here, we take $m_{\mathrm{DM}}=10$\,GeV, $m_{\pi'}=1\,\mbox{GeV}$ and $m_{\gamma'}=40\,\mbox{MeV}$.
The black solid and the dashed lines correspond 
to the spectra predicted from the $n'\bar{n}'$ and $p'\bar{n}'$ annihilation, respectively.
In the analysis, we ignore the contributions from the annihilation with large $(m,l)$ since the branching ratios of them are much suppressed. 
We stop taking the sum over $(m,l)$ if the size of contribution 
is less than $1\%$ of the total amount.

The figure shows that the ADM annihilation predict the continuous $\gamma$-ray spectrum peaked at the energy 
of ${\cal O}(m_{\rm{DM}}/10)$.
This is expected as the typical number of the dark pions for an annihilation is five, and the neutral dark pion
decays into two pairs of $e^+e^-$.

It should be reminded that the $\gamma$-ray emission from the ADM annihilation can happen at the present universe since the ADM
oscillation effectively happens at the late time scale.
The ADM signals can  therefore be tested by $\gamma$-ray telescope experiments from nearby sources, 
while evading the constraints from the observations of the cosmic microwave observations (see e.g. \cite{Elor:2015bho}).

\begin{figure}
	\centering
	\includegraphics[width=6.3cm,clip]{Spectra.pdf}
	\caption{
	The $\gamma$-ray spectrum at production from the $n'\bar{n}'$ (Solid line) and $p'\bar{n}'$ (Dashed line) annihilation. 
	We take $m_{\pi'}=1\,\mbox{GeV}$ and $m_{\gamma'}=40\,\mbox{MeV}$.
	}
	\label{fig:Spectra}
\end{figure}


The $\gamma$-ray flux from the dSphs for an energy bin from $E_{\rm{min}}$ to  $E_{\rm{max}}$ is calculated as
\begin{align}
\Phi
=
\int^{E_{\rm{max}}}_{E_{\rm{min}}}dEE
\int_{\Delta\Omega}\frac{d\Omega}{4\pi}
\int_{\rm{l.o.s.}} dl
\left(
n_{n'}n_{\bar{n}'}\langle{\sigma v}\rangle_{n'\bar{n}'}\frac{dN^{(n'\bar{n}')}_\gamma}{dE}
+
n_{p'}n_{\bar{n}'}\langle{\sigma v}\rangle_{p'\bar{n}'}\frac{dN^{(p'\bar{n}')}_\gamma}{dE}
\right)\ ,
\end{align}
where we perform the integrations over
a solid angle, $\Delta\Omega$, and the line-of-sight (l.o.s.).
Here $n_i$ and $\langle{\sigma v}\rangle_{ij}$ denote the number density of a particle $i$ at the dSphs
and the kinematically averaged cross section for $ij$ annihilation, respectively. $N^{(n'\bar{n}')}_\gamma$ and $N^{(p'\bar{n}')}_\gamma$ are the photon spectra from $n'\bar{n}'$ and $p'\bar{n}'$ annihilation which can be calculated from Eqs.\,\eqref{eq:dNnn} and \eqref{eq:dNpn}, respectively.


It should be noted that 
the total amount of the $\gamma$-ray flux can be large enough to be tested by the $\gamma$-ray searches on the dSphs
although the flux is suppressed by the factor,
\begin{align}
\frac{n_{\bar{n}'}}{n_{n'}}
&\simeq
\left(\frac{t_{0}}{t_{\rm{osc}}}\right)^2\nonumber\\
&\simeq
1.6\times 10^{-8}
\left(\frac{\Lambda'_{\rm{QCD}}}{2\,\mbox{GeV}}\right)^{12}
\left(\frac{\tilde{M}_C}{3\times 10^9\,\mbox{GeV}}\right)^{-8}
\left(\frac{M_R}{10^9\,\mbox{GeV}}\right)^{-2}\ .
\end{align}
where $t_0\simeq 4.3\times 10^{17}\,\rm{sec}$ is the age of the universe. 
This is because the thermally-averaged cross section can be large due to the strong interaction. 
In the following analysis, we take the annihilation cross sections to be
\begin{align}
\langle{\sigma v}\rangle_{n'\bar{n}'}
=
\langle{\sigma v}\rangle_{p'\bar{n}'}
=
\frac{4\pi}{m^2_{\rm{DM}}}\ ,
\label{eq:xsec}
\end{align}
to give rough estimation.
Such a large annihilation cross section multiplied by the relative velocity is supported 
by the cross section measurements of the non-relativistic nucleon and anti-nucleon 
annihilation~\cite{Armstrong:1987nu,Bertin:1997gn} 
(see also \cite{Huo:2015nwa,Lee:2015hma}).%
\footnote{The cross section multiplied by the relative velocity in Eq.\,\eqref{eq:xsec} 
is much smaller than the unitarity limit.

In the Appendix~\ref{sec:sommerfeld}, we discuss
the Sommerfeld enhancement effects by the exchange of the dark pions.
There, we find that the enhancement effects are not 
significant in the present setup.
}

\begin{figure}
	\centering
	\includegraphics[width=13cm,clip]{flux-draco.pdf}
	\caption{
	The predicted $\gamma$-ray flux from the Draco dSph and the $95\%$ C.L. upper bound obtained by the Fermi-LAT (green line). The black solid, dashed, and dotted lines correspond to the $\gamma$-ray flux when we take $\tilde{M}_C=10^9\,\mbox{GeV}$, $1.5\times 10^9\,\mbox{GeV}$, and $2\times 10^9\,\mbox{GeV}$, respectively. Here, we assume $m_{n'}=m_{p'}=m_{\rm{DM}}$ and fix $m_{\pi'}=1\,\mbox{GeV}$, $M_{R}=5\times 10^9\,\mbox{GeV}$ and $m_{\gamma'}=40\,\mbox{MeV}$. 
	}
	\label{fig:draco}
\end{figure}

 In Figure~\ref{fig:draco}, we show the predicted $\gamma$-ray flux from the Draco dSph. 
 The black solid, dashed, and dotted lines correspond to the $\gamma$-ray flux when we take $\tilde{M}_C=10^9\,\mbox{GeV}$, $1.5\times 10^9\,\mbox{GeV}$, and $2\times 10^9\,\mbox{GeV}$, respectively. Here, we assume $m_{n'}=m_{p'}=m_{\rm{DM}}=10\,\mbox{GeV}$ and fix $m_{\pi'}=1\,\mbox{GeV}$ and $m_{\gamma'}=40\,\mbox{MeV}$. 
 To obtain the predicted $\gamma$-ray spectrum, we use the $J$-factors estimated in \cite{Hayashi:2016kcy}
which takes into account the effects of the non-sphericity of the dSphs.%
\footnote{As for the $J$-factor of the Ursa Minor classical dSphs, we use the value given in
\cite{Geringer-Sameth:2014yza} as it is not analyzed in \cite{Hayashi:2016kcy}.
}
The green line corresponds to the upper bound ($95\%$ C.L.) on the $\gamma$-ray flux based on the  
6 years of Pass 8 data by the Fermi-LAT collaboration~\cite{Ackermann:2015zua}.
The figure shows that the $\gamma$-ray flux from the 
late-time annihilation becomes comparable to the 
upper limit on the observed flux for $\tilde{M}_C = \order{10^9}$\,GeV and $M_R = \order{10^{10}}$\,GeV, which corresponds to the oscillation time scale of 
$t_{\mathrm{osc}}=\order{10^{21}}\,$sec.
We discuss the constraints on the model parameters 
by the Fermi-LAT in subsection~\ref{sec:constraint}.

\subsection{Interstellar Electron/Positron Flux}
The Fermi-LAT observation does not constrain 
the late-time annihilation for $m_{\mathrm{DM}}\lesssim 3$\,GeV, 
since the Fermi-LAT is sensitive to the $\gamma$-ray with 
energy higher than $500$\,MeV.
For such a rather light ADM, 
the most stringent constraint is put by 
the observation of the interstellar $e^++e^-$ flux by the Voyager-1~\cite{Fisk:1976aw,Stone2013} (see also \cite{Boudaud:2016mos}).
In this subsection, we estimate the $e^++e^-$ flux from the late-time annihilation in the Milky Way.

The energy spectrum of $e^++e^-$ at production 
by the late-time ADM 
annihilation is obtained by replacing $d\tilde{N}_\gamma/dx_0$ in Eq.\,\eqref{eq:dNdx0} with
the $e^+/e^-$ spectrum in the dark photon rest frame,
\begin{align}
    \frac{d\tilde{N}_e}{dx_0} = \delta(x_0 - 1)\ .
\end{align}
Here, $x_0 = 2E_0/m_{\gamma'}$ with $E_0$ being 
the energy of either $e^-$ or $e^+$.
By repeating the same analysis in the previous section,
we can convert this spectrum to the one in the rest frame 
of the ADM annihilation.
In Figure~\ref{fig:SpectraEP}, we show the $e^+/e^-$ spectrum at production for 	$m_e\ll m_{\gamma'}\ll m_{\pi'}\ll m_{\mathrm{DM}}$. 

\begin{figure}
	\centering
	\includegraphics[width=6.3cm,clip]{SpectraEP.pdf}
	\caption{
	The $e^-/e^+$ spectrum predicted from the $n'\bar{n}'$ (Solid line) and $p'\bar{n}'$ (Dashed line) annihilation for 
	$m_e\ll m_{\gamma'}\ll m_{\pi'}\ll m_{\mathrm{DM}}$. 
	}
	\label{fig:SpectraEP}
\end{figure}

For a given $e^+/e^-$ spectra at production, the interstellar $e^++e^-$ flux at around the location of the Earth is 
given by \cite{Cirelli:2010xx,Buch:2015iya},%
\footnote{A typical propagation time of the 
cosmic ray to travel of $\order{1}$\,kpc is 
much shorter than the 
age of the universe.}
\begin{align}
\frac{d\Phi}{dE_e} = \frac{1}{4\pi b(E)}
\left(\frac{\rho_{\mathrm{DM}}}{m_{\mathrm{DM}}}\right)^2
\left(\frac{t_0}{t_{\mathrm{osc}}}\right)^2
\sum_{i=n'\bar{n}',p'\bar{n}'}\langle\sigma v \rangle_{i} \int_E^{E_{\rm DM}}dE_s\,
I(E,E_s)\,\frac{dN_{e\,i}}{dE}(E_s)\ .
\end{align}
Here, $\rho_{\mathrm{DM}}$ denotes a local dark matter 
density at around the location of the Earth, 
$I(E,E_s)$ is a Green function which encodes the propagation of 
$e^\pm$ from a source with a given energy $E_s$ to 
any energy $E$, and $b(E)$ is the $e^\pm$ energy loss function.%
\footnote{The Green function is dimensionless while  $b(E)$ has a unit of GeV/sec
which is typically 
$b(E)\simeq 10^{-(16-15)}$\,GeV/sec 
for $E=\order{10}$\,MeV to $\order{1}$\,GeV~\cite{Cirelli:2010xx,Buch:2015iya}.
}
\begin{figure}
	\centering
	\includegraphics[width=8cm,clip]{CR.pdf}
	\caption{The interstellar $e^++e^-$ flux in cosmic ray at around the location of the Earth from the late-time ADM annihilation. 
	The annihilation cross section and the osicllation time scale is taken to be 
	$(t_0/t_{\mathrm{osc}})^2\times\langle\sigma v\rangle = 1$\,pb ($t_{\mathrm{osc}}\gg t_0$).
	The green, red, and blue lines show the spectrum for $m_{\mathrm{DM}} = 3$\,MeV, $5$\,MeV, and $10$\,MeV, respectively. 
	The solid lines assume the MED propagation model,
	while the upper and the lower dotted lines assume 
	the MAX and the MIN propagation models, respectively.
	The NFW halo profile is used.
	}
	\label{fig:CR}
\end{figure}

In Figure~\ref{fig:CR}, we show the 
interstellar $e^++e^-$ flux at around the location of the Earth from the late-time ADM annihilation. 
Here, the annihilation cross section and the oscillation time scale is set to be $(t_0/t_{\mathrm{osc}})^2\times\langle\sigma v\rangle = 1$\,pb.
The Green function, $I(E,E_s)$, and the energy loss rate, $b(E)$, are those provided by \cite{Cirelli:2010xx,Buch:2015iya}.
In the figure, the solid lines assume the MED propagation model,
while the upper and the lower dotted lines assume 
the MAX and the MIN propagation models, respectively (see \cite{Donato:2003xg}).
The dark matter profile is assumed to be the NFW profile~\cite{Navarro:1996gj},%
\footnote{We numerically checked that the
spectra are not significantly changed even 
for a cored Burkert profile \cite{Burkert:1995yz},
though they are slightly suppressed.}
with the local dark matter density at around the Earth
to be $\rho_{\mathrm{DM}} = 0.3$\,GeV/cm$^3$.

In the figure, we also show the interstellar $e^++e^-$ spectrum observed by the Voyager-1~\cite{Fisk:1976aw,Stone2013}, where the data is taken from \cite{Maurin:2013lwa}. 
The figure shows that the $e^++e^-$ flux 
from the late-time ADM annihilation is much smaller
than the observed flux for $(t_{0}/t_{\mathrm{osc}})^2\times \langle \sigma v\rangle = \order{1}$\,pb.
We will summarize the constraints from the Voyger-1 
in the next subsection.

\subsection{Constraints on Parameter Space}
\label{sec:constraint}
\begin{figure}
	\centering
	\includegraphics[width=8cm,clip]{tosc_const.pdf}
	\caption{
	Constraints on the oscillation time scale. 
	The green and the yellow shaded regions are excluded by the $\gamma$-ray constraint (Fermi-LAT),
	and by the constraint on the $e^++e^-$ flux (Voyager-1), respectively. 
	We here assume $m_{n'}=m_{p'}=m_{\rm{DM}}$ and take $n_{n'}=n_{p'}$.
	We also assume $m_{\pi'}\ll m_{\mathrm{DM}}$ while fixing $m_{\gamma'}=40\,\mbox{MeV}$. 
	The red dotted line is a prospected lower limit by the $\gamma$-ray search 
	from the dSphs by e-ASTROGAM in one year of effective exposure.
	}
	\label{fig:tosc}
\end{figure}


As we have seen in the previous subsections, 
we can probe the time scale of the matter-antimatter oscillation by the $\gamma$-ray observation up to $t_{\mathrm{osc}}=\order{10^{21}}$\,sec for $m_{\mathrm{DM}}\simeq 10$\,GeV. 
This oscillation time scale corresponds to the effective annihilation cross section,%
\footnote{The effective cross section into the $\gamma$-ray is further suppressed by Eq.\,\eqref{eq:dNdx0}.}
\begin{align}
    \left(\frac{t_0}{t_{\rm osc}}\right)^2 \langle\sigma v\rangle \sim 10\,{\rm pb} \left(\frac{10\,{\rm GeV}}{m_{\rm DM}}\right)^2
    \left(\frac{10^{21}\rm {sec}}{t_{\rm osc}}\right)^2\ .
    \label{eq:effectiveCS}
    \end{align}
A lighter ADM can be also tested by the observation of the interstellar $e^++e^-$ flux.


In Figure~\ref{fig:tosc}, we show the constraints on the oscillation time scale from 
the observations by the Fermi-LAT and the Voyager-1.
Here, we assume $m_e\ll m_{\gamma'} \ll m_{\pi'}\ll m_{\mathrm{DM}}$ while we fix $m_{\gamma'}=40$\,MeV.%
\footnote{The constraints do not depend on $m_{\gamma'}$ significantly, as long as $m_e\ll m_{\gamma'} \ll m_{\pi'}\ll m_{\mathrm{DM}}$.}
The green region corresponds to the $95\%$ C.L. excluded region from the Fermi-LAT 
observations (see also \cite{Fermi-LAT:2016uux,Elor:2015bho}),
where we take into account the $\gamma$-ray fluxes
from the 8-classical dSphs.
The yellow shaded region corresponds to the $95\%$ C.L. excluded region from the Voyager-1 observation
for the MED propagation model with the NFW dark halo profile.
We see that, for $m_{\rm{DM}}\simeq 5\mbox{--} 10\,\mbox{GeV}$, the more stringent 
constraints are put by the Fermi-LAT observation, where the 
oscillation time scale shorter than $t_{\mathrm{osc}}\sim 10^{21}\,\mbox{sec}$ is excluded.
For a lighter mass region,
the Voyager-1 observation excludes the oscillation time scale shorter than 
$t_{\mathrm{osc}}\sim 10^{21-22}$\,sec.


\begin{figure}[t]
	\begin{minipage}{0.4\hsize}
	\centering
	\includegraphics[width=0.95\hsize,clip]{MRMC2.pdf}
    \end{minipage}
    \begin{minipage}{0.4\hsize}
	\centering
	\includegraphics[width=0.95\hsize,clip]{MRMC5.pdf}
    \end{minipage}
    \begin{minipage}{0.4\hsize}
	\centering
	\includegraphics[width=0.95\hsize,clip]{MRMC10.pdf}
    \end{minipage}
	\caption{
	Constraint on our ADM scenario for $m_{\mathrm{DM}}=2,5,10$\,GeV.
	The green and yellow shaded regions correspond the $95\%$ C.L. excluded regions by 
	the Fermi-LAT and the Voyager-1, respectively. 
	The lower gray region is excluded 
	in which the $B-L$ asymmetry is washed out (see Eq.\,\eqref{eq: lower bound on M_C}).
	Above the solid line, we require an on-shell $B-L$ portal sector (see Eq.\,\eqref{eq:scenario2}).
	We also show the prospected limits by e-ASTROGAM translated from Figure~\ref{fig:tosc}.}
	\label{fig:limit}
\end{figure}

In Figure~\ref{fig:limit},
we translate the constraints on the oscillation time scale to those on the 
parameters of the present model.
In the figure, we consider $m_{\mathrm{DM}}=2,5,10$\,GeV.
We also take  $\Lambda'_{\rm{QCD}}=2\,\mbox{GeV}\times(m_{\mathrm{DM}}/10\,\mathrm{GeV})$
to mimic QCD for each choice of the dark matter mass.
We also assume  $m_e \ll m_{\gamma'} \ll m_{\pi'}\ll m_{\mathrm{DM}}$.
The green and yellow shaded regions correspond the $95\%$ C.L. excluded regions by 
the Fermi-LAT and the Voyager-1, respectively. 
The lower gray region is excluded 
where the $B-L$ asymmetry is washed out (see Eq.\,\eqref{eq: lower bound on M_C}).
Above the solid line, we require an on-shell $B-L$ portal sector (see Eq.\,\eqref{eq:scenario2}).


We now see that the composite ADM scenario 
with $M_C = {\mathcal O}(10^9)$\,GeV can be tested by the $\gamma$-ray searches from the dSphs by the Fermi-LAT for $m_{\rm DM}\simeq 10$\,GeV.
Even for a lighter ADM scenario, we see that the region with $M_C = \order{10^8}$\,GeV
has been excluded by the Voyager-1 observation.
The resultant constraint is important in view of the fact that the parameter region with $M_C\sim 10^{9}$--$10^{10}$\,GeV
is highly motivated in the UV completion 
model based on $SU(4) \,(\supset SU(3)_{D}\times U(1)_D)$ gauge theory~\cite{Ibe:2018juk,Ibe:2018tex,Ibe:2019ena}.
In this UV completion, the tiny kinetic mixing of $\epsilon = 10^{-10}$--$10^{-9}$ which evades all the phenomenological constraints on the dark photon~\cite{Ibe:2018juk}
is achieved when the $SU(4)$ breaking scale 
is at around $10^{9}$--$10^{10}$\,GeV.
The $SU(4)$ breaking scale also leads 
to the colored dark Higgs mass in a similar range.
The $\gamma$-ray searches are already sensitive to such a well-motivated parameter region 
for $m_{\mathrm{DM}}\simeq 10$\,GeV.


Several comments are in order.
In our discussion, we consider only the $\gamma$-ray emitted by the FSR. 
This should be justified as the $\gamma$-rays made by the Synchrotron radiation and the inverse Compton scattering 
from the sub-GeV $e^+/e^-$ are very soft and below the Fermi-LAT sensitivity~\cite{Cirelli:2010xx}.
It should be also noted that the $\gamma$-ray signal from the galactic center does not lead to more stringent constraints,
despite the signal strength is higher than that from the dSphs.
This is because the $\gamma$-ray background is much higher 
for the galactic center, and hence, it is difficult to 
distinguish the continuous signal spectrum from the background spectrum.

Future $\gamma$-ray searches 
such as 
e-ASTROGAM~\cite{DeAngelis:2017gra,Rando:2019fzq},
SMILE~\cite{smile},
GRAINE~\cite{Aoki:2012nn},
and 
GRAMS~\cite{Aramaki:2019bpi}
projects will be important to test the model further.
It should be emphasized that those experiments 
are sensitive to the MeV $\gamma$-rays, and hence,
they are also able to test the models with $m_{\mathrm{DM}}=$\,a few GeV to which the Fermi-LAT loses 
sensitivity.
In Figure~\ref{fig:tosc}, we show the prospected lower limit on $t_{\mathrm{osc}}$ at 95\%CL 
by the  $\gamma$-ray search 
from the dSphs by e-ASTROGAM in one year of effective exposure.
In our analysis, we used the effective area and the prospected sensitivities for a $\gamma$-ray flux
from a point-like source at a high latitude (in Galactic coordinates) in \cite{DeAngelis:2017gra}.
The testable parameter region can be wider when the $J$-factors of the 
ultra-faint dSphs are determined more precisely by 
future spectroscopic observations such as the Prime Focus Spectrograph~\cite{Ellis:2012rn}.
For example, if the $J$-factor of Triangulum II converges to the central value in \cite{Hayashi:2016kcy},
i.e. $\log_{10} J \simeq 20$, the prospected lower limit on $t_{\mathrm{osc}}$
becomes higher for about a factor of $2^{1/2}$.



\section{Conclusions}
The composite ADM model is particularly motivated as it provides the DM mass of ${\cal O}(1)$\,GeV 
and a large annihilation cross section simultaneously.
In this paper, we discussed the indirect detection of the composite ADM where the portal operators 
of the $B-L$ asymmetry is generated in association with the seesaw mechanism.
In this model, the dark-neutron obtains a tiny Majorana mass,
and hence, ADM can pair-annihilate at later times.

As we have discussed, the late time annihilation of ADM results in multiple soft electrons/positrons and soft photons emitted as the FSR.
As a result, some parameter region of the composite ADM which 
is motivated by thermal leptogenesis and dark UV completion models
has been excluded by the Fermi-LAT and the Voyager-1 observations.
The obtained constraint is tighter than that from  the anti-neutrino flux  made by the 
decay of  ADM via the $B-L$ portal operator~\cite{Fukuda:2014xqa} (see Eq.\,\eqref{eq:lifetime}).
Future experiments which are 
sensitive to sub-GeV $\gamma$-rays such as 
e-ASTROGAM~\cite{DeAngelis:2017gra,Rando:2019fzq},
SMILE~\cite{smile},
GRAINE~\cite{Aoki:2012nn},
and 
GRAMS~\cite{Aramaki:2019bpi}
projects will be important to test the oscillating ADM model
further.

\section*{Acknowledgements}
\vspace{-.3cm}
MI thanks K.~Kohri for useful comments on the possible constraints on the model from the cosmic electron/positron rays.
This work is supported by JSPS KAKENHI Grant Numbers, 19K14701 (R. N.),  No. 15H05889, 
No. 16H03991, NO. 17H02878 
and No. 18H05542 (M. I.) and by World Premier International Research Center Initiative (WPI Initiative), MEXT, Japan. 
This work is also supported by the Advanced Leading Graduate Course for Photon Science (S. K.).
