Peaky Production of
Light Dark Photon Dark Matter

Let us first consider the background dynamics of the model. Throughout the present work, we assume that there is no vacuum expectation value of the vector field:¹¹1This assumption in fact depends on the scales of interest. Even if the classical background of $A_{\mu}^{\prime}$ is absent, the dark photon field can be generated from the quantum vacuum, as in our current consideration. If the produced field has a wavelength larger than a certain scale with a nonzero variance, then each local patch smaller than that stochastically experiences a coherent field of dark photon wave, which, strictly speaking, breaks spatial isotropy [45]. However, we are here interested in generation of the dark photon on scales smaller than the CMB ones, and therefore the vanishing $1$-point function is a relevant assumption, at least under the consideration of the isocurvature constraints which is the main observational obstacle for dark photon production of long wavelength modes.

$$\langle A^{\prime}_{\mu}\rangle=0\ .$$

(2.2)

Thanks to this assumption, and ensuring a posteriori that the energy density of the dark photon is a subdominant component, we can consider the flat Friedmann-Lemaître-Robertson-Walker (FLRW) metric,

$$ds^{2}=-dt^{2}+a(t)^{2}d\bm{x}^{2}=a(\tau)^{2}(-d\tau^{2}+d\bm{x}^{2})\ ,$$

(2.3)

as the background spacetime. Here, $a$ is the scale factor and $\tau$ denotes the conformal time, $d\tau=dt/a$. Throughout the paper, we assume that inflation proceeds as de Sitter essentially, and the time dependence of the scale factor is well approximated by $a\simeq\exp\left(Ht\right)\simeq-1/(H\tau)$, where $H\equiv\partial_{t}a/a$ is the Hubble rate during inflation.

For the kinetic function $I$ of $\sigma$, we adopt the following simple exponential form:

$$I(\sigma)=I_{0}\exp\left(\dfrac{\sigma}{\Lambda_{1}}\right)\,,$$

(2.4)

with some dimension-$1$ energy scale $\Lambda_{1}$. The normalization constant $I_{0}$ is assumed to be $I_{0}=1$ to realize $I(\sigma\rightarrow 0)=1$ after inflation ends when $\sigma$ settles into its potential minimum and recovers the canonical kinetic term of the dark photon field. The crucial parameter that dictates production of the dark photon, especially of the transverse modes, is the rate of the time variation of $I$ with respect to the number of e-foldings $N=\ln a$, i.e.,

$$n\equiv\dfrac{dI/dN}{I}=\frac{d\sigma/dN}{\Lambda_{1}}\ ,\qquad dN=Hdt\,.$$

(2.5)

The critical values of $n$ for dark photon production are $|n|=2$: for $|n|>2$, the associated dark electric or magnetic field, or both, goes through exponential amplifications [52].

With the form of $I$ fixed by (2.4), we take the freedom to choose the spectator potential $V(\sigma)$ so that the dark photon production occurs with nontrivial time dependence. In this regard, we consider a toy model,

$$V(\sigma)=\mu^{4}\tanh^{2}\left(\dfrac{\sigma}{\Lambda_{2}}\right)\ ,$$

(2.6)

where $\mu$ is a potential energy scale and $\Lambda_{2}$ is another energy scale which characterizes the steepness of the potential. This kind of functional form is sometimes considered in the context of the $\alpha$-attractor model [53], while in our current study we assume that another independent sector drives inflation. To characterize the energy density of the spectator field $\sigma$, playing no role in the inflationary expansion, we define the following ratio,

$$R\equiv\dfrac{\mu^{4}}{3M_{\rm Pl}^{2}H^{2}}\ll 1\ .$$

(2.7)

The Hubble slow-roll parameter $\epsilon_{H}\equiv-\dot{H}/H^{2}$ (dot represents derivative with respect to the physical time $t$) is given by

$$\epsilon_{H}=\epsilon_{\phi}+\epsilon_{\sigma}+\epsilon_{\gamma^{\prime}}\ll 1\ ,$$

(2.8)

where each contribution is defined as $\epsilon_{\alpha}\equiv\left(\rho_{\alpha}+p_{\alpha}\right)/(2M_{\rm Pl}^{2}H^{2})$ ($\alpha=\phi,\,\sigma,\,\gamma^{\prime}$), where $\rho_{\alpha}$ and $p_{\alpha}$ are the corresponding energy density and pressure, respectively, and reads

	$\displaystyle\epsilon_{\phi}$	$\displaystyle=\dfrac{\dot{\phi}^{2}}{2M_{\rm Pl}^{2}H^{2}}\ ,\quad\epsilon_{\sigma}=\dfrac{\dot{\sigma}^{2}}{2M_{\rm Pl}^{2}H^{2}}\ ,$		(2.9)
	$\displaystyle\epsilon_{\gamma^{\prime}}$	$\displaystyle=\dfrac{1}{M_{\rm Pl}^{2}H^{2}}\left\langle\dfrac{1}{3a^{2}}I^{2}\left(F^{\prime}_{0i}F^{\prime}_{0i}+\dfrac{1}{2a^{2}}F^{\prime}_{ij}F^{\prime}_{ij}\right)+\dfrac{1}{2}m_{\gamma^{\prime}}^{2}\left(A{{}^{\prime}}_{0}^{2}+\dfrac{1}{3a^{2}}A{{}^{\prime}}_{i}A{{}^{\prime}}_{i}\right)\right\rangle\ .$		(2.10)

Here, $\epsilon_{\gamma^{\prime}}$ denotes a contribution from the backreaction of the dark photon field and we assume that its magnitude is small in comparison with the other contributions. Note that, in (2.10), the expression assumes the physical time (not conformal time), and the bracket $\langle\bullet\rangle$ denotes the vacuum average. For our analytical convenience, we take the Hubble parameter as a constant value $H=\text{const}.$ to solve the dynamics of the spectator-dark photon system.

The equation of motion for the spectator field is

$$\dfrac{d^{2}\sigma}{dN^{2}}+3\dfrac{d\sigma}{dN}+\dfrac{V_{\sigma}}{H^{2}}=-\dfrac{II_{\sigma}}{2H^{2}}\langle F^{\prime}_{\mu\nu}F{{}^{\prime}}^{\mu\nu}\rangle\ ,$$

(2.11)

where the subscript $\sigma$ denotes derivative with respect to $\sigma$, and the time dependence of $H$ is assumed to be negligible ($|\partial_{N}H/H|\ll 3$). The right-hand side describes a backreaction effect on the motion of the spectator field. We consider the parameter space where this effect is negligible as we will discuss later. Then, assuming that the first term of the left-hand side in Eq. (2.11) is also negligible, we can use the slow-roll condition $3d\sigma/dN\simeq-V_{\sigma}/H^{2}$ and obtain from (2.5)

$$n\simeq-\dfrac{V_{\sigma}}{3H^{2}\Lambda_{1}}=-2cR\dfrac{M_{\rm Pl}^{2}}{\Lambda_{1}^{2}}\dfrac{\tanh(\tilde{\sigma})}{\cosh^{2}(\tilde{\sigma})}\ .$$

(2.12)

Here, we have defined dimensionless quantities, $c\equiv\Lambda_{1}/\Lambda_{2}$ and $\tilde{\sigma}\equiv\sigma/\Lambda_{2}$. Regarding the slow-roll solution of $\sigma$, we find that $\tilde{\sigma}$ and $n$ are described in terms of the Lambert W-function as

	$\displaystyle\tilde{\sigma}(N)$	$\displaystyle=\text{arcsinh}\left[\sqrt{W(y)}\right]\ ,$		(2.13)
	$\displaystyle n(N)$	$\displaystyle=-2cR\,\dfrac{M_{\rm Pl}^{2}}{\Lambda_{1}^{2}}\,\dfrac{\sqrt{W(y)}}{\left[1+W(y)\right]^{3/2}}\ ,$		(2.14)
	$\displaystyle y$	$\displaystyle\equiv-4R\,\frac{M_{\rm Pl}^{2}}{\Lambda_{2}^{2}}\,N+2\ln\left[\sinh(\tilde{\sigma}_{i})\right]+\sinh^{2}(\tilde{\sigma}_{i})\ ,$		(2.15)

where $\tilde{\sigma}_{i}=\tilde{\sigma}(0)>0$ is an initial value of $\tilde{\sigma}$. Inversely, using the definition of the W-function $W(ze^{z})=z$, the number of e-foldings is written as

$$N=\dfrac{\Lambda_{2}^{2}}{4RM_{\rm Pl}^{2}}\left[2\ln\left[\dfrac{\sinh(\tilde{\sigma}_{i})}{\sinh(\tilde{\sigma})}\right]+\sinh^{2}(\tilde{\sigma}_{i})-\sinh^{2}(\tilde{\sigma})\right]\ .$$

(2.16)

Note that, since we are requiring the weak coupling of the dark photon to the hidden sector, i.e. $I>1$, during inflation, we choose the branch $\sigma>0$, provided $\Lambda_{1,2}>0$.

We are interested in a situation where the magnitude of $|n|$ is initially small but becomes greater than 2 at some point during inflation. Namely, the dark photon production occurs at the intermediate stage of inflation, and its overproduction on large scales is avoidable. To realize such a situation, we consider the following scenario: the spectator field rolls down on the potential and its velocity experiences a peak value at an intermediate time of inflation (see Figure 1). At an initial stage, $\tilde{\sigma}$ is near the top of the potential and $|n|<2$ is satisfied. Then, the velocity of the scalar field gradually increases and gets maximized when it goes through an inflection point $\tilde{\sigma}=\tilde{\sigma}_{\rm IP}$ on the potential,

$$V_{\sigma\sigma}(\tilde{\sigma}_{\rm IP})=2\,\dfrac{\mu^{4}}{\Lambda_{2}^{2}}\,\dfrac{1-2\sinh^{2}(\tilde{\sigma}_{\rm IP})}{\cosh^{4}(\tilde{\sigma}_{\rm IP})}=0\qquad\longleftrightarrow\qquad\tilde{\sigma}_{\rm IP}\simeq 0.658\,,$$

(2.17)

at a certain number of e-foldings $N=N_{\rm IP}$, assuming $n(N_{\rm IP})<-2$. After that, the velocity decreases and finally the scalar field reaches the bottom of the potential. We define the number of e-foldings $N_{1}$ and $N_{2}$ at which $n(N_{1})=n(N_{2})=-2$ and $N_{1}<N_{\rm IP}<N_{2}$ are satisfied. Using Eq. (2.12), they are obtained as

$$N_{p}=\dfrac{1}{4c}\left[\dfrac{\tanh(\tilde{\sigma}_{p})}{\cosh^{2}(\tilde{\sigma}_{p})}\left(2\ln\left[\dfrac{\sinh(\tilde{\sigma}_{i})}{\sinh(\tilde{\sigma}_{p})}\right]+\sinh^{2}(\tilde{\sigma}_{i})-\sinh^{2}(\tilde{\sigma}_{p})\right)\right]\quad(p=1,2)\ ,$$

(2.18)

where $\tilde{\sigma}_{p}=\tilde{\sigma}(N_{p})$. Therefore, the parameter $c$, characterizing the slope of the potential, partially determines a duration of e-folds for which the particle production occurs. To obtain $N_{p}=\mathcal{O}(10)$, we may choose $c=\mathcal{O}(10^{-2})$ for $\tilde{\sigma}_{i}=\mathcal{O}(1)$.

Figure 2 shows the time evolution of $\tilde{\sigma}$ and that of the corresponding $n$. We set the inflation end at $N_{\rm end}=60$. The solid and dashed curves represent the numerical and slow-roll approximate solutions (2.13) and (2.14), respectively. One can find that the slow-roll solution is well-fitted to the exact one. At an initial stage of inflation, $|n|$ is smaller than 2 and the vector field is not amplified. At a certain time, $|n|$ becomes greater than 2 and the vector field starts to be amplified. In the next section, we will show that the dark photon production occurs at an intermediate stage of inflation. As a result, the sourced power spectrum of the dark photon exhibits a peaky structure on small scales, reflecting the peaky evolution of $n$.

Let us analyze the evolution of the dark photon field kinetically coupled to the scalar field $\sigma$. Since the vector field possesses a nonzero mass, its spatial components are decomposed into transverse modes and a longitudinal mode:

$$A^{\prime}_{i}=A^{\prime T}_{i}+\partial_{i}\chi\ ,\ \qquad\partial_{i}A^{\prime T}_{i}=0\ .$$

(2.19)

Integrating out the non-dynamical component $A^{\prime}_{0}$ by

$\displaystyle A_{0}^{\prime}=\frac{-I^{2}\partial^{2}}{-I^{2}\partial^{2}+a^{2}m_{\gamma^{\prime}}^{2}}\,\partial_{\tau}\chi\;,$

(2.20)

the quadratic action of the transverse modes and that of the longitudinal mode are respectively obtained as

	$\displaystyle S_{T}$	$\displaystyle=\dfrac{1}{2}\int d\tau\,d\bm{x}\Big{[}I^{2}(\partial_{\tau}A^{\prime T}_{i}\partial_{\tau}A^{\prime T}_{i}-\partial_{i}A^{\prime T}_{j}\partial_{i}A^{\prime T}_{j})-a^{2}m_{\gamma^{\prime}}^{2}A^{\prime T}_{i}A^{\prime T}_{i}\Big{]}\ ,$		(2.21)
	$\displaystyle S_{L}$	$\displaystyle=\dfrac{1}{2}\int d\tau\,d\bm{x}\,a^{2}m_{\gamma^{\prime}}^{2}\left[\partial_{\tau}\chi\left(\dfrac{-I^{2}\partial^{2}}{-I^{2}\partial^{2}+a^{2}m_{\gamma^{\prime}}^{2}}\partial_{\tau}\chi\right)-\partial_{i}\chi\partial_{i}\chi\right]\ ,$		(2.22)

with $\partial^{2}\equiv\partial_{i}\partial_{i}$. Modulations of the kinetic couplings in these actions are crucial for the evolution of the vector field. To see this, we decompose the transverse modes ($2$ degrees of freedom) and the longitudinal one ($1$ d.o.f.) in terms of their Fourier modes,

	$\displaystyle IA^{\prime T}_{i}(\tau,\bm{x})$	$\displaystyle=\int\dfrac{d\bm{k}}{(2\pi)^{3}}\left(\hat{V}^{X}_{\bm{k}}(\tau)e^{X}_{i}(\hat{\bm{k}})+i\hat{V}^{Y}_{\bm{k}}(\tau)e^{Y}_{i}(\hat{\bm{k}})\right)e^{i\bm{k}\cdot\bm{x}}\ ,$		(2.23)
	$\displaystyle\chi(\tau,\bm{x})$	$\displaystyle=\int\dfrac{d\bm{k}}{(2\pi)^{3}}\dfrac{\hat{X}_{\bm{k}}(\tau)}{z_{k}}e^{i\bm{k}\cdot\bm{x}}\ ,\qquad z_{k}\equiv\dfrac{am_{\gamma^{\prime}}Ik}{\sqrt{I^{2}k^{2}+a^{2}m_{\gamma^{\prime}}^{2}}}\ ,$		(2.24)

where $k\equiv|\bm{k}|$, superscripts $X$ and $Y$ label the polarization states, and the orthonormal polarization vectors take the forms

$$e^{X}_{i}(\hat{\bm{k}})=(\cos\theta\cos\phi,\ \cos\theta\sin\phi,\ -\sin\theta)\ ,\quad e^{Y}_{i}(\hat{\bm{k}})=(-\sin\phi,\ \cos\phi,\ 0)\ ,$$

(2.25)

with the unit wave vector $\hat{\bm{k}}=(\sin\theta\cos\phi,\ \sin\theta\sin\phi,\ \cos\theta)$. They obey the properties: $k_{i}e^{\sigma}_{i}(\hat{\bm{k}})=0\ (\sigma=X,Y),\ e^{X}_{i}(\hat{\bm{k}})e^{X}_{i}(\hat{\bm{k}})=e^{Y}_{i}(\hat{\bm{k}})e^{Y}_{i}(\hat{\bm{k}})=1,\ e^{X}_{i}(\hat{\bm{k}})e^{Y}_{i}(\hat{\bm{k}})=0$, $e_{i}^{\sigma\,*}(\hat{\bm{k}})=e_{i}^{\sigma}(\hat{\bm{k}}),\ e_{i}^{X}(-\hat{\bm{k}})=e_{i}^{X}(\hat{\bm{k}})$, $e_{i}^{Y}(-\hat{\bm{k}})=-e_{i}^{Y}(\hat{\bm{k}})$.²²2Note that we define the parity change $\hat{k}\to-\hat{k}$ by the simultaneous operations $\theta\to\pi-\theta$ and $\phi\to\phi+\pi$. The canonical variables $\hat{V}^{\sigma}_{\bm{k}}$ and $\hat{X}_{\bm{k}}$ satisfy the reality conditions,

$$\hat{V}_{\bm{k}}^{X\,\dagger}=\hat{V}_{-\bm{k}}^{X}\;,\qquad\hat{V}_{\bm{k}}^{Y\,\dagger}=\hat{V}_{-\bm{k}}^{Y}\;,\qquad\hat{X}_{\bm{k}}^{\dagger}=\hat{X}_{-\bm{k}}\;.$$

(2.26)

Note that we put $i$ in front of the $\hat{V}^{Y}_{\bm{k}}$ term in (2.23) so that $\hat{V}^{Y}_{\bm{k}}$ respects the above property of Hermitian conjugate. In terms of these variables, the actions (2.21) and (2.22) are rewritten as

	$\displaystyle S_{T}$	$\displaystyle=\dfrac{1}{2}\sum_{\sigma=X,Y}\int d\tau\,\dfrac{d\bm{k}}{(2\pi)^{3}}\left[\partial_{\tau}\hat{V}^{\sigma\dagger}_{\bm{k}}\partial_{\tau}\hat{V}^{\sigma}_{\bm{k}}-\left(k^{2}-\dfrac{\partial_{\tau}^{2}I}{I}+\dfrac{a^{2}m_{\gamma^{\prime}}^{2}}{I^{2}}\right)\hat{V}^{\sigma\dagger}_{\bm{k}}\hat{V}^{\sigma}_{\bm{k}}\right]\ ,$		(2.27)
	$\displaystyle S_{L}$	$\displaystyle=\dfrac{1}{2}\int d\tau\,\dfrac{d\bm{k}}{(2\pi)^{3}}\left[\partial_{\tau}\hat{X}^{\dagger}_{\bm{k}}\partial_{\tau}\hat{X}_{\bm{k}}-\left(k^{2}-\dfrac{\partial_{\tau}^{2}z_{k}}{z_{k}}+\dfrac{a^{2}m_{\gamma^{\prime}}^{2}}{I^{2}}\right)\hat{X}^{\dagger}_{\bm{k}}\hat{X}_{\bm{k}}\right]\ .$		(2.28)

Then, decomposing $\hat{V}^{\sigma}_{\bm{k}}$, $\hat{X}_{\bm{k}}$ into creation/annihilation operators,

	$\displaystyle\hat{V}^{\sigma}_{\bm{k}}$	$\displaystyle=V^{\sigma}_{k}\hat{a}^{\sigma}_{\bm{k}}+V^{\sigma*}_{k}\hat{a}^{\sigma\dagger}_{-\bm{k}}\ ,\qquad[\hat{a}^{\sigma}_{\bm{k}},\hat{a}^{\sigma^{\prime}\dagger}_{-\bm{k^{\prime}}}]=(2\pi)^{3}\delta^{\sigma\sigma^{\prime}}\delta(\bm{k}+\bm{k^{\prime}})\ ,$		(2.29)
	$\displaystyle\hat{X}_{\bm{k}}$	$\displaystyle=X_{k}\hat{b}_{\bm{k}}+X^{*}_{k}\hat{b}^{\dagger}_{-\bm{k}}\ ,\qquad[\hat{b}_{\bm{k}},\hat{b}^{\dagger}_{-\bm{k^{\prime}}}]=(2\pi)^{3}\delta(\bm{k}+\bm{k^{\prime}})\ ,$		(2.30)

with the vacuum $|0\rangle$ defined by $\hat{a}_{\bm{k}}^{\sigma}|0\rangle=\hat{b}_{\bm{k}}|0\rangle=0$, we obtain the equations of motion for the mode functions as

	$\displaystyle\partial_{\tau}^{2}V_{k}+\left(k^{2}-\dfrac{\partial_{\tau}^{2}I}{I}+\dfrac{a^{2}m_{\gamma^{\prime}}^{2}}{I^{2}}\right)V_{k}=0\ ,$		(2.31)
	$\displaystyle\partial_{\tau}^{2}X_{k}+\left(k^{2}-\dfrac{\partial_{\tau}^{2}z_{k}}{z_{k}}+\dfrac{a^{2}m_{\gamma^{\prime}}^{2}}{I^{2}}\right)X_{k}=0\ ,$		(2.32)

where the index $\sigma\,(=X,Y)$ of the transverse mode functions has been omitted since the two modes obey the same equation. Historically, for analytical convenience, many of the studies on inflationary kinetic coupling models have assumed a special functional form of $I$ to make $n$ constant in time. Such a case, with a negligible mass term, is briefly summarized in Appendix A. In the case of small mass, the longitudinal mode obeys the equation of motion as in (A.11) and does not experience an exponential enhancement, in contrast to the transverse modes when $|n|>2$ is realized. Thus the contribution to the dark photon production from the longitudinal mode is always subdominant in the parameter space of our interest and shall be disregarded in the following consideration.

In order to find the solution of the transverse modes with a mild time dependence of $n(\tau)$, we assume a negligible mass of the dark photon, and then Eq. (2.31) is rewritten as

$$\partial_{\tau}^{2}V_{k}+\left(k^{2}-\dfrac{n(\tau)(n(\tau)+1)+dn/dN}{\tau^{2}}\right)V_{k}=0\ ,$$

(2.33)

with a constant $H$. Neglecting a small correction of the velocity term $dn/dN$, we obtain a differential equation of the same form as Eq. (A.1). Although in the case of time-dependent $n(\tau)$ there would be no analytically closed form of the solution to eq. (2.33) in general, there exists a useful technique called uniform approximation [54, 55, 56], which we employ here. By introducing variables

$$g(\tau)\equiv\dfrac{\nu^{2}(\tau)}{\tau^{2}}-k^{2}\ ,\qquad f(\tau)\equiv\begin{cases}\displaystyle-\left[\dfrac{3}{2}\int_{\tau}^{\tau_{*}}d\tilde{\tau}\sqrt{-g(\tilde{\tau})}\right]^{2/3}\;,&\quad\tau<\tau_{*}\;,\vspace{1mm}\\ \displaystyle\left[\dfrac{3}{2}\int_{\tau_{*}}^{\tau}d\tilde{\tau}\sqrt{g(\tilde{\tau})}\right]^{2/3}\;,&\quad\tau_{*}<\tau\;,\end{cases}$$

(2.34)

with $\nu^{2}\equiv(1+4n(n+1))/4$ and the turning point $\tau_{*}<0$ defined by $g(\tau_{*})=0$, the solution of Eq. (2.33) is well described by the following linear combination of Airy functions:

$$V^{\rm UA}_{k}\equiv A_{k}\left(\dfrac{f}{g}\right)^{1/4}\text{Ai}(f)+B_{k}\left(\dfrac{f}{g}\right)^{1/4}\text{Bi}(f)\ ,$$

(2.35)

where the superscript “UA” is to remind (the leading order of) the uniform approximation. To utilize $V_{k}^{\rm UA}$ as an approximate solution to (2.33), the coefficients $A_{k}$ and $B_{k}$ are chosen to realize the adiabatic initial condition in the sub-horizon regime and found to be

$$A_{k}=iB_{k}\ ,\qquad B_{k}=\sqrt{\dfrac{\pi}{2}}e^{i\theta}\ ,$$

(2.36)

where $\theta$ is an overall phase factor irrelevant in the following discussion. Then, using the asymptotic behavior of Airy functions, the solution on the super-horizon regime is approximately given by

$\displaystyle V^{\rm UA}_{k}$

$\displaystyle\simeq\dfrac{B_{k}}{\sqrt{\pi}g^{1/4}}\exp\left(\dfrac{2}{3}f^{3/2}\right)\simeq e^{i\theta}\dfrac{\sqrt{-k\tau}}{\sqrt{2\nu k}}\exp\left(\int^{\tau}_{\tau_{*}}d\tilde{\tau}\sqrt{\frac{\nu^{2}(\tilde{\tau})}{\tilde{\tau}^{2}}-k^{2}}\right)$

(2.37)

for $-k\tau\ll|\nu(\tau)|$. One caution to note is that this leading-order solution of the uniform approximation captures the correct spectral behavior, but its amplitude can be slightly off. Even in the case of constant $n$, or equivalently constant $\nu$, denoted by $\nu_{0}$, the comparison to the known exact solution using the Hankel function (see e.g. [30]), $V_{k}^{\rm exact}$, reveals the difference

$\displaystyle V_{k}^{\rm UA}\simeq\frac{e^{i\theta}}{\sqrt{k}}\,e^{-\nu_{0}}\left(\frac{2\nu_{0}}{-k\tau}\right)^{\nu_{0}-\frac{1}{2}}\;,\quad V_{k}^{\rm exact}\simeq\frac{e^{i\theta}}{\sqrt{k}}\,\frac{2^{\nu_{0}-1}\Gamma(\nu_{0})}{\sqrt{\pi}\left(-k\tau\right)^{\nu_{0}-\frac{1}{2}}}\;,\quad\bigg{|}\frac{V_{k}^{\rm UA}}{V_{k}^{\rm exact}}\bigg{|}\simeq\frac{\sqrt{2\pi}\,\nu_{0}^{\nu_{0}-\frac{1}{2}}}{e^{\nu_{0}}\Gamma(\nu_{0})}\simeq 0.95\;,$

(2.38)

where the last numerical value is evaluated for the case $n_{0}=-2$, or $\nu_{0}=3/2$. This is due to the error of the truncation at the leading order. The uniform approximation can accommodate the calculations of higher-order terms in an iterative manner for an arbitrary $\nu(\tau)$ [57]. For constant $\nu=\nu_{0}$, for example, the first sub-leading correction gives $V_{k}^{\rm UA(2)}=V_{k}^{\rm UA}/(12\nu_{0})$ [55], for which the difference from the exact result becomes $|(V_{k}^{\rm UA}+V_{k}^{\rm UA(2)})/V_{k}^{\rm exact}|\simeq 0.9999$ for $\nu_{0}=3/2$. Nonetheless, our leading-order solution (2.35), and its super-horizon limit (2.37), gather the correct time evolution and the spectral feature, as can be seen in Fig. 3. We thus use the leading order and admit the $\sim 5\,\%$ error in the amplitude in the following considerations.

Let us first derive the relic abundance of the dark photon DM,

$$\Omega_{\gamma^{\prime}}=\dfrac{\langle\rho_{\gamma^{\prime}}\rangle_{t=t_{0}}}{\rho(t_{0})}\ ,$$

(4.1)

where the average density of the dark photon is evaluated at the present time, denoted by subscript $0$, and $\rho(t_{0})$ is the critical density of our present universe: $\rho(t_{0})=3M_{\rm Pl}^{2}H_{0}^{2}$. To ensure the produced dark photon becomes non-relativistic some time after inflation ends at $t=t_{\rm end}$, we compare the dark photon mass $m_{\gamma^{\prime}}$ and a physical, time-dependent momentum scale with a comoving wave number at which the spectrum is peaked, $q_{\rm peak}(t)\equiv k_{\rm peak}/a(t)$. Using $k_{\rm peak}=a(t_{\rm peak})H(t_{\rm peak})$, where the time $t_{\rm peak}$ is defined by this relation and is taken to be during inflation, and defining a time duration of e-folds $\Delta N\equiv\log(a(t_{\rm end})/a(t_{\rm peak}))$, the physical momentum scale $q_{\rm peak}(t)$ is written as

$$q_{\rm peak}(t)=H(t_{\rm peak})e^{-\Delta N}\dfrac{a(t_{\rm end})}{a(t_{\rm reh})}\dfrac{a(t_{\rm reh})}{a(t)}\ ,$$

(4.2)

where $t_{\rm reh}$ denotes the time when the reheating completes. For simplicity, we assume an instantaneous reheating, $a(t_{\rm end})=a(t_{\rm reh})$, where the inflationary energy scale is related to the reheating temperature $T_{\rm reh}$ as

$$\rho(t_{\rm end})=3M_{\rm Pl}^{2}H(t_{\rm end})^{2}=\dfrac{\pi^{2}}{30}g_{*}(t_{\rm reh})T_{\rm reh}^{4}\ .$$

(4.3)

Here, $g_{*}$ is the number of relativistic degrees of freedom. The momentum scale of the dark photon at the reheating period is then evaluated as

$$q_{\rm peak}(t_{\rm reh})\sim 0.1\,\text{GeV}\left(\dfrac{H_{\rm inf}}{10^{12}\ \text{GeV}}\right)\left(\dfrac{e^{30}}{e^{\Delta N}}\right)\ ,$$

(4.4)

where $H_{\rm inf}$ denotes the inflationary Hubble scale. In our scenario, we assume a high-scale, almost de-Sitter inflation and $H(t_{\rm peak})\simeq H(t_{\rm end})\simeq H_{\rm inf}$. Therefore, the dark photon with mass $m_{\gamma^{\prime}}\ll\text{eV}$ in our interest is relativistic right after inflation ends. Then, the current energy density of the dark photon evolves from the end of inflation as

$$\langle\rho_{\gamma^{\prime}}\rangle_{t=t_{0}}=\langle\rho_{\gamma^{\prime}}\rangle_{t=t_{\rm end}}\left(\dfrac{a_{\rm reh}}{a(t_{\rm NR})}\right)^{4}\left(\dfrac{a(t_{\rm NR})}{a_{0}}\right)^{3}\ ,$$

(4.5)

where $t_{\rm NR}$ is the time when the dark photon becomes non-relativistic, determined by

$$q_{\rm peak}(t_{\rm NR})=m_{\gamma^{\prime}}\ .$$

(4.6)

By using the entropy conservation law, the present abundance of the dark photon DM is obtained as

$$\Omega_{\gamma^{\prime}}=\left(\dfrac{\pi^{2}}{90}\right)^{3/2}\dfrac{g_{S}(t_{0})g_{*}(t_{\rm reh})^{3/2}}{g_{S}(t_{\rm reh})}\dfrac{T_{\rm reh}^{3}}{M_{\rm Pl}^{3}}\dfrac{m_{\gamma^{\prime}}T_{0}^{3}}{3M_{\rm Pl}^{2}H_{0}^{2}}\dfrac{\langle\rho_{\gamma^{\prime}}\rangle_{t=t_{\rm end}}}{H_{\rm inf}^{4}}\dfrac{1}{-k_{\rm peak}\tau_{\rm end}}\ .$$

(4.7)

Therefore, the mass of the dark photon is evaluated as

$$m_{\gamma^{\prime}}\simeq 8.2\times 10^{-11}\,{\rm eV}\left(\frac{\Omega_{\gamma^{\prime}}h^{2}}{0.14}\right)\left(\dfrac{10^{7}}{\langle\rho_{\gamma^{\prime}}\rangle_{t=t_{\rm end}}/H_{\rm inf}^{4}}\right)\left(\dfrac{e^{30}}{e^{\Delta N}}\right)\left(\dfrac{10^{15}\text{GeV}}{T_{\rm reh}}\right)^{3}\ ,$$

(4.8)

where we have used $H_{0}=2.133\times 10^{-42}\,h\,\text{GeV}$, $T_{0}=2.725\,\text{K}=2.348\times 10^{-13}\,\text{GeV}$, $M_{\rm Pl}=2.435\times 10^{18}\text{GeV}$ and $g_{S}(t_{0})=3.91$, and assumed $g_{*}(t_{\rm reh})=g_{S}(t_{\rm reh})=106.75$. From this expression, we can observe that, in order for the dark photon that explains the total dark matter abundance to have a smaller mass, either its energy density $\rho_{\gamma^{\prime}}$ during inflation , $\Delta N$ or the reheating temperature $T_{\rm reh}$ must be larger. The choice of these parameters in (4.8) is rather an optimistic one, in favor of small mass, and thus it gives a rough lower bound for the dark photon mass $m_{\gamma^{\prime}}$ in our model. Yet we remark that it is sensitive to the values of $\Delta N$ and $T_{\rm reh}$, and a mass smaller than $\sim 10^{-10}\,{\rm eV}$ by a few orders of magnitude is still feasible.