CHEIDI: HEIDI.fr

(* *************************************************************************************
                        HEIDI Higgs with a compact extra dimension.

   The model is intended for use on top of the Standard Model where it replaces the
   Higgs sector. As the theory contains a variable number of modes below the cutoff
   (the masses of which cannot be expressed analytically), the HEIDI parameters have to
   be set before creating the Feynman rules and model files. Those parameters are:

   Heidi$v      : higgs VeV. Optional, default is the value.
   Heidi$mh     : higgs Mass. Optional, default is the SM value
   Heidi$cs     : compactification scale
   Heidi$mb     : bulk mass
   Heidi$g2     : 5D trilinear coupling squared
   Heidi$cutoff : cutoff scale
   Heidi$nmodes : if Heidi$cutoff is undefined, Heidi$nmodes determines the number of
                  modes in the spectrum (starting from index 0)

   Masses, widths and couplings are calculated statically. DO NOT CHANGE THOSE OR THE
   ELECTROWEAK PARAMETERS AFTER CREATING THE MODEL.

   This model file is a descendant of Claude Duhr's SMScalars model.

*)


(* Make sure everything we need is at our disposal *)
Needs["Heidi`"];
Block[{error},
   error[s_] := (Print["HEIDI error: please set " <> s]; Throw[Null]);
   If[!NumericQ[Heidi$v], error["Heidi$v"]];
   If[!NumericQ[Heidi$mh], error["Heidi$mh"]];
   If[!NumericQ[Heidi$cs], error["Heidi$cs"]];
   If[!NumericQ[Heidi$mb], error["Heidi$mb"]];
   If[!NumericQ[Heidi$g2], error["Heidi$g2"]];
   If[FeynmanGauge,
      Print["HEIDI error: the model currently only supports unitarity gauge"
         <> "; please set FeynmanGauge = False."];
      Throw[Null];
   ];
];

(* Init HEIDI *)
InitHeidi[Heidi$v, Heidi$mh, Heidi$cs, Heidi$mb, Heidi$g2];

(* Determine the number of modes we have been requested to include *)
If[NumericQ[Heidi$cutoff],
   Heidi$nmodes = FixedPoint[If[HeidiMass[#] < Heidi$cutoff, # + 1, #]&, 0];
];
If[!NumericQ[Heidi$nmodes],
   Print["HEIDI error: please set either Heidi$cutoff or Heidi$nmodes"];
   Throw[Null];
];

(* Model metadata *)
M$ModelName = "HEIDI";
M$Information = {Authors -> "C. Speckner",
   Institutions -> "Universitaet Freiburg",
   Emails -> "Christian.Speckner@physik.uni-freiburg.de",
   Date -> "30.07.2010",
   Version  ->  "0.1",
   URLs  -> "None",
   References -> (
      ""
         <> "\"Exploring the golden channel for HEIDI models using an "
         <> "interface between WHIZARD and FeynRules\", "
         <> "N. D. Christensen, C. Duhr, B. Fuks, J. Reuter, C. Speckner, "
         <> "arXiv:1010.3251")
};

(* The flavor index for the HEIDI modes *)
IndexRange[Index[Heidi]] = Range[Heidi$nmodes];
IndexStyle[Heidi, k];

M$Parameters = {
(* The Higgs wavefunction is a tensor parameter, the entries of which are *
 * stored numerically as external parameters.                             *)
   xi == {
      Indices -> {Index[Heidi]},
      ComplexParameter -> False,
      Value -> Table[Rule[xi[n], HeidiWavefunction[n - 1, 0]],
         {n, 1, Heidi$nmodes}],
      Description -> "Higgs wave function",
      ParameterType -> External,
      AllowSummation -> True
   },
(* K factor for the effective ggh coupling *)
   kggh == {
      Value -> 1,
      ParamterType -> External,
      Description -> "ggh K factor"
   }
};


(* Enter the HEIDI scalars *)
M$ClassesDescription =  {
   S[4] == {
      ClassName -> Hh,
      ClassMembers -> Table[Symbol["Hh" <> ToString[n]], {n, 1, Heidi$nmodes}],
      SelfConjugate -> True,
      Indices -> {Index[Heidi]},
      FlavorIndex -> Heidi,
      Mass -> Table[{Symbol["mhh" <> ToString[n]], HeidiMass[n - 1]}, {n, 1, Heidi$nmodes}],
      Width -> Table[{Symbol["whh" <> ToString[n]], HeidiWidth[n - 1]}, {n, 1, Heidi$nmodes}],
      PDG -> Table[990000 + n, {n, 1, Heidi$nmodes}]
  }
}

(* Construct HEIDI Lagrangian *) 
LHEIDI := (LSM /. {muH -> 0, \[Lambda] -> 0}) /. H -> (Hh[n]*xi[n]);
LHEIDI::usage = "LHEIDI: HEIDI lagrangian";

(* Optional: the effective ggh_i lagrangian *)
LHEIDIgg[mode_] := Block[{fun, piece, tt, mm, hh},
   Print["adding effective scalar-gluon-gluon vertices " <> If[mode == "heavytop",
      "in heavy top limit", "with dynamical loop factor"]];
   tt[m_] := 4*MT^2 / m^2;
   fun[tau_] := If[mode == "heavytop", 2/3,
      tau * (1 + (1 - tau)* If[NumericalValue[tau] > 1,
         ArcSin[Sqrt[1 / tau]]^2
      ,
         -1/4 * (Log[(1 + Sqrt[1 + tau])/(1 - Sqrt[1 - tau])]
            - I*Pi)^2
      ])
   ];
   mm[n_] := Symbol["mhh" <> ToString[n]];
   hh[n_] := Symbol["Hh" <> ToString[n]];
   piece[n_] :=
      Sqrt[kggh] * gs^2/32/Pi^2/v * fun[tt[mm[n]]] * xi[n] *
      hh[n] * (del[G[mu, a], nu] - del[G[nu, a], mu])^2;
   Plus @@ (piece /@ Range[Heidi$nmodes])
];
LHEIDIgg[] := LHEIDIgg[""];
LHEIDIgg::usage = "LHEIDIgg [mode]: HEIDI effective scalar-gluon-gluon couplings
   mode : optional, set to \"heavytop\" to use heavy top limit";

(* Protect our symbols. *)
Protect[LHEIDI, LHEIDIgg];

(* Spit out final warnings. *)
Print["HEIDI initialized. Please doublecheck that the runtime parametrization of "
   <> "the electroweak sector is consistent with the parameters you used for HEIDI "
   <> "before running any code."];