\subsection*{#1}}

\newcommand{\Abst}[1]{\,#1}

\newcommand{k_{\rm B}}{k_{\rm B}}

\newcommand{m_{\rm e}}{m_{\rm e}}

\newcommand{N_{\rm e}}{N_{\rm e}}

\newcommand{n_{\rm p}}{n_{\rm p}}

\newcommand{n_{\rm b}}{n_{\rm b}}

\newcommand{T_{\rm e}}{T_{\rm e}}

\newcommand{T_{\gamma}}{T_{\gamma}}

\newcommand{\theta_{\rm e}}{\theta_{\rm e}}

\newcommand{\theta_{\gamma}}{\theta_{\gamma}}

\newcommand{m_{\rm p}}{m_{\rm p}}

\newcommand{\sigma_{\rm T}}{\sigma_{\rm T}}

\newcommand{\sigma_{\rm K-N}}{\sigma_{\rm K-N}}

\newcommand{r_{\rm c}}{r_{\rm c}}

\newcommand{n_{\rm Pl}}{n_{\rm Pl}}

\newcommand{n_{\rm BE}}{n_{\rm BE}}

\newcommand{n_{\rm W}}{n_{\rm W}}

\newcommand{\rho_{\rm Pl}}{\rho_{\rm Pl}}

\newcommand{\rho_{\rm BE}}{\rho_{\rm BE}}

\newcommand{N_{\rm Pl}}{N_{\rm Pl}}

\newcommand{N_{\rm BE}}{N_{\rm BE}}

\newcommand{\vek} [1]{\mbox{\boldmath${#1}$\unboldmath}}

\newcommand{\matr} [1]{\mbox{$\hat{#1}$}}

\newcommand{\SP}[2]{\mbox{$\vek{#1}\cdot\vek{#2}$}}

\newcommand{\partial}{\partial}

\newcommand{\pAb}[2]{\frac{\displaystyle\partial #1}{\displaystyle\partial #2}}
\newcommand{\pAbc}[3]{\left.\frac{\displaystyle\partial #1}{\displaystyle\partial #2}\right|_{#3}}

\newcommand{\PAb}[3]{\frac{\displaystyle\partial^{#3} #1}{\displaystyle\partial {#2}^{#3}}}

\newcommand{\tAb}[2]{\frac{\displaystyle d #1}{\displaystyle d #2}}

\newcommand{\TAb}[3]{\frac{\displaystyle d^{#3} #1}{\displaystyle d {#2}^{#3}}}

\newcommand{\Abl}[2]{\frac{{\rm d} #1}{{\rm d} #2}}

\newcommand{\AKc}[1]{{\rm C}_{#1}}

\newcommand{\AKs}[1]{{\rm S}_{#1}}

\newcommand{\sgn}[1]{\,{\rm sgn}\,{#1}}

\newcommand{\Mitw}[1]{\left<#1\right>}

\newcommand{\Vv}[1]{\mbox{${#1}$}}

\newcommand{\Vvek}[2]{\mbox{${#1}^{#2}$}}

\newcommand{\VvN}[1]{\mbox{${#1}_{0}$}}

\newcommand{\VvNN}[1]{\mbox{${#1}_{0}^2$}}

\newcommand{\VSP}[2]{\mbox{$\Vv{#1}\Vv{#2}$}}

\newcommand{\STV}[1]{\mbox{$\tilde{#1}$}}

\newcommand{\figref}[1]{Abb.\ref{#1}}

\newcommand{\tabref}[1]{Tab.\ref{#1}}

\newcommand{\kapref}[1]{Kap.\ref{#1}}

\newcommand{\deltaD}[1]{\mbox{$\delta( {#1} )$}}

\newcommand{\Aut}[1]{\it{#1}\rm }

\newcommand{\AutP}[1]{\it{#1}\rm }

\newcommand{\AutPet}[1]{\it{#1 et al.}\rm }

\newcommand{\AutPno}[1]{{#1}\rm }

\newcommand{\AutPnoet}[1]{{#1 et al.}\rm }

\newcommand{\Vvekm}[1]{\mbox{$\Vvek{#1}{\mu}$}}

\newcommand{\Vvekcv}[2]{\mbox{${#1}_{#2}$}}

\newcommand{\Vvekcvm}[1]{\mbox{$\Vvekcv{#1}{\mu}$}}

\newcommand{\pot}[2]{#1 \times 10^{#2}}
\newcommand{\poterr}[3]{(#1\pm#2)\times10^{#3}}

\newcommand{n_{\gamma}}{n_{\gamma}}

\newcommand{t_{\rm C}}{t_{\rm C}}

\newcommand{t_{\rm K}}{t_{\rm K}}

\newcommand{t_{\gamma \rm e}}{t_{\gamma \rm e}}

\newcommand{t_{\rm e\gamma}}{t_{\rm e\gamma}}

\newcommand{t_{\rm exp}}{t_{\rm exp}}

\newcommand{\Theta_{2.7}}{\Theta_{2.7}}
\newcommand{\theta_{z}}{\theta_{z}}
\newcommand{Y_{\rm p}}{Y_{\rm p}}

\newcommand{\Omega_{\rm b}h^2}{\Omega_{\rm b}h^2}

\newcommand{\Omega_{\rm m}h^2}{\Omega_{\rm m}h^2}

\newcommand{z_{\rm eq}}{z_{\rm eq}}

\section{Introduction}
\label{sec:intro}
In modern cosmology, our detailed understanding of the cosmic microwave background (CMB) has opened the floodgates to precision tests of the $\Lambda$CDM\xspace model and different flavours of new physics beyond. This has been made possible with precise measurements of the CMB anisotropies using space- and ground-based experiments \citep{wmap9results, ACTPol, Planck2018over}. Next generation CMB experiments furthermore promise to uncover unparalleled details of the background radiation, allowing us to probe even deeper into the underpinnings of cosmology \citep{Matsumura2014, CMBS42016, SOWP2018}.

In spite of the great successes for the $\Lambda$CDM model, many extensions have been considered.
Additions to the standard model of cosmology have included {\it modifications of the effective numbers of relativistic species} and {\it sterile neutrinos} \citep{PlanckNeutrino, Nollett2011, BattyeNeutrinos, Lesgourgues2014, Abazajian2015}, {\it dark matter annihilation} \citep{Padmanabhan2005, Galli2009, Huetsi2009, Chluba2010a}, {\it primordial magnetic fields} \citep{Sethi2005, Shaw2010PMF, Kunze2014, Chluba2015PMF, Planck2016PMF, Paoletti2019, Jedamzik2019, Jedamzik2020} and {\it variations of fundamental constants} \citep{Avelino2001,Battye2001, Galli2011,Menegoni2012,Planck2015var_alp, Hart2017, Hart2020a}. So far, no {\it significant deviation} from the $\Lambda$CDM\xspace model has been identified.
However, several tensions have been discussed \citep[see][for review]{Perivolaropoulos2021, Abdalla2022tensions}. Of these, the \emph{Hubble tension}, a discrepancy of the average expansion rate between the early and late Universe, seems to persist and grow in significance, suggesting modifications to fundamental physics might be required \citep{Bernal2016, Verde2019, Divalentino2020Intertwined}.

Indeed, several of the aforementioned new physics examples have been applied in an attempt to alleviate the notorious Hubble tension \citep[see][for a comparative study]{Schoeneberg2021}.
One viable solution is connected to the presence of an oscillating scalar field, such as those relating to \emph{ultra-light axions} (henceforth ULA), as has been considered with the consequences for the CMB anisotropies outlined in \citet{Poulin2018}. The dynamical dilution of this field in the pre-recombination era leads to an \emph{early dark energy} phenomenon, changing the expansion rate before last scattering \citep{Poulin2019}. This scalar field can be treated like an effective dark energy fluid, an approximation which has been rigorously tested against the direct field evolution \citep{Smith2020}. However, the jury is still out on whether or not this can be reconciled with large-scale structure and BAO data at lower redshifts \citep{Ivanov2020EDE,Smith2020LSS, Simon2022, Cruz2022}. Several similar dark energy theories have built on these ideas \citep[e.g.,][]{Lin2019, Alestas2020, Hill2020EDE,McDonough2021, Karwal2022, Wang2022,Kojima2022}, and most recently, even a link to possible detection of birefringence was drawn \citep{Murai2022}.

Another promising addition to the standard model of cosmology that could alleviate the Hubble tension involves the variations of fundamental constants \citep{Hart2020a}. In the interactions between matter and radiation, the main fundamental constants that garner interesting physical insights are the fine structure constant ($\alpha_{\rm EM}$) and the effective electron mass ($m_{\rm e}$).
At low redshifts ($z\lesssim 2$), the fine-structure constant has been tested with many astrophysical probes such as quasar absorption spectra \citep{Bize2003, Murphy2017, Kotus2017, Levshakov2019, Wilczynska2020VFC}, white dwarves \citep{Hu2020VFC} and black holes \citep{Hees2020VFC}.
More recently, variations in the electron-proton mass ratio have also been studied using quasar spectra \citep{Levshakov2020_me}. These works all indicate consistency with the standard value known from local lab experiements.

At higher redshifts, the changes in the 21cm line radiation that arise due to variations in the fine structure constant have been forecasted for future experiments \citep{Lopez2020}. Similarly, the changes to the light element abundances arising from Big Bang Nucleosynthesis (BBN) have been tested for fundamental constant variations \citep{Avelino2001,Coc2013,Alvey2019}.
However, the effects of varying fundamental constants (VFCs) on the CMB anisotropies indicate an interesting avenue in connection with the Hubble tension. With the most recent {\sc Planck} data, the specific dependencies of these constants during recombination have shown unique imprints during hydrogen and helium recombination \citep{Hart2017}. Specifically, the variations from the recombination epoch lead to a significant geometric degeneracy between $m_{\rm e}$ and $H_0$ which can alleviate the Hubble tension \citep{Hart2020a}.
Several reviews have been published on the motivation and various methods of detecting VFCs \citep{Uzan2003, Uzan2011, Martins2017review}, and it is important to ask if there are indeed new methods for shedding light on early VFCs.

The study of models that alter the expansion history of the universe can be carried out with many of the aforementioned cosmological probes (e.g., CMB anisotropies, weak lensing, 21cm). As explored recently, primordial $\mu$-type spectral distortions of the CMB may also provide information on the expansion rate  \citep{Lucca2020b}. However, a particular distortion that directly probes different periods of cosmological time is the \emph{cosmological recombination radiation} (CRR) from $z\simeq 1000-8000$ \citep{Sunyaev2009}. Predicted by $\Lambda$CDM\xspace, this distortion arises as the CMB photon field departs from thermal equilibrium due to the transitions within hydrogen and helium atoms as well as the continuum during recombination \citep{Zeldovich68,Peebles68}. This manifests as a unique spectral signal in the CMB spectrum \citep{Dubrovich1975, RybickiDell94}, which can now be accurately computed using {\tt CosmoSpec} \citep{Chluba2016CosmoSpec}. Given the superposition of hydrogen and helium lines, the spectral changes caused by variations in cosmological parameters can be constrained with futuristic spectrometers (e.g., \emph{PRISM}, \emph{Voyage 2050}\xspace) \citep{PRISM2013WP,PRISM2013WPII,Vince2015, Mayuri2015, Chluba2019Voyage,Hart2020c}.
Similarly, new physics can affect the dynamics of the recombination process and thus leave unique imprints in the CRR \citep{Jose2008, Chluba2008T0, Chluba2008c}.

In this paper, we will outline the ways that the CRR can probe the effects from early dark energy theories and variations of the fundamental constants $\alpha_{\rm EM}$ and $m_{\rm e}$. In Section~\ref{sec:ede}, we introduce the approach used in the previous ULA constraints papers \citep{Poulin2018, Poulin2019} and revisit the main effects on the background expansion history. We then study how early dark energy can affect the ionisation history and consequently, impact the recombination lines in unique ways.
We briefly discuss the detectability of these variations for different models using rudimentary signal-noise predictions and then investigate more complete parameter correlations with a Fisher matrix analysis.
Our estimates are meant to give a first rough feeling about the observability of these effects; however, a rigorous analysis in combination with CMB anisotropy constraints is left for a future investigation.

In Section~\ref{sec:vfc}, we show the differences in the CRR caused by VFCs. We explain how these variations can be related to the features discussed in \citet{Hart2017} with a particular emphasis on the modifications to the recombination process caused by these changes.
We provide a comparative study indicating the impact of including CMB anisotropy results with future spectrometers. This leads to discussion on the possible solutions to the Hubble tension involving $m_{\rm e}$ \citep{Hart2020a,Hart2021}.

\section{Early dark energy}\label{sec:ede}
The equations of motion for the ULA can be reconstructed using the effective fluid approximation as has been validated in comparison to the full scalar field evolution \citep{Smith2020}. The field dynamics lead to an evolving energy density,
\begin{equation}
    \Omega_\phi\left(z\right) = \frac{2\Omega_\phi\left(z_{\rm\,c}\right)}{1+\left[\left(1+z_{\rm\,c}\right)/\left(1+z\right)\right]^{\,3\left(1+w_n\right)}}\;,
\end{equation}
with an equation of state,
\begin{equation}
    1+w_\phi(z) = \frac{1+w_n}{1+\left[\left(1+z\right)/\left(1+z_{\rm\,c}\right)\right]^{\,3\left(1+w_n\right)}}\;.
\end{equation}
Here $z_{\rm\,c}$ signifies the redshift when the field becomes dynamical and $n$, the order of the oscillating potential for a ULA field \cite[][for more details]{Poulin2018}, determines the dilution rate of the energy density according to $w_{n} = \left(n-1\right)/\left(n+1\right)$.

The energy density for several potential orders, $n$, are illustrated in Fig.~\ref{fig:ede}. To quantify the amplitude of the early dark energy density, we will use the parametrisation $f_{\rm\, EDE} = \Omega_\phi(z_{\rm\,c})/\Omega_{\rm tot}(z_{\rm\,c})$, following previous ULA papers. At early times ($z\gg z_{\rm\,c}$), the early dark energy behaves like a cosmological constant due to the \emph{Hubble friction} term of the evolving field. As the dark energy fluid becomes dynamical, it decays according to $\Omega_{\phi}\propto\left(1+z\right)^{3\left(1+w_n\right)}$. For $n=2,3$ this means the early dark energy density drops off like radiation or an ultra-relativistic species respectively. In the extreme case that $n\rightarrow\infty$, the field energy density rapidly decays $\propto a^{-6}$. This corresponds to the field energy being totally dominated by kinetic term, similar to a scalar field term\footnote{This was explored with detailed changes to the background cosmology in \citet{Karwal2016}.}.

\begin{figure}
    \centering
    \includegraphics[width=.98\columnwidth]{figures/ede_omega_n_2.pdf}
    \\
    \includegraphics[width=.98\columnwidth]{figures/ede_omega_n_3.pdf}
    \\
    \includegraphics[width=.98\columnwidth]{figures/ede_omega_n_500.pdf}
    \\
    \caption{The energy density for EDE under the fluid-approximation for $f_{\rm\, EDE}=0.3$. Curves for 3 different critical redshifts $z_{\rm\,c}$ are shown for $n=2$ ({\it top}), $n=3$ ({\it middle}) and $n=\infty$ ({\it bottom}) panels, respectively. The matter and relativistic density evolution are shown for comparison as well as the matter-radiation equality epoch.}
    \label{fig:ede}
\end{figure}

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{figures/hz.pdf}
    \\
    \caption{Examples of the evolution of the Hubble parameter $H(z)$ as a function of the different indices $n=\{2,3,\infty\}$ at $z_{\rm\,c}=3000$. The parameter growth with comparison against $\Lambda$CDM\xspace (gold, dashed) is shown in the top panel, while the difference of the models, $\Delta H$, are shown in the bottom panel. Overlaid are the rough neutral hydrogen, neutral and ionised helium recombination eras \citep[see][for more details]{Sunyaev2009}. The Hubble parameter here is measured in ${\rm km \,s^{-1} Mpc}$, as usual.}
    \label{fig:hz}
\end{figure}

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{figures/Xe_zc_1500.pdf}
    \includegraphics[width=\columnwidth]{figures/Xe_zc_6000.pdf}
    \\
    \caption{The ionisation history, $x_{\rm e}$ for added early dark energy contributions§ with $z_{\rm\,c} = 1500$ (\emph{top}) and $z_{\rm\,c} = 6000$ (\emph{bottom}) compared to the $\Lambda$CDM\xspace case (gold/dashed). The different order of potentials $n\in\{2,3,500\}$ is shown by the colour scheme \{black, purple, orange\}, respectively. The critical redshift, $z=z_{\rm\,c}$ is also shown by a vertical black line.}
    \label{fig:ede_xe}
\end{figure}

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{figures/DXe_zc_1500.pdf}
    \includegraphics[width=\columnwidth]{figures/DXe_zc_6000.pdf}
    \\
    \caption{Variations in the ionization history arising from $f_{\rm\, EDE} = 0.3$ and dynamical redshifts $z_{\rm\,c}=1500$  (\emph{top}) and $z_{\rm\,c}=6000$  (\emph{bottom}). For both cases, the different dilution rates for $n\in\{2,\,3,\,500\}$ are shown. Rough guidelines for the different recombination eras are also included.}
    \label{fig:ede_dxe}
\end{figure}

Here will we discuss the impact of an early dark energy contribution to the recombination lines so it is important to emphasise that we will only consider the effects on the background physics. In this case we do not consider the full changes to perturbations carried out in previous analyses, rather we consider the changes to the recombination process due to the altered expansion history shown in Fig.~\ref{fig:ede}. Furthermore, we ensure that flatness is conserved for these models (i.e., $\Omega_\Lambda\rightarrow \Omega_\Lambda + \Omega_{\phi,0}$).
\changeJ{This condition should not directly affect the recombination lines, as $\Omega_\Lambda$ only becomes important at late times.}
The recombination calculations are carried out using {\tt CosmoRec} \citep{Chluba2010b} with the {\tt CosmoSpec} module \citep{Chluba2016CosmoSpec} to compute the CRR. This allows us to cleanly propagate all the effects on the hydrogen and helium recombination dynamics.

Current constraints on $f_{\rm\, EDE}$ derived from CMB and large-scale structure measurements imply $f_{\rm\, EDE}\simeq 0.05-0.1$ for critical redshifts $z_{\rm\,c}\simeq 3000-4000$ [and fixed $n\simeq 3$] \citep{Hill2020EDE, Simon2022, Cruz2022}. To better illustrate the effects on the CRR, we will use larger values of $f_{\rm\, EDE}\simeq 0.3-0.5$ and also widen the range of critical redshifts explored for varying values of $n$ as stated. The value of $f_{\rm\, EDE}$ will mostly lead to an overall rescaling of the corresponding signals, while both $z_{\rm\,c}$ and $n$ affect the shape of the distortion responses, as we illustrate below.

\subsection{Ionization history effects}\label{sec:ede_xe}
The modifications to the background dynamics propagate to the ionisation history. The changes to the Hubble rate $H(z)$ alter the total energy density of the Universe as a function of time and affect the mapping to redshift (see Fig.~\ref{fig:hz}). Consequently, recombination is delayed as can be seen from the ionisation history ($x_{\rm e}$) variations shown in Fig.~\ref{fig:ede_xe}. These are compared to the case for standard $\Lambda$CDM\xspace (gold/dashed). For all dilution rates, the variations in $x_{\rm e}$ are larger for a lower critical redshift, $z_{\rm\,c} = 1500$. Due to the smaller relative helium fractions \citep[and helium feedback processes at $z\simeq 2000$,][]{Chluba2009c, Chluba2012HeRec}, the relative changes in the ionisation history are more pronounced during hydrogen recombination, which is expected to lead to greater variations in the Thomson visibility function (since the Thomson visibility is larger around the last scattering epoch at $z\simeq 1100$).

However the residuals of the ionisation history shown in Fig.~\ref{fig:ede_dxe} reveal the bigger picture for the earlier epochs of recombination. For both cases ($z_{\rm\,c} = 1500$ and $z_{\rm\,c} = 6000$), this corresponds to a positive residual in the ionisation history ($\Delta \xe /\xe >0$). When the critical redshift (the time when the field becomes dynamical and the energy density dilutes) is shifted from the hydrogen recombination era ($z_{\rm\,c} = 1500$) to the doubly-ionised helium ($\text{HeIII}\rightarrow\text{HeII}$) recombination era ($z_{\rm\,c} = 6000$), modifications in the epochs of helium recombination become visible. Specifically, a noticeable variation around $z\simeq 6000$ arises in the free electron fraction, while the changes during hydrogen and neutral helium recombination show an increased sensitivity to the value of $n$. This highlights that the CRR can in principle be used as a probe of EDE. In particular for models with $z_{\rm\,c}\gtrsim 5000$ this could nicely complement probes based on the CMB anisotropies alone, which already tightly constrain scenarios with $z_{\rm\,c}\lesssim3000$ \citep{Simon2022}.

\subsection{Propagating changes to spectral distortions}\label{sec:ede_dist}
The changes to the ionization history that arise from early dark energy model extensions can be propagated into the deviations of the CRR using {\tt CosmoSpec}. In this section, we isolate some of the key features that have been modified in the CRR by considering an EDE species added into the expansion rate.

\subsubsection {Hydrogen recombination lines}
\label{sec:hEDE}

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{figures/hydrogen1000.pdf}
    \\
    \caption{Hydrogen recombination lines for EDE models fixed with a dynamical time of $z_{\rm\,c} = 1000$ for a variety of dilutions [$n=\{2,3,\infty\}$]. These are compared against the $\Lambda$CDM\xspace benchmark (\emph{black}) over the frequency of interest $\nu \simeq 30-2000$ GHz.}
    \label{fig:hydrogenEDEzc1000}
\end{figure}

In Fig.~\ref{fig:hydrogenEDEzc1000}, the hydrogen recombination lines are shown for $\Lambda$CDM\xspace (\emph{black}) and against the added EDE fluid with an `late' dynamical time of $z_{\rm\,c}=1000$. For illustration, we vary the dilution rate ($\leftrightarrow n$) at fixed $f_{\rm\, EDE}=0.5$.
There is very little movement in the line positions; however, the amplitudes of most spectral features is suppressed, including the Balmer-$\alpha$ line ($\nu\simeq120$ GHz) and Paschen-$\alpha$ line ($\nu\simeq350$ GHz). In addition, the lines are broadened due to the more extended duration of the recombination process.
At $\nu \gtrsim 1800\,{\rm GHz}$, we can furthermore see an increased blue-wing of the HI Lyman-$\alpha$ line while the level of the two-photon continuum ($\nu \simeq 800-1500\,{\rm GHz}$) is reduced, indicating a delay in the recombination process.
Since the hydrogen line emission process occurs at $z\simeq 1400$, for the chosen example the variations in the lines appear invariant with the speed of dilution (affected by $n$) for this model.

In Fig.~\ref{fig:hydrogenEDEN2peaks}, we focus our attention on the Balmer-$\alpha$ and Paschen-$\alpha$ lines.
Here, the ULA treatment of EDE is calculated with a variety of dynamical times and a radiation-like dilution ($n=2$).
For $n=2$, we can see that even the changes of the CRR are mostly independent of the chosen dynamical time, though there is a very small ($\lesssim 5\%$) effect for $z_{\rm\,c}=5000$, which is marginally closer to the $\Lambda$CDM\xspace case. %
The main conclusion from both Figs.~\ref{fig:hydrogenEDEzc1000} and~\ref{fig:hydrogenEDEN2peaks} is that the variations for slower dynamical EDE models closer to the surface of last scattering are qualitatively indifferent for hydrogen recombination lines.

We present the contrary to this in Fig.~\ref{fig:hydrogenN500weak}, where we focus on a pre-recombination dynamical time, $z_{\rm\,c}=6000$. For the different dilution models, the changes compared to $\Lambda$CDM\xspace are heavily hindered, with the kinetic ULA example ($n\rightarrow\infty$; \emph{orange}) being almost identical to the $\Lambda$CDM\xspace case. For a dynamical field that dilutes exponentially fast, the variations do not seed in the hydrogen lines at all. Since the ULA contributions to the energy density dilute away very quickly, the hydrogen recombination process is shielded from these modifications to the expansion rate. In comparison, and with reference to the species dilution curves in Fig.~\ref{fig:ede}; the denser EDE models such as radiation-like ($n=2$) and ultra-relativistic ($n=3$) tail off more slowly, since the dilution of the field is much harder. Consequently, the net impact on the expansion rate is larger and the interplay across the hydrogen recombination lines for an early dynamical time such as $z_{\rm\,c}=5000$ becomes more apparent, as presented in Fig.~\ref{fig:hydrogenN500weak}. This shows that the CRR is sensitive to the dilution rate of EDE models with $z_{\rm\,c}\gtrsim 3000$, as already anticipated from Fig.~\ref{fig:ede_dxe}.

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{figures/hydrogenN2_peaks.pdf}
    \\
    \caption{Example of effects on the Balmer-$\alpha$ (H$\alpha$) and Paschen-$\alpha$ lines caused by EDE modifications. Here the ULA field dilutes like cold-dark matter however is shown for a variety of dynamical redshifts, all showing suppression of the emission peaks.}
    \label{fig:hydrogenEDEN2peaks}
\end{figure}

\begin{figure}
    \centering
    \includegraphics[width=\columnwidth]{figures/hI_n500_weak.pdf}
    \\
    \caption{Hydrogen lines over the same frequency range as shown in Fig.~\ref{fig:hydrogenEDEzc1000}, however focussing on a dynamical time of $z_{\rm\,c} = 5000$. The variations arising from the EDE models are much smaller in comparison to the model described by the aforementioned figure.}
    \label{fig:hydrogenN500weak}
\end{figure}

\subsubsection{Helium recombination lines}\label{sec:heEDE}

As emphasised in Sec.~\ref{sec:ede_xe}, the helium recombination lines will be more susceptible to changes in the expansion rate encroaching on earlier epochs ($z>2000$). However given the lower impact of helium recombination, both within the lines and the optical depth of CMB photons arising during the decoupling era, we expect this effect to be unique but smaller.
One of the most distinctive features in helium recombination is the backwash of photons arising from feedback processes between He-I and He-II. In Fig.~\ref{fig:heIabsorption}, the absorption trough for helium where $\nu\simeq270$ GHz is shown. By adding the diluting field that emulates EDE, the trough shifts to higher frequencies and dips to a weaker signal ($\Delta I_\nu\simeq0.002$ Jy/sr). Furthermore, the profile defining this `absorption' is much sharper for an EDE model. The earlier dynamical time starts to impact singly-ionized helium recombination (He-II$\rightarrow$He-I) as the Hubble flow accelerates before this epoch. Consequently, the energy rates involved are naturally dampened by the increased expansion making the effects from the absorption sharper than in the $\Lambda$CDM\xspace case.

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/heI_absorption_n2n3.pdf}
    \caption{The absorption trough in the \changeJ{He-I radiation} at $\nu\simeq260$ GHz. The matter and radiation-like dilution models are shown for a fixed dynamical time associated from before He-I recombination ($z_{\rm\,c} = 3000$).}
    \label{fig:heIabsorption}
\end{figure}

A wider example of this absorption feature can be seen if we look at the full \changeJ{He-I} recombination spectra between the frequencies that are relevant for future space missions\footnote{In the  \emph{PIXIE} \citep{Kogut2011PIXIE} and \emph{Voyage 2050}\xspace \citep{Chluba2019Voyage} mission concepts, the frequency bands were defined for $\nu_{\rm min} = 30$ GHz and $\nu_{\rm max} = 3000$ GHz.}. This is highlighted in Fig.~\ref{fig:heIbb}, where we have also shown the impact of the full distortion against the case where bound-bound transitions are considered only (\emph{faded, dashed}). Removing the logarithmic $y$-axis, the rich structure of the helium emission and absorption features are clearer to see; however, the impact of the EDE changes does not reveal added information when we omit the bound-free absorption components in the distortions (free-free typically affects much lower frequencies). The contrast between amplification of the peaks and suppression of the troughs in the wider frequency information for helium also suggests that the expansion's impact on the net energy transitions is not trivial, due to the interlinked nature of He-I recombination to both the later hydrogen recombination epoch and the earlier doubly-ionised helium region.

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/heI_bb_n23.pdf}
    \caption{Helium-I recombination lines from Fig.~\ref{fig:heIabsorption} with wider frequencies to show the impacts on the $\Lambda$CDM\xspace case and the bound-bound (BB) spectra when a similar EDE model is added (e.g., $z_{\rm\,c}=3000$,$n=\{2,3\}$).}
    \label{fig:heIbb}
\end{figure}

Finally, we focus solely on the doubly-ionised helium ($\text{HeIII}\rightarrow\text{HeII}$) region, which occurs significantly earlier ($z\simeq 6000$). In Fig.~\ref{fig:heII_shifts}, the absorption between the two helium species is shown ($300-500$ GHz) for $\Lambda$CDM\xspace against two dilution examples. The line at $\nu\simeq380$ GHz represents a pivot for the added EDE models. At lower frequencies, this feature in helium is suppressed whereas it is amplified for higher frequencies. The full recombination spectrum for doubly-ionised helium is shown for context in Fig.~\ref{fig:heII_all}, however this was shown for $z_{\rm\,c}=6000$ since it is more directly related to the doubly helium recombination era. The variations according to early dark energy all show distinct changes to the $\Lambda$CDM\xspace CRR; however, the different dilution models do not seem to create appreciable changes in the spectra for the chosen value of $z_{\rm\,c}$. Whether the field dilutes as radiation or ultra-radiation, the spectra look very similar.
\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/heII_shifts.pdf}
    \caption{The effects on the He-II distortion arising from an EDE model with very early dynamical time ($z_{\rm\,c}=10000$) for $300$ GHz $<\nu<500$ GHz. The vertical line at $\nu=378$ GHz acts as a pivot for the variations when EDE is introduced. This leads to amplification of the signal at $\nu>378$ GHz and suppression at $\nu<378$ GHz. There is subtle smaller pivoting in the accompanying helium resonance at $\nu\simeq320$ GHz.}
    \label{fig:heII_shifts}
\end{figure}

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/heII_all.pdf}
    \caption{The full distortion arising from doubly-ionised helium recombination between the proposed \emph{PIXIE} bandwidths ($30$ GHz $<\nu<3000$ GHz). Here the different EDE models are shown with various dilution speeds as we have previously discussed ($n=\{2,3,500\}$) for a critical redshift $z_{\rm\,c}=6000$.}
    \label{fig:heII_all}
\end{figure}
These changes arise because, as with previous cases, the dynamic changes for these models are harder to distinguish for earlier redshifts (i.e., during He-III$\rightarrow$He-II recombination). However, the knee-pivot discussed in Fig.~\ref{fig:heII_shifts} at $\nu\simeq400$ GHz is an isolated large change which can be associated with continuum processes. These continuum features \citep[discussed in more detail in][]{Switzer2007I, Jose2008} can be amplified if the electrons recombine into the bound states. This is more emphasised with states at higher energies, that are simultaneously closer to the ground state. There is also interplay with the fine structure lines in the helium atom \citep[see][for more details]{Chluba2012HeRec}.
\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/h_all.pdf}
    \caption{Comparison for an EDE model that has different slopes for $z_{\rm\,c}=5000$ across the different atomic species. Hydrogen (\emph{top}) is compared against helium I and II (\emph{middle,bottom}) respectively. For better clarity, this has been plotted with $f_{\rm\, EDE} = 0.8$.}
    \label{fig:h_all}
\end{figure}
\subsubsection{Combined effect from all the atomic species at $z_{\rm\,c}=5000$}\label{sec:total5000}
To illustrate the total impact of various models when the early dark energy starts diluting right in the middle of the two helium recombination eras, we show the effects on the hydrogen and two helium energy spectra in Fig.~\ref{fig:h_all}. For the purposes of this figure, we have increased the relative amount of EDE to $f_{\rm\, EDE}=0.8$. While the hydrogen lines are showing signs of spectral deviations according to the dilution slope, $n$; the helium atomic variants are showing this less prominently. There are difference tracers in the lower frequency ($\nu\simeq200$ GHz), particularly around the He-I absorption trough; however, the effects are less apparent in the helium plateau at $\nu\simeq800$ GHz and there is negligible spectral variance in the doubly-ionised helium lines.
This is anticipated, since the lever arm is shortest for the He-II spectrum and biggest for the H-I radiation, rendering the $n$-sensitivity largest for H-I.
These variations are hinting that the EDE mechanism can begin to act as a pivot: highlighting not only changes in different magnitudes across the epochs, but more importantly that the impact of the sloping parameter $n$ is non-negligible as you tune your critical redshift to earlier times.

\subsection{Basic forecast for future spectrometers}\label{sec:forecastEDE}
In this section, we use a basic signal-noise ratio (SNR) test to validate the most responsive EDE models affecting the recombination lines. From modelling the various combinations of the dilution speed and critical redshift for the EDE theories (i.e., $n$ and $z_{\rm\,c}$), we can quantify the models with the largest net change to the $\Lambda$CDM\xspace reference using a basic signal-to-noise (diagonal Fisher matrix) testing case \citep[statistical use cases explained in][]{Sellentin2016,Bhandari2021}. Since the early dark energy fraction parameter ($f_{\rm\, EDE}$) is the amplitude of the signal and is fairly model-independent (larger values of $f_{\rm\, EDE}$ imply larger responses) then this is the parameter for which we model an effective signal-to-noise ratio to identify exceptional models.

For this simple case study, we have isolated the frequency bins using the proposals for channels outlined for \emph{SuperPIXIE}, where we consider 3 instruments: low-frequency ($10$ GHz $<\nu<40$ GHz; $\Delta\nu = 2.4$ GHz), mid-frequency ($20$ GHz $<\nu<600$ GHz; $\Delta\nu = 19.2$ GHz) and high-frequency ($400$ GHz $<\nu<6000$ GHz; $\Delta\nu = 57.6$ GHz). The full details of this spectrometer setup can be found in many forecast proposal papers \citep{Kogut2011PIXIE,Abitbol2017, Chluba2019Voyage, Hart2020c}. The signal-to-noise is isolated in each bin and then the summation RMS value is taken for every given model. After evaluating each of the known dilution models [$n = \{2,3,\infty\}$] for a wide range of critical redshifts [$z_{\rm\,c} = \{1000,1500,2000,3000,5000,6000,10000,20000\}$], the full response of each model can be seen in Fig.~\ref{fig:edeMatrix}. Here the specifications for a spectrometer following the \emph{Voyage 2050}\xspace outline \citep{Hart2020c} have been used to calculate an effective SNR. Note that this process omits correlated parameters but we can start to implement a basic Fisher matrix analysis and see how the parameter dependencies vary for different EDE models.

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/edeVoyageSNR.pdf}
    \caption{Signal-to-noise forecasts for \emph{Voyage 2050}\xspace for various EDE models defined by $n$ and $z_{\rm\,c}$. The largest SNR values corresponds to the brightest colours in the colorbar. Similar SNR distributions can be found for \emph{SuperPIXIE} and \emph{Voyage 2050+}\xspace by multiplying the SNR by the appropriate change in sensitivity. For this particular study, $f_{\rm\, EDE} = 0.05$. }
    \label{fig:edeMatrix}
\end{figure}
The SNR values for these models have been visually summarised in Fig.~\ref{fig:edeMatrix}. The top five models have been quantitatively presented as a list in Table~\ref{tab:snrEDE}. As initially predicted from the large hydrogen variations in Sec.~\ref{sec:ede_dist}, the strongest models for a \emph{Voyage 2050}\xspace like mission have a slope consistent with a radiation-like dilution ($n=2$). Constraining the EDE model requires a \emph{`sweet spot'} solution where the field dilutes quickly enough that the net change in dynamics across recombination is large, while not requiring such a high initial value that may affect the initial conditions. These initial conditions will ultimately affect the background as recombination lines begin to form and the non-thermal processes between electrons and photons emerge. More specifically, the variations favour a dilution that begins in the HeI recombination era ($z\simeq2000,3000$) since this provides the largest amount of variation within the recombination radiation for the smallest additional factor of EDE. The most constraining detectable model according to this comparison is the model is a radiation-like dilution with $z_{\rm\,c} = 2000$ with an SNR of $2.02$ for \emph{Voyage 2050}\xspace, as shown in Table~\ref{tab:snrEDE} ($M_A$). Note that the results shown for the \emph{SuperPIXIE}\xspace and \emph{Voyage 2050+}\xspace configurations are multiplied by constant factors of $0.2$ and $10$ respectively, as their noise curves are defined.
Surprisingly, for radiation-like dilution, a redshift switch of $z_{\rm\,c}=5000$ has shown similar SNR in Table~\ref{tab:snrEDE} (Model $M_E$), suggesting that some of the pre-HeI recombination models can survive. Large variations from models like these will still appear during hydrogen recombination; however, as shown in Fig.~\ref{fig:edeMatrix}, the sharper dilution models quickly deteriorate to low SNRs once the EDE phase-transition redshift is pushed to earlier times.

In summary, we can see that for radiation-like dilution ($n=2$), a wide range of critical redshifts may be probed, even reaching deep into the primordial universe at $z\gtrsim 10^4$. For steeper dilution, sensitivity to $z_{\rm\,c}\gtrsim 6000$ quickly drops, and lower values of $z_{\rm\,c}$ give preferable responses.

\subsection{Using Fisher matrix analysis to estimate the errors}
\label{sec:fisherEDE}
From the leading order estimates on the strength of the EDE amplitude parameters, we can compute more realistic responses between the $\Lambda$CDM\xspace parameters and $f_{\rm\, EDE}$ for various models.
In this section, we obtain Fisher forecasts for the 3 highest SNR models found in Table.~\ref{tab:snrEDE}. We discuss the correlations from model to model and refer to the main features of the spectral distortion variations highlighted in Sec.~\ref{sec:ede_dist}.\footnote{It is also important to note that we will revisit this methodology in Sec.~\ref{sec:vfc} where we apply a simpler version of the Fisher forecast to constraining fundamental constants with the recombination lines.}
One of the simplest statistical measures that we can use to test for parameter correlations is the Fisher information matrix that defines the covariances at the peak of the likelihood. The Fisher matrix is defined in many ways however for spectral distortions, the matrix is defined by,
\begin{equation}
    F_{ij} = \sum _{\nu \nu'}\frac{\partial \Delta I_{\nu}}{\partial p_i}\; {\bf{\Sigma}}^{-1}_{\nu \nu'} \; \frac{\partial \Delta I_{\nu'}}{\partial p_j},
\end{equation}
where $\{p_{i},p_{j}\}$ are the parameters in the correlation study that correspond to the Fisher matrix element $F_{ij}$. Here $\Sigma_{\nu\nu'}$ is the covariance matrix for the $\nu\times\nu'$ frequency bands. In this analysis, our covariance matrix is going to be made from the total signal coming from the fiducial $\Delta I_\nu$ spectra and the noise spectra discussed in Sec.~\ref{sec:forecastEDE}.
Note that this is the same formalism of the Fisher matrix that was used in the previous paper constraining cosmological parameters with the recombination lines \citep{Hart2020c}. In Sec.~\ref{sec:resultsF}, we will present the covariances for the 3 most constrainable models (shown in \,Table.~\ref{tab:snrEDE}).
\changeL{The EDE amplitude $f_{\rm\, EDE}$ is bound by a hard-prior where $f_{\rm\, EDE}>0$. When sampling the Fisher matrices required for Fig.~\ref{fig:edeContours2}, we have not considered the physical limits such as those imposed by scalar parameters. Hard priors will be included in evaluating the posterior, designed for a full MCMC analysis. The methodology assumes perfectly Gaussian likelihoods and does not appropriately evaluate hard-priors for non-negative parameters (such as $f_{\rm\, EDE}$ and $N_{\rm eff}$). The errors calculated in this forecast are designed to be order-magnitude estimates for comparison and assessment against other probes.}

\begin{table}
    \centering
    \begin{tabular}{c|c|c|c|c|c}
        \hline\hline
        Model & $n$ & $z_{\rm\,c}$ & SNR  & SNR  & SNR \\
        & & & (\emph{SuperPIXIE}) & (\emph{V2050}) & (\emph{V2050+})\\
        \hline
        $M_A$ & $2$ & $2000$ & $0.40$ & $2.02$ & $20.20$ \\
        $M_B$ & $2$ & $3000$ & $0.40$ & $1.98$ & $19.78$ \\
        $M_C$ & $2$ & $1500$ & $0.35$ & $1.76$ & $17.59$ \\
        $M_D$ & $3$ & $2000$ & $0.34$ & $1.72$ & $17.22$ \\
        $M_E$ & $2$ & $5000$ & $0.34$ & $1.70$ & $17.06$ \\
        \hline
    \end{tabular}
    \caption{The signal-to-noise ratios for the 5 most constrained models of EDE using the recombination lines (referred to as $M_i$), for 3 different model configurations. Note that \emph{Voyage 2050+}\xspace (\emph{V2050+}) is 10 times more sensitive than \emph{Voyage 2050}\xspace (\emph{V2050}). Here the dilution index $n$ and the dynamical redshift $z_{\rm\,c}$ of each of these models are given along with the SNR (signal-noise ratio) value. \changeJ{We used $f_{\rm EDE} = 0.05$ as reference value.}}
    \label{tab:snrEDE}
\end{table}

\subsubsection{Contour results for the Fisher matrix}\label{sec:resultsF}
Results from the Fisher matrix can tell us about the underlying parameter degeneracies in this EDE model. Using the Fisher matrix with a generation of random Gaussian samples, one can visualise the degeneracies using the same contours as in MCMC\footnote{The generation of Gaussian samples and subsequent plots were done with the {\tt GetDist} package \citep{GetDist}.}. The results for the 3 most promising models described in Sec.~\ref{sec:forecastEDE} are shown in Fig.~\ref{fig:edeContours2}. Note that the EDE amplitude parameter has been set to $f_{\rm EDE} = 0.0$ as the fiducial pivot point for the calculation. The stability of the step size for the Fisher matrix treatment was akin to the previous paper, as well as the diagonal Gaussian likelihood approach to the Fisher matrix \citep{Hart2020c}.

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/edeContours2}
    \\
    \caption{Estimated probability contours with a {\it Voyage 2050+} style mission for the EDE amplitude at $z=z_{\rm\,c}$, $f_{\rm EDE}$. The fiducial values for the standard parameters arise from \emph{Planck}\xspace and the fiducial value for the EDE amplitude is set to 0. The different colours refer to the three most promising models, whereas the black line indicates $\Lambda$CDM\xspace.}
    \label{fig:edeContours2}
    \vspace{-3mm}
\end{figure}

All three models in Fig.~\ref{fig:edeContours2} have some degeneracies with $\omega_{\rm c}$ and $N_{\rm eff}$. Interestingly, the higher values of $z_{\rm\,c}$ in this particular configuration get gradually more constrained. The contours for the $z_{\rm\,c}=3000$ model have the smallest contours where the $1\sigma$ limit $\sigma_f\approx0.3$ assuming {\it Voyage 2050+}.
The marginalised errors attained with \emph{Planck}\xspace for similar models are $\simeq 4$ times smaller \citep{Poulin2018, Hill2020EDE}. Not only do the contours for $n=2$ models show signs of deformation which can unknowingly bloat contours in simplistic analyses such as the Fisher presented here; it was shown in our previous paper that the expansion rate quantities would require a further $\simeq 50$ times greater sensitivity than \emph{Voyage 2050}\xspace to get the desired parameter constraints that would complement \emph{Planck}\xspace. For these particular models, with a $50 \times$\emph{Voyage 2050}\xspace sensitivity, we could hope to constrain $f_{\rm EDE}\lsim0.05$ using only a CMB spectrometer. However this also neglects the involvement of foregrounds, albeit their broader spectral shape mitigates how badly they affect cosmological parameter constraints \citep[more details in][]{Hart2020c}.
Comparatively, for a \emph{Voyage 2050+}\xspace style mission, the model that begins to dilute during the peak of hydrogen recombination ($z_{\rm\,c} = 1500$) is only at $f_{\rm EDE} \lesssim 1$. This becomes the case as the contributions correlate directly with a matter-dilution effect on the CRR which primarily dominates hydrogen recombination: the fractional cold dark matter density today $\omega_c$.

\changeJ{We note that at this stage we had no access to MCMC chains relating to the analysis of EDE models with {\sc Planck}. We therefore could not explore robustly how the addition of external data sets could help breaking parameter degeneracies. For the standard CRR analysis, we saw this to yield significant gains for $N_{\rm eff}$ and $Y_{\rm p}$ \citep{Hart2020c}. Below we will illustrate how in the case of VFCs the addition of {\sc Planck} priors indeed significantly improves matters. In a similar vein, we expect significant gains for EDE models when combining CRR measurements with external data set; however, a more detailed assessment is beyond the scope of this work.}

\vspace{-3mm}
\section{Varying fundamental constants}\label{sec:vfc}
One of the promising extensions to the standard $\Lambda$CDM\xspace paradigm is the addition of variations to fundamental constants as we have alluded to in Sec.~\ref{sec:intro}. Specifically, recombination physics relies on the couplings between charged particles (electrons) and photons; therefore, the main constants that will affect this interaction are the fine structure constant ($\alpha_{\rm EM}$) and the effective electron mass ($m_{\rm e}$). In this section, we will recap the effect to the  recombination process arising from VFCs and succinctly discuss the updates of {\tt CosmoSpec} since the last analysis of the CRR. These updates proved vital for the accurate calculation of fine structure variations within the recombination epoch. Finally, we present the resulting spectral distortion changes from these parameters and use those in a Fisher forecast akin to Sec.~\ref{sec:fisherEDE} to get some zeroth level detectability estimates.

\subsection{Recap of VFCs in recombination}
\label{sec:recap_vfc}
Fundamental constants have a wide impact on various parts of cosmology, however the physics of recombination can be directly traced by variations in two constants: the fine structure constant $\alpha_{\rm EM}$ and the effective electron mass $m_{\rm e}$. These both affect many aspects of the atomic physics picture that directly impact the evolution of the free electrons during recombination. For reference, we present the summary table for state variables affected by these constants \citep{Scoccola2009,Hart2017}:
\begin{align}
\label{eq:vfcScale}
\begin{split}
\sigma_{\rm T} \propto \alpha_{\rm EM}^2 m_{\rm e}^{-2}
\qquad A_{2\gamma} &\propto \alpha_{\rm EM}^8 m_{\rm e}
\qquad P_{\rm S} A_{1\gamma} \propto \alpha_{\rm EM}^{6}m_{\rm e}^{3}
\\
\alpha_{\rm rec} \propto \alpha_{\rm EM}^2 m_{\rm e}^{-2}
\qquad \beta_{\rm phot} &\propto \alpha_{\rm EM}^5 m_{\rm e}
\qquad T_{\rm eff} \propto \alpha_{\rm EM}^{-2}m_{\rm e}^{-1}.
\end{split}
\end{align}
We can extend this to the scalings to effective rate coefficients that are very important for full calculations to recombination as,
\begin{align}
\mathcal{A}_i\left(T_\gamma, T_{\rm e}\right) &\rightarrow \alpha_{\rm EM}^2m_{\rm e}^{-1}\,\mathcal{A}_i\left(\alpha_{\rm EM}^{-2}m_{\rm e}^{-1}\,T_\gamma, \alpha_{\rm EM}^{-2}m_{\rm e}^{-1}\,T_{\rm e}\right), \\
\mathcal{B}_i\left(T_\gamma, T_{\rm e}\right) &\rightarrow\alpha_{\rm EM}^5m_{\rm e}\,\mathcal{B}_i\left(\alpha_{\rm EM}^{-2}m_{\rm e}^{-1}\,T_\gamma, \alpha_{\rm EM}^{-2}m_{\rm e}^{-1}\,T_{\rm e}\right), \\
\mathcal{R}_{ij}\left(T_\gamma\right) &\rightarrow \alpha_{\rm EM}^5m_{\rm e}\,\mathcal{R}_{ij}\left(\alpha_{\rm EM}^{-2}m_{\rm e}^{-1}\,T_\gamma\right).
\end{align}
Here the effective recombination and photoionisation rates are $\mathcal{A}_i$ and $\mathcal{B}_i$ respectively, whereas $\mathcal{R}_{ij}$ represents the transitions between excited states \citep{Yacine2010}. For the spectral conductances \citep{Yacine2013RecSpec, Chluba2016CosmoSpec} that are required to calculate the CRR, the coefficients $\mathcal{G}_{n' n}^{X}$ are scaled by the same factor as the transition rates $\mathcal{R}_{ij}$. It is important to note that for these variations, we have not assumed a particular model and therefore assume no knowledge of an external field that could potentially manipulate the underlying background cosmology. This in turn means that we have not considered any modifications to the Hubble flow $H(z)$ arising from such fields. \footnote{This is a potential direction for more complex fundamental constant variations and may even couple to quintessence-like fields that resemble the EDE discussed in Sec.~\ref{sec:ede} \citep{Calabrese2011}.}

\vspace{-3mm}
\subsection{Modifications to {\tt CosmoSpec}}
\label{sec:mods_CosmoSpec}

\subsubsection{Rescaling effective conductances}
\label{sec:Gnu}
To include the effect of varying fundamental constants (VFCs) on the CRR, we follow the description presented in Appendix~B of \citet{Chluba2016CosmoSpec} and revisited in Sec.~\ref{sec:recap_vfc}. A few important differences with respect to the original version of {\tt CosmoSpec} are:
\begin{itemize}
    \item the \ion{He}{II} spectrum is now computed using rescaled conductances of \ion{H}{I}. This reduces the storage of data and we confirmed the results carefully by direct computation.
    \item the redshift range over which the conductances are tabulated was extended, as VFCs can allow recombination to occur at higher and lower temperatures than in the standard scenario.
    \item the effect of electron scattering is included for modified scattering cross section, with rescaled $y$-parameter obtained as $y'\propto \sigma_{\rm T}/m_{\rm e} \simeq (\alpha_{\rm EM}'/\alpha_{\rm EM})^2\,(m_{\rm e}'/m_{\rm e})^{-3}\,y$ and updated recombination history.
    \item similarly the effect of free-free absorption is modelled with the optical depth scaled as $\tau_{\rm ff}'\simeq (\alpha_{\rm EM}'/\alpha_{\rm EM})^3\,(m_{\rm e}'/m_{\rm e})^{-1.5}\,\tau_{\rm ff}$ and modified recombination history.
\end{itemize}
Aside from these rather straightforward modifications we also improved the analytic treatment of photon escape from the main \ion{He}{I} resonances, as explained next.

\subsubsection{Treatment of \ion{H}{I} absorption during \ion{He}{I} recombination}
\label{sec:DPesc}
One of the important corrections to the helium recombination history at $z\simeq 1700-2000$ is the effect of neutral hydrogen continuum absorption, which leads to a significant acceleration of the recombination process \citep{Kholupenko2007, Switzer2007I, Jose2008}. To approximately model this process, one can compute the correction to the Sobolev escape probabilities of the main singlet and triplet resonances as \citep[see Appendix~B  of][]{Jose2008}:
\begin{align}
  \label{eq:DPesc_chi_1D}
  \Delta P^{1\rm D}_{\rm esc}
  &\approx\!\!\int_0^1 \!\!\!{\,\rm d}\chi\!
  \left\{1-\expf{-{\tau_{\rm  S}}(1-\chi)}-\kappa(\chi)
  \left[1-\expf{-[{\tau_{\rm  S}}+{\tilde{\tau}_{\rm  c}}(\chi)](1-\chi)}\right]\right\},
\end{align}
where $\chi=\int_{-\infty}^{{x_{\rm D}}} \phi(y, a){\,\rm d} y$ is the integral over the Voigt-profile of the resonance, $\phi({{x_{\rm D}}}, a)$, with Voigt-parameter $a$; ${{x_{\rm D}}}$ is the distance from line center at frequency $\nu_0$ expressed in Doppler-widths; ${\tau_{\rm  S}}$ is the Sobolev optical depth of the line. We furthermore defined the \ion{H}{I} continuum opacity variables
\begin{subequations}
\begin{align}
{\tilde{\tau}_{\rm  c}}(\chi)&=\frac{c\, N^{\rm H}_{\rm 1s} \sigma^{\rm H}_{\rm 1s}(\nu)}{H}
\,\frac{\Delta\nu_{\rm D}}{\nu}\,\frac{1}{\phi({{x_{\rm D}}}, a)},
\\[1mm]
\kappa(\chi)&=\frac{\tau_{\rm S}}{{\tau_{\rm  S}}+{\tilde{\tau}_{\rm  c}}(\chi)},
\end{align}
\end{subequations}
where one should think of $\nu$ and ${{x_{\rm D}}}$ as functions of $\chi$. The Doppler width of the line is given by
\begin{align}
\frac{\Delta\nu_{\rm D}}{\nu_0}&=\sqrt{\frac{2 kT_{\rm e}}{{m_{\rm He}} c^2}}
\approx
\pot{1.7}{-5} \left[\frac{(1+z)}{2500}\right]^{1/2},
\end{align}
and directly depends on the electron temperature, with the typical value given for $T_{\rm e}=T_{\rm CMB}$.

It was shown that Eq.~\eqref{eq:DPesc_chi_1D} provides a good first approximation to the escape probability corrections of the \ion{He}{I} singlet and triplet resonances. This can be used to model the main effect of \ion{H}{I} continuum absorption on the recombination process and further corrections can then be added using the radiative transfer module of {\tt CosmoRec/CosmoSpec} \citep{Chluba2012HeRec}.

Instead of computing the integral in Eq.~\eqref{eq:DPesc_chi_1D} repeatedly, in {\tt CosmoRec/CosmoSpec} the escape probability correction is pre-tabulated as a function of ${\tau_{\rm  S}}$, $\eta_{\rm  c}=c\, N^{\rm H}_{\rm 1s} \sigma^{\rm H}_{\rm 1s}(\nu_0)/H$ and $T_{\rm e}$. This eases the computation for the standard recombination history. However, to include the effect of varying fundamental constants, the strategy has to be slightly changed. The parameter ${\tau_{\rm  S}}$ can be used as before, as it is computed internally for the rescaled atomic properties. The other dependencies on atomic properties and the electron temperature enter through the mapping of ${{x_{\rm D}}}={{x_{\rm D}}}(\chi, a)$, which is a function of the Voigt-parameter $a$, and also $\nu/\nu_0=1+{{x_{\rm D}}}(\chi, a)\,\Delta\nu_{\rm D}/\nu_0\equiv f(\chi, a, T_{\rm e})$. Overall, this suggests that $\Delta P^{1\rm D}_{\rm esc}=\Delta P^{1\rm D}_{\rm esc}({\tau_{\rm  S}}, \eta_{\rm  c}, T_{\rm e}, a)$.
By writing
\begin{align}
{\tilde{\tau}_{\rm  c}}
&=
\eta_{\rm  c}\,
\frac{\Delta\nu_{\rm D}}{\nu_0}
\frac{1}{\phi(\chi, a)}\,\frac{\nu_0\, \sigma^{\rm H}_{\rm 1s}(\nu)}{\nu \, \sigma^{\rm H}_{\rm 1s}(\nu_0)}
\end{align}
we can realize that the explicit dependence of ${\tilde{\tau}_{\rm  c}}$ on $\Delta\nu_{\rm D}/\nu_0$ can be captured by using ${\tilde{\eta}_{\rm  c}}=\eta_{\rm  c}\,\Delta\nu_{\rm D}/\nu_0$ as a new parameter for tabulation. Since in the optically-thick regime most of the contributions to the escape probability come from around $\nu\simeq \nu_0$, this in fact means that the explicit dependence on $\nu/\nu_0$ can be omitted. Consequently, a 3D table in ${\tau_{\rm  S}}$, ${\tilde{\eta}_{\rm  c}}$ and $a$ should provide a good representation for the main scaling of $\Delta P^{1\rm D}_{\rm esc}({\tau_{\rm  S}}, \eta_{\rm  c}, T_{\rm e}, a)\approx \Delta P^{1\rm D}_{\rm esc}({\tau_{\rm  S}}, {\tilde{\eta}_{\rm  c}}, a)$.

In practice, we can simply go back to the previous tabulation scheme but keep the electron temperature instead of the Voigt-parameter $a$. Since $a=A_{21}/[4\pi \Delta\nu_{\rm D}]\propto A_{21}/[\nu_0
\sqrt{T_{\rm e}}]$, we have
\begin{align}
a'(T_{\rm e})=a(T_{\rm e}) \, \frac{A'_{21}}{A_{21}}\frac{\nu_0}{
\nu'_0}\equiv a(f_{\rm V} T_{\rm e})
\end{align}
with $f_{\rm V}=\left(A_{21}\,\nu'_0/[A'_{21} \nu_0]\right)^2=(\alpha_{\rm EM}'/\alpha_{\rm EM})^{-6}$, where we used the fundamental constant scaling $A_{21}/\nu_0\propto \alpha_{\rm EM}^3$. This relation allows one to trade $a$ in terms of $T_{\rm e}$.
Preparing the tables over
${\tau_{\rm  S}}$, $\eta_{\rm  c}$ and $T_{\rm e}$, we can then obtain scaled versions as
\begin{align}
P^{1\rm D}_{\rm esc}({\tau_{\rm  S}}, \eta_{\rm  c}, T_{\rm e})
\rightarrow
P^{1\rm D}_{\rm esc}({\tau_{\rm  S}}', f_{\eta_{\rm  c}} \eta_{\rm  c}', f_{\rm V} T_{\rm e})
\end{align}
where $f_\eta=1/\sqrt{f_{\rm V}}$ restores the temperature dependence of ${\tilde{\eta}_{\rm  c}}$ and ${\tau_{\rm  S}}'$ and $\eta_{\rm  c}'$ are evaluated using the scaled atomic variables.
For the tables we use $T_{\rm e}=300-10^5$~K and $\eta_{\rm  c}=10^{-8}-10^8$. For the singlet $\HeIlevel{2}{1}{P}{1}-\HeIlevel{1}{1}{S}{0}$ line, we use ${\tau_{\rm  S}}=10^{-5}-10^{10}$, while for the triplet $\HeIlevel{2}{3}{P}{1}-\HeIlevel{1}{1}{S}{0}$ resonance ${\tau_{\rm  S}}=10^{-5}-10^{3}$ suffices.\footnote{This is motivated by values that are encountered for the standard recombination problem but then scaled for a range of $\simeq 20\%$ $\alpha$-variations.}
For ${\tau_{\rm  S}}\leq 10^{-5}$, linear extrapolation is applied. The grid is setup in log-space with a density of 10 to 15 points per decade. A four-point Lagrange polynomial interpolation is performed in each dimension.
Indeed we find that this procedure works extremely well (usually to better than 0.1\% precision) over a wide range of the parameters.

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/DPesc_S_alp_1.pdf}
    \\
    \includegraphics[width=\linewidth]{figures/DPesc_T_alp_1.pdf}
    \\
    \caption{Dependence of $P^{1\rm D}_{\rm esc}$ for the singlet and triplet lines on $T_{\rm e}$, $\eta_{\rm  c}$ and ${\tau_{\rm  S}}$ for the standard value of $\alpha_{\rm EM}$ (see text for discussion). The results were obtained by direct integration, but are accurately represented by our tabulation scheme.}
    \label{fig:DPesc_LCDM}
    \vspace{-3mm}
\end{figure}
In Fig.~\ref{fig:DPesc_LCDM}, we illustrate the dependence of $P^{1\rm D}_{\rm esc}$ for the singlet and triplet lines on $T_{\rm e}$, $\eta_{\rm  c}$ and ${\tau_{\rm  S}}$ for the standard value of $\alpha_{\rm EM}$. The temperature dependence decreases with $\eta_{\rm  c}$, as expected from the fact for $\eta_{\rm  c}\rightarrow 0$ the dependence on $a=a(T_{\rm e})$ and $T_{\rm e}$ drops out of the expression in Eq.~\eqref{eq:DPesc_chi_1D}. For the triplet line, even when $\eta_{\rm  c}\simeq 1$ a significant dependence on $T_{\rm e}$ can be observed, while for the singlet case the corresponding curves become practically independent of $T_{\rm e}$. As we will see below (Fig.~\ref{fig:DPesc_var_alpha}), the escape probability of the triplet line is practically independent of $a$, such that the main temperature dependence only enters through ${\tilde{\eta}_{\rm  c}}=\eta_{\rm  c}\,\Delta\nu_{\rm D}/\nu_0$, while for the singlet line also $a$ matters.

At ${\tau_{\rm  S}}\ll1$, a quasi-linear scaling with ${\tau_{\rm  S}}$ is found (see Fig.~\ref{fig:DPesc_LCDM}), which directly follows when performing a Taylor-series expansion of  Eq.~\eqref{eq:DPesc_chi_1D} to linear order in ${\tau_{\rm  S}}$:
\begin{align}
  \label{eq:DPesc_chi_1D_corr}
  \Delta P^{1\rm D}_{\rm esc}
  &\approx {\tau_{\rm  S}} \int_0^1 \!\!\!{\,\rm d}\chi\!
  \left\{1-\chi-
  \frac{1-\expf{-{\tilde{\tau}_{\rm  c}}(\chi)[1-\chi]}}{{\tilde{\tau}_{\rm  c}}(\chi)}\right\}
\end{align}
and justifies the adopted extrapolation procedure for ${\tau_{\rm  S}}\rightarrow 0$.

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/DPesc_S_var_alp.pdf}
    \\
    \includegraphics[width=\linewidth]{figures/DPesc_T_var_alp.pdf}
    \\
    \caption{Dependence of $P^{1\rm D}_{\rm esc}$ for the singlet and triplet lines on $\alpha_{\rm EM}'/\alpha_{\rm EM}$, $\eta_{\rm  c}$ and ${\tau_{\rm  S}}$ for the standard value of $\alpha_{\rm EM}$. The values of ${\tau_{\rm  S}}$ and $\eta_{\rm  c}$ are directly evaluated for the scaled atomic species and therefore do not directly modify the escape probability. In this sense, all the visible changes are due to changes of the Voigt-parameter $a'/a\propto (\alpha_{\rm EM}'/\alpha_{\rm EM})^{-3}$. The results were obtained by direct integration, but are accurately represented by our tabulation scheme.}
    \label{fig:DPesc_var_alpha}
\end{figure}
To illustrate the dependence of the escape probability on the fundamental constants, we note that the escape integral has no explicit dependence on $m_{\rm e}$, but only indirectly from the changes of ${\tau_{\rm  S}}$ and $\eta_{\rm  c}$ at fixed physical parameters like temperature and density. The only explicit dependence enters due to modifications of $a$. This is shown in Fig.~\eqref{fig:DPesc_var_alpha}. The singlet line exhibits noticeable changes for all shown values of $\eta_{\rm  c}$, while the triplet line is practically independent of $a$.
This arises from the much smaller typical value of $a$ for the triplet line and associated dominance of the Doppler core in the escape problem \citep[see Fig.~B.1 of][for some illustrations of the integrand]{Jose2008}, which makes the escape problem for the triplet line less dependent on $a$.
For the singlet line we find $\Delta P^{1\rm D}_{\rm esc}\propto (\alpha_{\rm EM}'/\alpha_{\rm EM})^{1.5}$ for the shown cases.
Our tabulation scheme described above captures all dependencies accurately.

\subsection{Effects on the recombination lines}\label{sec:linesVFC}
The impact of varying the fine structure constant on the cosmological recombination lines is shown in Fig.~\ref{fig:vfcCrr}. As in Sec.~\ref{sec:DPesc}, the fiducial value is compared to the two cases with $\Delta\alpha_{\rm EM}/\alpha_{\rm EM} = \pm \,0.1$. In Fig.~\ref{fig:vfcCrr}, the modified recombination lines are shown for HI \emph{(top)}, HeI \emph{(middle)} and HeII \emph{(bottom)}.
For hydrogen, the net impact of a larger $\alpha_{\rm EM}$ is an amplification on the spectral features. We find that the high frequency peaks roughly change as $\propto\alpha_{\rm EM}^2$. Given the primary recombination mechanism described in Sec.~\ref{sec:recap_vfc}, with an acceleration of recombination expected for increasing $\alpha_{\rm EM}$, this makes sense; however, the detailed effects are quite subtle, where the impact on recombination is a complicated combination of line enhancements versus broadening. The latter leads to a redistribution of photons across the CRR and hence reduction of emission in some bands, with the width being directly linked to the relative duration of the recombination process. We also note that at low frequencies (not shown here), we do not see any significant change in the amplitude of the distortion, aside from some modifications from the free-free absorption process. Given that at low frequencies no spectral features are visible, this shows that a net redistribution of the emission in different bands occurs while leaving the total number of photons added nearly constant. Thus visible changes to the CRR only appear where individual lines can be identified.

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/alpha_HI.pdf}
    \\[-0.5mm]
    \includegraphics[width=\linewidth]{figures/alpha_HeI.pdf}
    \\[-0.5mm]
    \includegraphics[width=\linewidth]{figures/alpha_HeII.pdf}
    \\[-0.5mm]
    \caption{Changes in the cosmological recombination lines due to variations in $\alpha_{\rm EM}$ for the values $\alpha_{\rm EM}/\alpha_{\rm EM, 0} = \{0.9,1.0,1.1\}$. The dashed curve is shown for the $\Lambda$CDM\xspace case ($\alpha_{\rm EM}/\alpha_{\rm EM, 0} = 1.0$). This has been shown for the hydrogen lines \emph{(top)}, singly ionised helium \emph{(middle)} and doubly ionised helium \emph{(bottom)}.}
    \label{fig:vfcCrr}
\end{figure}
For the helium lines, additional complications arise from the changes to the importance of electrons scattering. With increasing $\alpha_{\rm EM}$ the effective scattering $y$-parameter increases and so does the smearing of the lines. This is indeed visible in the HeII spectrum. We can notice that the HeII lines remain almost constant in amplitude but mainly change their width. This highlights that for the total spectrum it is hard to understand how the effects propagate into the final CRR just from the simple scaling of variables given above.

An increase in $\alpha_{\rm EM}$ also causes a small net drift of the main lines for all species of distortion to lower frequencies. However, the effect is much smaller than what would be guessed from the energy scaling of the transition frequencies, $\nu_{ij}\propto \alpha_{\rm EM}^2$. In reality, the recombination process also occurs at roughly $\alpha_{\rm EM}^2$ times higher redshifts, leaving the ratio $\nu_{ij}/(1+z_{\rm em})$ roughly constant for individual transitions. Nevertheless, the positions of the lines are tracers of when the recombination process happened. Interestingly, the shifts are more noticeable for the HI and HeI contributions than for HeII, for which line broadening effects seem to be more pronounced.

For larger $\alpha_{\rm EM}$ we can also see added spectral peaks emerging. Examples are at $\nu\simeq80{\rm GHz}$ and $\nu\simeq160{\rm GHz}$ in the HI contribution. This has to do with the way the emission from various transitions overlaps and compensate each other, partially canceling or interfering constructively.
For HeII (lower part of Fig.~\ref{fig:vfcCrr}), larger values of $\alpha_{\rm EM}$ more strongly smear the overall shape, with peaks being compressed and troughs being raised. In this case, we also notice a deterioration of the peak feature at $\nu\simeq450{\rm GHz}$, indicating modifications in the relative importance of various transitions, here related to the HeII Balmer lines.

When all three atomic species are combined into a single distortion, the result is as shown in Fig.~\ref{fig:vfcCrrTotal}. Here the more intricate variations due to helium are smeared out by the larger magnitude signal from hydrogen, as was the case for early dark energy. However, this still provides richer structure on the inter-level emission bumps throughout the profile. This should generally allow one to distinguish $\alpha_{\rm EM}$ variations in the CRR from other parameters using high resolution spectrometers.

For comparison, in Fig.~\ref{fig:vfcCrrTotal} we also show the same spectra but for changes to $m_{\rm e}$. Although we have doubled the variation to $\Deltam_{\rm e}/m_{\rm e, 0}=\pm0.2$, one can see that the changes of the total CRR are a lot smaller. By decreasing the electron mass one can observe a small overall enhancements of the spectral features.
Therefore, $m_{\rm e}$ will be harder to detect with a CMB spectrometer, with the largest response being visible at about $\simeq 1$~THz.
We confirmed numerically that this has to do with the precise interplay of the emissivities and energies scaling, leaving the main effect to be a small change in the width of the spectral features.
This directly contrasts with the influence that $m_{\rm e}$ has on the CMB anisotropies as previously discussed in the literature \citep{Hart2017,Hart2020a}. However, here the observable is related to the Thomson visibility function, which has a different non-canceling dependence on $m_{\rm e}$.

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/alpha_Total.pdf}
    \includegraphics[width=\linewidth]{figures/me_Total.pdf}
    \caption{The total impact on the recombination lines from the variations shown in Fig.~\ref{fig:vfcCrr} for $\Delta\alpha_{\rm EM}/\alpha_{\rm EM, 0}\pm0.1$ \emph{(top)}. The total variations due to changes in $\Deltam_{\rm e}/m_{\rm e, 0}\pm0.2$ are shown for comparison \emph{(bottom)}.}
    \label{fig:vfcCrrTotal}
\end{figure}

\subsubsection{Relationship between $\alpha_{\rm EM}$ and $T_{\rm CMB}$}\label{sec:alphaT}
The fine structure constant has very clear degeneracies with the monopole temperature of the CMB, $T_{\rm CMB}$ due to the similar effects they have on the last scattering surface position. In Fig.~\ref{fig:dlnI_alpha_T}, we show the relative derivatives for $\alpha_{\rm EM}$ \emph{(purple)} and $T_{\rm CMB}$ \emph{(orange)}. The reference CRR spectrum for hydrogen and helium combined is shown as a dotted line. The large-scale features of the derivatives spectrum align in most cases: for example the positive-negative envelope in $\alpha_{\rm EM}$ at $\nu\simeq1800$ GHz, close to the edge of the Lyman-continuum, is mirrored at a higher amplitude in $T_{\rm CMB}$. Similarly the mirror in the derivatives is also present between $\alpha_{\rm EM}$ and $T_{\rm CMB}$ in the predominantly bound-bound transitions\footnote{Though the free-bound radiation is present, it is void of the features that we typically associate with bound-bound transitions at $\nu\lsim100$ GHz.} at $\nu\lesssim 100$ GHz.
The remainder of the spectra, in the mid-range of the proposed bands by \emph{Voyage 2050}\xspace, show distinct differences between the responses from the two parameters. It was known from previous studies that $\alpha_{\rm EM}$ and $T_{\rm CMB}$ had some degeneracy breaking effects due to $\alpha_{\rm EM}$ affecting the recombination physics while not affecting the global expansion rate (under the current assumptions). In comparison, $T_{\rm CMB}$ affects the positioning of the \emph{matter-radiation equality} epoch as well. This leads to a wider impact at higher redshifts \cite[see][for more details]{Hart2020a}. With the full deviations implicitly found in the interaction rates and the conductance variations (see Sec.~\ref{sec:recap_vfc}), the changes that will primarily affect the Lyman, Balmer and higher order series will be present in the $\alpha_{\rm EM}$ derivatives and less obvious in the $T_{\rm CMB}$ responses. As originally shown in \citet{Chluba2008T0}, the generic effect from changing $T_{\rm CMB}$ is to shift the peaks of the CRR back and forth. Therefore, the interaction rates have a more complex effect from variations in $\alpha_{\rm EM}$.

\subsubsection{Relationship between $m_{\rm e}$ and $N_{\rm eff}$}\label{sec:meNeff}
The covariance that can be highlighted otherwise is the relationship between $m_{\rm e}$ and $N_{\rm eff}$, the relativistic degrees of freedom. The weighted derivatives are shown in Fig.~\ref{fig:dlnI_me_Neff}. Specifically the negative and positive spectra for $m_{\rm e}$ and $N_{\rm eff}$ almost perfectly emulate each other; however, there are small structural changes in the inter-spectral gaps. This is due to the more complicated atomic variations from $m_{\rm e}$ compared with the broader effects on the expansion rate associated with $N_{\rm eff}$. When looking at the finer structure of the spectral responses, differences can be ascertained as the resolution and sensitivity increase. For example, at $\nu\simeq 1500$ GHz in Fig.~\ref{fig:dlnI_me_Neff}, the broader changes in the helium feedback lines for $m_{\rm e}$ are different to $N_{\rm eff}$. Therefore, we would expect that a sufficiently high sensitivity would distinguish electron mass from changes to the relativistic degrees of freedom.

\changeJ{We also note that the changes from $N_{\rm eff}$ are numerically harder to compute, given that very large changes are required to see any responses \citep{Hart2020c}. This is the origin of some numerical imperfections, visible as kinks in Fig.~\ref{fig:dlnI_me_Neff}, however, they do not affect the main conclusions, in particular once external priors are added.}

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/dlnDIv_alphaT.pdf}
    \\
    \caption{Weighted derivatives for changes in the recombination lines for variations of $\alpha_{\rm EM}$ vs. variations in $T_{\rm CMB}$. The negative derivatives are indicated by dashed lines. Both are compared against the $\Lambda$CDM reference spectra  to see how the derivatives compare against the full spectral signal. Similarly, we have used the $\ln p$ weighting for the denominator as this will more appropriately compare to the fiducial CRR.}
    \label{fig:dlnI_alpha_T}
\end{figure}

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/dlnDIv_meNeff.pdf}
    \\
    \caption{Recombination radiation derivatives for the electron mass rescaling against the $N_{\rm eff}$ derivatives originally shown in \citet{Hart2020c} (multiplied by a factor of 10 for readability). These have been weighted by the variations in parameter $\delta p/p$. }
    \label{fig:dlnI_me_Neff}
\end{figure}

\subsection{Possible detections for $\alpha_{\rm EM}$ and $m_{\rm e}$ with the CRR}\label{sec:vfcFisher}
In this section, we use the spectral variations in the CRR caused by a change in $\alpha_{\rm EM}$ and $m_{\rm e}$ to test the constraining strength for future spectrometer missions. This involves a Fisher forecast where one additional parameter ($\alpha_{\rm EM}$, $m_{\rm e}$) modifies the standard $\Lambda$CDM\xspace matrix. This will be shown for the following variety of assumptions for the priors on the standard parameters: no priors and priors from \emph{Planck}\xspace 2018.
The contours in this section were generated using {\tt GetDist} as in Sec.~\ref{sec:forecastEDE} using Gaussian random samples. However in this section, we use the previously generated \emph{Planck}\xspace chains to add inverse covariances to the Fisher matrix \citep{Hart2020a}.

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/alphaVFCTest.pdf}
    \\
    \caption{Rudimentary contour plot for the explored parameters from \citet{Hart2020c} with $\alpha_{\rm EM}$ added into the analysis. The contours are compared for the futuristic configuration \emph{Voyage 2050++}\xspace with and without the added fundamental constants.}
    \label{fig:alphaFisherCRR}
\end{figure}
\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/meVFCTest.pdf}
    \\
    \caption{The same forecasted posterior contours as in Fig.~\ref{fig:alphaFisherCRR} however for variations in electron rest mass $m_{\rm e}$. Note that the configurations and varying parameter sets are the same as well.}
    \label{fig:meFisherCRR}
\end{figure}

\subsubsection{Case 1: Spectrometer only}\label{sec:specOnly}
The most simple case for this is where we forecast the detectability with a futuristic spectrometer such as \emph{Voyage 2050++}\xspace. For this we use the same noise profiles as described in Sec.~\ref{sec:forecastEDE} for \emph{Voyage 2050}\xspace. Specifically, we will use the standard Voyage setup (\emph{Voyage 2050}\xspace) and $50\times$ higher \changeJ{sensitivity} (\emph{Voyage 2050++}\xspace) to highlight the potential of future CMB spectrometers as a cosmological probe.

The comparison for \emph{Voyage 2050++}\xspace for $\Lambda$CDM\xspace and adding $\alpha_{\rm EM}$ is shown in the contours of Fig.~\ref{fig:alphaFisherCRR}.
When varying $\alpha_{\rm EM}$ the contours show a general broadening once the fine structure constant is added to the analysis. This softens the direct degeneracy between $\omega_{\rm c}$ and $N_{\rm eff}$ that we have previously discussed in this paper and previous works \citep{Hart2020c}, transferring it to correlations with $\alpha_{\rm EM}$. However, the contours are not altered to a great amount in general; this is owed to the largely distinctive and unique spectral patterns associated with $\alpha_{\rm EM}$ that we have discussed in Sec.~\ref{sec:linesVFC}. The standard deviations for this configuration with $\alpha_{\rm EM}$ are shown in  Table~\ref{tab:alphaCRR}. From the data in this table, and the previous findings from \citet{Hart2020c}, we can see that $\omega_{\rm c}$ and $N_{\rm eff}$ are far too weakly constrained by \emph{Voyage 2050}\xspace alone, requiring higher sensitivity and complimentary probes to reach competitive results.\footnote{We ask the reader to bear in mind that the errors for $N_{\rm eff}$ in Table~\ref{tab:alphaCRR} and Table~\ref{tab:meCRR} are indeed large compared to modern probes; however, they are quoted here to underpin the potential of CRR oriented spectrometers when the sensitivity can be reached.} As $\alpha_{\rm EM}$ is included in the analysis, we can see from Table~\ref{tab:alphaCRR} that the errors on the parameters do not expand, save for a $\simeq 19\%$ increase for $\omega_{\rm c}$.

\begin{table}
    \centering
    \begin{tabular}{c|c|c|c}
    \hline\hline
        Parameter & \emph{Voyage 2050}\xspace & \emph{Voyage 2050}\xspace & \emph{Voyage 2050++}\xspace \\
         & $\Lambda$CDM\xspace & + $\alpha_{\rm EM}$ & + $\alpha_{\rm EM}$ \\
        \hline
        $\omega_{\rm b}$ & $0.00057$ & $0.00058$ & $0.00002$ \\
        $\omega_{\rm c}$ & $0.186$ & $0.222$ & $0.004$ \\
        $Y_{\rm p}$ & $0.0478$ & $0.0481$ & $0.0010$ \\
        $N_{\rm eff}$ & $20.40$ & $21.30$ & $0.43$ \\
        $\alpha_{\rm EM}$ & \textemdash & $0.059$ & $0.001$ \\
        \hline\hline
    \end{tabular}
    \caption{Initial covariance forecasts from \emph{Voyage 2050}\xspace and \emph{Voyage 2050++}\xspace spectrometers when the parameters listed here are modified, for varying $\alpha_{\rm EM}$. In this case, the \emph{Voyage 2050++}\xspace configuration is $50\times$ more sensitive than the standard \emph{Voyage 2050}\xspace.}
    \label{tab:alphaCRR}
\end{table}

\begin{table}
    \centering
    \begin{tabular}{c|c|c|c}
    \hline\hline
        Parameter & \emph{Voyage 2050}\xspace & \emph{Voyage 2050}\xspace & \emph{Voyage 2050++}\xspace \\
         & $\Lambda$CDM\xspace & + $m_{\rm e}$ & + $m_{\rm e}$ \\
        \hline
        $\omega_{\rm b}$ & $0.00057$ & $0.00059$ & $0.00001$ \\
        $\omega_{\rm c}$ & $0.186$ & $0.309$ & $0.006$ \\
        $Y_{\rm p}$ & $0.0478$ & $0.0481$ & $0.0010$ \\
        $N_{\rm eff}$ & $20.40$ & $27.90$ & $0.56$ \\
        $m_{\rm e}$ & \textemdash & $1.53$ & $0.030$ \\
        \hline\hline
    \end{tabular}
    \caption{Electron mass $m_{\rm e}$ forecasted covariances, $\sigma_{ij}$, from \emph{Voyage 2050}\xspace and \emph{Voyage 2050++}\xspace spectrometers when these parameters are modified. The configurations are the same as the ones in Table~\ref{tab:alphaCRR}.}
    \label{tab:meCRR}
\end{table}

For $m_{\rm e}$, the situation is slightly more complicated as the contours in Fig.~\ref{fig:meFisherCRR} show. Both $\omega_{\rm c}$ and $N_{\rm eff}$ distributions are broadened and the effect from $m_{\rm e}$ all but removes the $\omega_{\rm c}-N_{\rm eff}$ degeneracy. In this case, the co-varying power is passed to the $m_{\rm e}$-$N_{\rm eff}$ relation, as the $m_{\rm e}$ variations partially emulate an acceleration in the radiation era. This is corroborated by the results in Table~\ref{tab:meCRR} where the error on $\omega_{\rm c}$ increased by $66\%$ and $N_{\rm eff}$ by $36\%$ for \emph{Voyage 2050}\xspace.

We note that none of the other $\Lambda$CDM\xspace parameters matter for the CRR and that the value of $T_{\rm CMB}$ will be measured to very high precision using the blackbody part. Of course the precludes cosmologies with varying temperature-redshift relation, which would run into several other issues of course \citep{Chluba2014TRR}.

\subsubsection{Case 2: Spectrometer with Planck 2018}\label{sec:specPlanck}

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/alphaVFC_updated.pdf}
    \\
    \caption{Forecasted contours for variations in $\alpha_{\rm EM}$ that correspond to the results in Table~\ref{tab:alphaForecast}. The emphasised parameter set $\{\omega_{\rm b},\omega_{\rm c},\alpha_{\rm EM},100\,\theta_{\rm MC},n_{\rm s}\}$ show the most important variations when the fine structure constant is considered with spectral distortions. \changeL{Here the parameter drift in $\alpha_{\rm EM}$ arises from assuming a fiducial picture of $\alpha_{\rm EM}/\alpha_{\rm EM, 0}=1.0$ for the CRR case.}}
    \label{fig:alphaFullContour}
\end{figure}
From the spectrometer only results, the next stage is adding the \emph{Planck}\xspace covariance into the analysis. Specifically, we can do this by adding the Fisher matrix $F_{ij}$ to the inverse covariance matrix from \emph{Planck}\xspace ${\bf: \Sigma}_{ij}$. Note that in this case the parameters coming from \emph{Planck}\xspace ($\theta_{\rm MC},\tau,A_{\rm s},n_{\rm s}$) have a zero contribution in the CRR Fisher matrix; we are simply investigating how the parameter errors change as the covariances are influenced.

Adding different mission configurations with \emph{Planck}\xspace data for variations of $\alpha_{\rm EM}$ are shown in Fig.~\ref{fig:alphaFullContour}. The standard $\Lambda$CDM\xspace parameters associated with the the power spectrum amplitude and reionisation era, $A_s$ and $\tau$ have been omitted from the plots as the addition of CRR does not alter their values significantly. The configuration for the highest sensitivity is \emph{Voyage 2050++}\xspace since this is the level of sensitivity in our setup where the errors start to markedly improve compared to \emph{Planck}\xspace. In our forecasts with $\alpha_{\rm EM}$ and $m_{\rm e}$, we have only included the two energy densities: $\omega_{\rm b}$ and $\omega_{\rm c}$ as added free parameters since we have not included chains where $Y_{\rm p}$ and $N_{\rm eff}$ co-vary with the fundamental constants. Responses in $H_0$ are so small in the CRR they have not been included here \cite[see][for more details]{Hart2020c}.

Firstly, the \emph{Planck}\xspace contours are shown as yellow dashed lines, with the addition of \emph{Voyage 2050}\xspace data in dark purple. The contours for \emph{Planck}\xspace are almost identical whether one includes \emph{Voyage 2050}\xspace or omits it in the case of $\alpha_{\rm EM}$. However, there are some small shifts in the maximum likelihood positions. This is reflected in Table~\ref{tab:alphaForecast} where the values of the errors do not change for this configuration either. When we look at a higher precision \emph{Voyage 2050++}\xspace configuration, the errors start to diminish as shown in Fig.~\ref{fig:alphaFullContour}. Including $\alpha_{\rm EM}$ with \emph{Planck}\xspace+\emph{Voyage 2050++}\xspace forecasts, the finer spectral shape of $\alpha_{\rm EM}$ is complimented by the much higher sensitivity and the error on $\alpha_{\rm EM}$ is significantly reduced. \changeJ{Some gains can even be seen for other parameters like $n_{\rm s}$ due to a reduction of the $\alpha_{\rm EM}$-$100\,\theta_{\rm MC}$ contours.}

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/alphaVFC_2d.pdf}
    \\
    \caption{Single contours from Fig.~\ref{fig:alphaFullContour} specifically focused on the $\alpha_{\rm EM}$-$100\,\theta_{\rm MC}$ degeneracy to showcase the much better forecasting that \emph{Voyage 2050++}\xspace gives as well as the similar degeneracy line; note here that both \emph{Voyage 2050}\xspace and \emph{Voyage 2050++}\xspace have been folded in with \emph{Planck}\xspace errors.}
    \label{fig:alpha2d}
\end{figure}

The degeneracy of $\alpha_{\rm EM}$ with $100\,\theta_{\rm MC}$ is not removed; however, since the parameter degeneracies with $\theta$ are not altered by other parameters, the correlations remain similar (correlation between the parameters $\xi_{\theta,\alpha} = 0.99 \rightarrow0.91$). This has been highlighted in Fig.~\ref{fig:alpha2d}. When the high precision \emph{Voyage 2050++}\xspace configuration is considered, even with latent degeneracies, the error on $\alpha_{\rm EM}$ is diminished to $\sigma_{\alpha_{\rm EM}} \simeq 0.0005$, which is $\simeq 5$ times smaller than \emph{Planck}\xspace alone and also another factor of $\simeq 2$ better than the \emph{Voyage 2050++}\xspace spectrometer alone. Note that for CMB anisotropy measurements, the addition of BAO data did not change the error on $\alpha_{\rm EM}$ in previous studies \citep{Hart2020a}.

\begin{figure}
    \centering
    \includegraphics[width=\linewidth]{figures/meVFC_updated.pdf}
    \caption{Forecasted contours for variations in $m_{\rm e}$ that correspond to the results in Table~\ref{tab:alphaForecast}. The emphasised parameter set $\{\omega_{\rm b},\omega_{\rm c},m_{\rm e},100\,\theta_{\rm MC},n_{\rm s}\}$ show the most important variations when the fine structure constant is considered with spectral distortions.\changeL{As in Fig.~\ref{fig:alphaFullContour}, the shifts in $m_{\rm e}$ emulate the drift expected between \emph{Planck}\xspace and the CRR configurations associated with the change in errors and the CRR fiducial case being $m_{\rm e}/m_{\rm e, 0}=1.0$.}}
    \label{fig:meFullContour}
\end{figure}
Conversely, the influence of $m_{\rm e}$ has a more substantial effect. In \emph{Planck}\xspace 2018, the errors for $m_{\rm e}$ are much larger since the marginalised result is $m_{\rm e}/m_{\rm e, 0} = 0.888\pm0.059$ as shown in Table~\ref{tab:meForecast}.
In particular, the non-zero off-set with respect to the standard electron mass helped alleviate the Hubble tension \citep{Hart2020a}.
However, already when we combine $\emph{Voyage 2050}\xspace$ with \emph{Planck}\xspace, we see the error in $m_{\rm e}$ to cascade down and move to the fiducial value: the density covariances with $m_{\rm e}$ are limited by the spectrometer. The error in $m_{\rm e}$ is already $\simeq 3$ times smaller in this case, allowing one to contest the non-zero shift in the value for $m_{\rm e}$ seen from \emph{Planck}\xspace alone.
\changeL{Primarily, the error on $m_{\rm e}$ in \emph{Planck}\xspace is caused by the geometric degeneracy between $m_{\rm e}$ and $100\,\theta_{\rm MC}/H_0$. Since the changes in $\Delta I_\nu$ caused by variations of $m_{\rm e}$ are distinct from any associated with the Hubble constant, this separates the variations of the two parameters and begins to remove the correlation.}
When we further look at \emph{Planck}\xspace with \emph{Voyage 2050++}\xspace, the error drops by a another factor of $\simeq7$. Evidently, this very futuristic version of the \emph{Voyage 2050}\xspace spectrometer will give unprecedented precision on $m_{\rm e}$ which may help rule out a possible VFC origin of the Hubble tension.
\changeJ{We note that adding \emph{Voyage 2050+}\xspace to \emph{Planck}\xspace (not shown here) did not lead to significant additional improvements on $m_{\rm e}$ since the gains on other CMB parameters did eat up some of the gains from the increased spectrometer sensitivity.}

\begin{table}
    \centering
    \begin{tabular}{c|c|c|c}
     \hline
     Parameters & \emph{Planck}\xspace only  & \emph{Voyage 2050}\xspace & \emph{Voyage 2050++}\xspace \\
      & ($\mu_i \pm \sigma_i$) & \& \emph{Planck}\xspace & \& \emph{Planck}\xspace\\
     \hline\hline
       $\omega_b$ & $0.02236\pm0.00015$ & $0.00014$ & $0.00001$ \\
       $\omega_c$ & $0.1201\pm0.0014$ & $0.0014$ & $0.0010$ \\
       $\alpha_{\rm EM}$ & $1.0005\pm0.0024$ & $0.0024$ & $0.0005$ \\
       $100\, \theta_{\rm MC}$ & $1.0416\pm0.0034$ & $0.0034$ & $0.0007$ \\
       $\tau$ & $0.0540\pm0.0075$ & $0.0075$ & $0.0070$ \\
       $\ln(10^{10} A_s)$ & $3.043\pm0.015$ & $0.015$ & $0.014$ \\
       $n_s$ & $0.9637\pm0.0070$ & $0.0069$ & $0.0041$ \\
       \hline\hline
    \end{tabular}
    \caption{Forecasted standard deviations found for different mission configurations combining CMB anisotropies (\emph{Planck}\xspace) and spectral distortions (\emph{Voyage 2050}\xspace, \emph{Voyage 2050++}\xspace). Note that here the isolated \emph{Planck}\xspace parameters $\left\{100\,\theta_{\rm MC},\tau,\ln\left(10^{10} A_{\rm s}\right),n_{\rm s}\right\}$ are added into the Fisher and then recalculated with the CRR influence. The '\emph{Planck}\xspace only' parameters have been quoted by their marginalised values as well as their standard deviations.}
    \label{tab:alphaForecast}
\end{table}

\begin{table}
    \centering
    \begin{tabular}{c|c|c|c}
     \hline
     Parameters & \emph{Planck}\xspace only & Voyage 2050  & \emph{Voyage 2050++}\xspace  \\
      & ($\mu_i \pm \sigma_i$) & \& \emph{Planck}\xspace & \& \emph{Planck}\xspace \\
     \hline\hline
       $\omega_b$ & $0.0199^{+0.0012}_{-0.0014}$ & $0.00044$ & $0.00001$ \\
       $\omega_c$ & $0.1058\pm0.0076$ & $0.0029$ & $0.0007$ \\
       $m_{\rm e}$ & $0.888\pm0.059$ & $0.020$ & $0.0026$ \\
       $100\, \theta_{\rm MC}$ & $0.958\pm0.045$ & $0.0155$ & $0.0022$ \\
       $\tau$ & $0.0512\pm0.0077$ & $0.0074$ & $0.0069$ \\
       $\ln(10^{10} A_s)$ & $3.029\pm0.017$ & $0.015$ & $0.014$ \\
       $n_s$ & $0.9640\pm0.0040$ & $0.0040$ & $0.0034$ \\
       \hline\hline
    \end{tabular}
    \caption{Similar error forecasts for $m_{\rm e}$ following the format of Table~\ref{tab:alphaForecast}. As before, the configurations are the same and the 'derived' parameters from \emph{Planck}\xspace are not directly probed, rather they are modified by the CRR impact.}
    \label{tab:meForecast}
\end{table}

\changeJ{We also mention that in our forecasts we have not allowed for time-dependent variations of the fundamental constants. This had interesting effects on the CMB anisotropies and could be independently constrained \citep{Hart2017, Hart2021}. Here, we have three emission eras from each of the individual atomic species \citep{Sunyaev2009}. This means one could expect the interesting inter-atomic interplay to manifest in the CRR and even be constrainable with future CMB spectrometers. However, a more detailed analysis of this problem is beyond the scope of this paper.}

\section{Conclusions}\label{sec:conc}
In this work, we studied the effects of EDE and VFC on the CRR. We illustrated how the various model parameters affect the CRR and provided simple forecasts for the expected sensitivities of various CMB spectrometer concepts.

The effect of EDE solely enters through changes of the expansion rate in the pre-recombination era. The associated effects on the three recombination phases depend on the details of the underlying model-parameters. For radiation-like dilution, all three recombination contributions can be modified, while for significantly faster dilution, the responses remain more localized in redshift. \changeJ{These effects in principle allow probing EDE models with spectrometers comparable to \emph{Voyage 2050++}\xspace. Since we do not have access to chains from {\sc Planck} for EDE model, we have not explored how a combination with a spectrometer could improve the constraints. However, we do anticipate significant gains but leave a more detailed analysis to future work.}

For VFC models, the recombination physics is directly affected, leaving distinct responses in the CRR that in principle again allow testing various phases in the pre-recombination era. The leading order effects for variations of $\alpha_{\rm EM}$ are a change in the amplitude of the CRR and smaller shift in the position of the recombination lines. The latter effect appears degenerate to changes in the value of the CMB monopole temperature; however, as we show here, the responses are distinct, in principle allowing to distinguish the two (see Fig.~\ref{fig:vfcCrrTotal}) directly with the CRR.
For variations of $m_{\rm e}$, we find the responses to be much smaller with significant cancellations between various effects (see lower panel of Fig.~\ref{fig:vfcCrrTotal}).

\changeJ{Our simple forecasts show that the CRR provides the principle possibility to test EDE and VFC models. However, futuristic spectrometer sensitivities are required to derive independent but competitive constraints. By combining with CMB anisotropy measurements, significant improvements can be found. For example, \emph{Voyage 2050}\xspace with {\sc Planck} could improve the allowed error on variations of $m_{\rm e}$ by a factor of more than $\simeq 3$ and also remove a large part of the geometric degeneracy allowed by {\sc Planck} alone (see Fig.~\ref{fig:meFullContour}). This might shed new light on the origin of the Hubble tension, allowing us to rule out a VFC cause.}

\changeJ{A more comprehensive forecast, that combines the CRR responses with other cosmological probes could improve the forecasts, potentially allowing to distinguish between various scenarios.
The addition of foregrounds will also be important. However, in contrast to $\mu$ and $y$ type distortions \citep{Abitbol2017}, the CRR does not suffer as strongly from foregrounds \citep{Hart2020c}, such that the main conclusions should not change as much.
In addition, EDE models can simultaneously create VFC effects, which then potentially enhance the sensitivity of the CRR to the underlying physics model.
An exploration of these possibilities is left to future work.}

\section{Data Availability}\label{sec:data}
The recombination lines were simulated and modelled using {\tt CosmoSpec} along with modifications for non-standard physics. The modelling of the Fisher matrix was created with the repository {\tt vfcFisher}, which will be made public upon the publication of this article\footnote{\url{http://www.github.com/cosmologyluke/vfcFisher}}.

\section*{Acknowledgements}
The authors would like to thank Alan Heavens for useful discussions surrounding the philosophy of Fisher matrices and likelihoods. All contour plots in this paper were made using the computing package {\tt GetDist}. This work was supported by the ERC Consolidator Grant CMBSPEC (No. 725456) as part of the European Union’s Horizon 2020 research and innovation program. LH would like to acknowledge the support of TNEI Services Ltd.

\bibliographystyle{mn2e}