Skip to content

File scalefactorinitializer.h

File List > code_source > cosmolattice > include > CosmoInterface > initializers > scalefactorinitializer.h

Go to the documentation of this file

#ifndef COSMOINTERFACE_INITIALIZERS_SCALEFACTORINITIALIZER_H
#define COSMOINTERFACE_INITIALIZERS_SCALEFACTORINITIALIZER_H

/* This file is part of CosmoLattice, available at www.cosmolattice.net .
   Copyright Daniel G. Figueroa, Adrien Florio, Francisco Torrenti and Wessel Valkenburg.
   Released under the MIT license, see LICENSE.md. */

// File info: Main contributor(s): Daniel G. Figueroa, Adrien Florio, Francisco Torrenti,  Year: 2020

#include "CosmoInterface/common.h"
#include "TempLat/util/rangeiteration/for_in_range.h"
#include "TempLat/lattice/algebra/operators/power.h"
#include "CosmoInterface/evolvers/evolvertype.h"
#include "CosmoInterface/runparameters.h"
#include "CosmoInterface/definitions/energies.h"
#include "CosmoInterface/definitions/nonminimalcoupling.h"

namespace TempLat
{
  class ScaleFactorInitializer
  {
  public:
    // Put public methods here. These should change very little over time.
    ScaleFactorInitializer() = delete;

    // @label:scalefactorinitializer_class
    template <class Model, typename T> static void initializeScaleFactor(Model &model, RunParameters<T> &rPar)
    {
      //For model swith defects, only a fixed background is allowed by default
      if (!rPar.fixedBackground && Model::DefectsModel)
        throw(RunParametersInconsistent("Running a defects model with self-consistent expansion is not tested, and features such as (extra)fattening may not work correctly. If you really want to run this option, comment out this exception in scalefactorinitializar.h"));

      // If fixed background, the initial Hubble parameter H0 is given by the user
      if (rPar.fixedBackground) {
        auto fbe = FixedBackgroundExpansion(model, rPar);
        model.aI = fbe(0.);
        model.aDotI = fbe.dot(0.);
        if constexpr (Model::IsNonMinimallyCoupled) model.RI = fbe.R(0.);
        // H0 is in GeV, so we transform to program variables
      }

      // If self-consistent expansion, the initial Hubble parameter must be computed
      // via the 1st Friedmann eqn:
      else {
        T kin = 0; // initial kinetic energy
        ForLoop(i, 0, Model::Ns - 1,
                kin += Energies::kineticS(model, pow<2>(model.piS0(i) / model.fStar / model.omegaStar)););
        ForLoop(i, 0, Model::NCs - 1,
                kin += Energies::kineticCS(model, norm2(model.piCS0(i) / model.fStar / model.omegaStar)););
        ForLoop(i, 0, Model::NSU2Doublet - 1,
                kin +=
                Energies::kineticSU2Doublet(model, norm2(model.piSU2Doublet0(i) / model.fStar / model.omegaStar)););
        // Note: piS0(i), piCS0(i), piSU2Doublet0(i) are given in GeV^2,
        // so we divide by fStar*omegaStar to transform to program variables

        auto rhoMC = kin + model.pot0;

        model.aI = rPar.a0;
        model.aDotI = model.fStar / Model::MPl * sqrt(rhoMC / 3.0); // 1st Friedmann eqn

        if constexpr (Model::IsNonMinimallyCoupled) {
          auto A = 1.0 - Total(i, 0, Model::Ns - 1,
              IfElse(Model::NonMinimalCouplings::couples(i, 0_c),
                pow<2>(1.0 / Model::MPl) * model.xis(i, 0_c) * pow<2>(model.fldS0(i)), 0));
          auto B = -2.0 * Total(i, 0, Model::Ns - 1,
              IfElse(Model::NonMinimalCouplings::couples(i, 0_c),
                pow<2>(1.0 / Model::MPl) / model.omegaStar * model.xis(i, 0_c) * model.fldS0(i) * model.piS0(i), 0));
          auto C = -pow<2>(model.fStar / Model::MPl) / 3.0 * rhoMC;

          auto Delta = pow<2>(B) - 4 * A * C;
          auto H1 = (-B + sqrt(Delta)) / 2.0 / A;

          model.aI = rPar.a0;
          model.aDotI = H1;
          model.piAI = model.aDotI;
          model.RI = NonMinimalCoupling::R(model);
        }

        // Note: Initially, the gradients are 0, so the scale factors receives
        //  only contributions from the kinetic energies and the potential.
      }
    }
    // @endlabel
  };
} // namespace TempLat

#endif