Classes in C++

The Cosmo class

The mimes::Cosmo< LD > class is used to interpolate various quantities of the plasma. It is defined inside the header file MiMeS/src/Cosmo/Cosmo.hpp. The constructor is

    template< class LD >
    mimes::Cosmo< LD >(std::string cosmo_PATH, LD minT=0, LD maxT=mimes::Cosmo< LD >::mP)

The argument cosmo_PATH is the path of the data file that contains $T$ (in GeV), $h_{eff}$, $g_{eff}$, with increasing $T$. The parameters minT and maxT are minimum and maximum interpolation temperatures. These temperatures are just limits, and the action interpolation is done between the closest temperatures in the data file. Moreover, beyond the interpolation temperatures, both $h_{eff}$ and $g_{eff}$ are assumed to be constants.

Interpolation of the RDOF, allows us to define various quantities related to the plasma. These quantities are given as the member functions:

1. template< class LD > LD mimes::Cosmo< LD >::heff(LD T): $h_{eff}$ as a function of $T$.
2. template< class LD > LD mimes::Cosmo< LD >::geff(LD T): $g_{eff}$ as a function of $T$.
3. template< class LD > LD mimes::Cosmo< LD >::dheffdT(LD T): $dh_{eff}/dT$ as a function of $T$.
4. template< class LD > LD mimes::Cosmo< LD >::dgeffdT(LD T): $dg_{eff}/dT$ as a function of $T$.
5. template< class LD > LD mimes::Cosmo< LD >::dh(LD T): $\delta_h = 1 + \frac{1}{3} \frac{d\log h_{eff}}{d\log T}$ as a function of $T$.
6. template< class LD > LD mimes::Cosmo< LD >::s(LD T): The entropy density of the plasma as a function of $T$.
7. template< class LD > LD mimes::Cosmo< LD >::rhoR(LD T): The energy density of the plasma as a function of $T$.
8. template< class LD > LD mimes::Cosmo< LD >::Hubble(LD T): The Hubble parameter assuming radiation dominated expansion as a function of $T$.

Moreover, there are several cosmological quantities are given as members variables:¹

1. template< class LD > constexpr static LD mimes::Cosmo< LD >::T0: CMB temperature today in GeV.
2. template< class LD > constexpr static LD mimes::Cosmo< LD >::h_hub: Dimensionless Hubble constant.
3. template< class LD > constexpr static LD mimes::Cosmo< LD >::rho_crit: Critical density in GeV$^3$.
4. template< class LD > constexpr static LD mimes::Cosmo< LD >::relicDM_obs: Central value of the measured DM relic abundance.
5. template< class LD > constexpr static LD mimes::Cosmo< LD >::mP: Planck mass in GeV.

The AxionMass class

The mimes::AxionMass<LD> class is responsible for the definition of the axion mass. The header file of this class is \mimes/src/AxionMass/AxionMass.hpp.

    template< class LD >
    mimes::AxionMass(std::string chi_PATH, LD minT=0, LD maxT=mimes::Cosmo::mP)

The first argument, chi_PATH, is the path to a data file that contains two columns; $T$ (in GeV) and $\chi$ (in $GeV^4$), with increasing $T$. The arguments minT and maxT are the interpolation limits. These limits are used in order to stop the interpolation in the closest temperatures that exist in the data file. That is the actual interpolation limits are $T_{\min}\geq$minT and $T_{\max}\leq$maxT. Beyond these limits, by default, the axion mass is assumed to be constant. However, this can be changed by using the member functions

    void set_ma2_MIN(std::function<LD(LD,LD)> ma2_MIN)
    void set_ma2_MAX(std::function<LD(LD,LD)> ma2_MAX)

Here, ma2_MIN and ma2_MAX are functors that define the axion mass squared beyond the interpolation limits. In order to ensure that the axion mass is continuous, usually we need $T_{\min}$, $T_{\max}$, $\chi(T_{ min})$, and $\chi(T_{ max})$. These values can be obtained using the member functions

1. template< class LD > LD mimes::AxionMass::getTMin(): This function returns the minimum interpolation temperature, $T_{ min}$.
2. template< class LD > LD mimes::AxionMass::getTMax(): This function returns the maximum interpolation temperature, $T_{ max}$.
3. template< class LD > LD mimes::AxionMass::getChiMin(): This function returns $\chi(T_{ min})$.
4. template< class LD > LD mimes::AxionMass::getChiMax(): This function returns $\chi(T_{ max})$.

    template< class LD >
    mimes::AxionMass(std::function<LD(LD,LD)> ma2)

Here, the only argument is the axion mass squared, $\tilde{m}_a(T)$, defined as a callable object.

    template< class LD >    LD mimes::AxionMass::ma2(LD T, LD fa)

We should note that ma2 is a public std::function<LD(LD,LD)> member variable. Therefore, it can be assigned using the assignment operator. However, in order to change its definition, we can also use the following member function:

    template< class LD > void mimes::AxionMass::set_ma2(std::function<LD(LD,LD)> ma2)

The Axion class

    template< class LD, const int Solver, class Method >
    mimes::Axion< LD, Solver, Method >(LD theta_i, LD fa, LD umax, LD TSTOP, 
                    LD ratio_ini, unsigned int N_convergence_max, LD convergence_lim, 
                    std::string inputFile, AxionMass *axionMass, LD initial_step_size=1e-2, 
                    LD minimum_step_size=1e-8, LD maximum_step_size=1e-2, LD absolute_tolerance=1e-8, 
                    LD relative_tolerance=1e-8, LD beta=0.9, LD fac_max=1.2, LD fac_min=0.8, 
                    unsigned int maximum_No_steps=10000000)

1. theta_i: Initial angle, $\theta_i$.
2. fa: The PQ scale.
3. umax : If $u>$ umax the integration stops (remember that $u=\log(a/a_i)$). Typically, this should be a large number ($\sim 1000$), in order to avoid stopping the integration before the axion begins to evolve adiabatically.
4. TSTOP: If the temperature drops below this, integration stops. In most cases this should be around $10^{-4}$ GeV, in order to be sure that any entropy injection has stopped before integration stops should not be violated).
5. ratio_ini: Integration starts when $3H/ m_a(T) $ ratio_ini (the exact point depends on the file `` inputFile", which we will see later).
6. N_convergence_max and convergence_lim: Integration stops when the relative difference between two consecutive peaks is less than convergence_lim for N_convergence_max consecutive peaks. This is the point beyond which adiabatic evolution is assumed.
7. inputFile : Relative (or absolute) path to a file that describes the cosmology. the columns should be: $u$ $T$ (in GeV) $\log H$, sorted so that $u$ increases.² It is important to remember that MiMeS assumes that the entropy injection has stopped before the lowest temperature given in inputFile. Since MiMeS is unable to guess the cosmology beyond what is given in this file, the user has to make sure that there are data between the initial temperature (which corresponds to ratio_ini), and TSTOP.
8. axionMass: An instance of the mimes::AxionMass class, passed by pointer.
9. initial_stepsize (optional): Initial step the solver takes.
10. maximum_stepsize (optional): This limits the step-size to an upper limit.
11. minimum_stepsize (optional): This limits the step-size to a lower limit.
12. absolute_tolerance (optional): Absolute tolerance of the RK solver.
13. relative_tolerance (optional): Relative tolerance of the RK solver. Generally, both absolute and relative tolerances should be $10^{-8}$. In some cases, however, one may need more accurate result (if fa is extremely high, the oscillations happen violently, and the system destabilizes). In any case, if the tolerances are below $10^{-8}$, LD should be long double. MiMeS by default uses long double variables, in order to change it see the options available in Compile-time choices.
14. beta (optional): Controls how agreesive the adaptation is. Generally, it should be around but less than 1.
15. fac_max, fac_min (optional): The stepsize does not increase more than fac_max, and less than fac_min. This ensures a better stability. Ideally, fac_max$=\infty$ and fac_min$=0$, but in reality one must tweak them in order to avoid instabilities.
16. maximum_No_steps (optional): Maximum steps the solver can take. Quits if this number is reached even if integration is not finished.

    template< class LD, const int Solver, class Method >
    void mimes::Axion< LD, Solver, Method >::solveAxion()

Once this function finishes, the results are stored in several member variables.

    template< class LD, const int Solver, class Method > 
    std::vector< std::vector< LD > > mimes::Axion< LD, Solver, Method >::points

    template< class LD, const int Solver, class Method > 
    std::vector< std::vector< LD > > mimes::Axion< LD, Solver, Method >::peaks

Note that these points are computed using linear interpolation between two integration points with a change in the sign of $\zeta$.

    template< class LD, const int Solver, class Method > 
    std::vector< LD > mimes::Axion< LD, Solver, Method >::dtheta
    
    template< class LD, const int Solver, class Method > 
    std::vector< LD > mimes::Axion< LD, Solver, Method >::dzeta

    template< class LD, const int Solver, class Method >
    LD mimes::Axion< LD, Solver, Method >::T_osc
    
    template< class LD, const int Solver, class Method >
    LD mimes::Axion< LD, Solver, Method >::a_osc

    template< class LD, const int Solver, class Method >
    LD mimes::Axion< LD, Solver, Method >::theta_osc

    template
    LD mimes::Axion< LD, Solver, Method >::gamma

    template< class LD, const int Solver, class Method >
    LD mimes::Axion< LD, Solver, Method >::relic

    template< class LD, const int Solver, class Method >
    void mimes::Axion< LD, Solver, Method >::setTheta_i(LD theta_i)

We should note that running this function all variables are cleared. So we lose all information about the last time axionSolve() ran.

    template< class LD, const int Solver, class Method >
    void mimes::Axion< LD, Solver, Method >::restart()

Again, this function clears all member variables. So it should be used with caution.

    template< class LD, const int Solver, class Method >
    static mimes::Cosmo< LD > mimes::Axion< LD, Solver, Method >::plasma

This variable can be used without an instance of the mimes::Axion< LD,Solver,Method > class.