# Models with Thresholds in SPheno

## General

In the SPheno Output of SARAH it is possible to include thresholds in the RGE evaluation at which heavy degrees of freedoms are removed from the spectrum. The implementation depends if the thresholds are also associated with the breaking of an additional gauge symmetry or now.

## Thresholds without gauge symmetry breaking

### Assumptions

If all scales have the same gauge structure, it is possible for SARAH to derive the RGEs for all scales from the RGEs for the highest scale by performing the following steps:

• For those fields which should be integrated out during the run new variables $n_{gen}(\Phi_i)$ are introduced, which define the number of generation of the heavy field $\Phi_i$. All gauge group constants like the the Dynkin index summed over chiral superfields, $S(R)$, are expressed as function of $n_{gen}(\Phi_i)$. These $n_{gen}(\Phi_i)$ are dynamically adjusted, when the energy scale crosses a threshold.
• It is also necessary to set the couplings which involve heavy fields to zero when a threshold is crossed. For example, the Yukawa type coupling of the form$Y^{ij} \Phi_i \phi_j H$ involves three generations of the heavy field$\Phi$. At the threshold of$\Phi_k$, the$k$-th row of$Y$ is set to zero. That happens similarly for all other superpotential and soft-breaking parameters.
• It is assumed, that the masses of the scalar and fermionic component of a heavy superfield are the same, i.e. the masses are much larger than the soft-breaking masses. Furthermore, it is assumed that the masses are given by a bilinear superpotential term.

SARAH included automatically the threshold corrections to gauge couplings and gaugino masses in the SPheno output.

### Implementation

To include thresholds without gauge symmetry breaking, the following steps have to be performed:

1. The heavy fields must be deleted in the SARAH model definition:

 DeleteFields = {...};

This ensures, that the these particle are not take into account for the calculation of decays or loop corrections at the SUSY scale.

2. The thresholds have to be defined in SPheno.m :

 Thresholds = {{Scale1, {HeavyFields1}},
{Scale2, ... }};

For all scales an entry in the array Thresholds has to be added. Each entry defines first the threshold scale, at second position a list with the heavy superfields is given. Also specific generations for a superfield can be given.

It is possible to define boundary conditions at each threshold scale for running up and down separately:

 BoundaryConditionsUp[[x]] = { ...};
BoundaryConditionsDown[[x]] = { ...};
##### Example

The simplest example for a SUSY model with thresholds is the seesaw I version of the MSSM. The changes compared to the MSSM are as follows:

1. In the model file, the right handed neutrino superfield is defined and additional superpotential couplings are added
SuperFields[[8]] = {v, 3, vR,          0, 1,  1, RpM};
SuperPotential =  ... + Yv v.l.Hu + Mv/2 v.v +WOp/2 l.Hu.l.Hu;


Note the presence of the Weinberg operator WOp

2. The right-handed fields are deleted from the SUSY scale Lagrangian
DeleteParticles={v};

3. The thresholds in SPheno.m are defined
Thresholds={
{Abs[MvIN[1,1]],{v[1]}},
{Abs[MvIN[2,2]],{v[2]}},
{Abs[MvIN[3,3]],{v[3]}}
};

4. When running down, boundary conditions to calculate the shift in the Weinberg operator are added for each threshold scale
BoundaryConditionsUp=Table[{},{Length[Thresholds]}];
BoundaryConditionsDown=Table[{},{Length[Thresholds]}];

BoundaryConditionsDown[[1]]={ {WOp[index1,index2],WOp[index1,index2] - Yv[1,index1] Yv[1,index2]/MassOfv[1]} };
BoundaryConditionsDown[[2]]={ {WOp[index1,index2],WOp[index1,index2] - Yv[2,index1] Yv[2,index2]/MassOfv[2]} };
BoundaryConditionsDown[[3]]={ {WOp[index1,index2], - Yv[3,index1] Yv[3,index2]/MassOfv[3]} };


## Thresholds with gauge symmetry breaking

If the gauge structure at the different scales are different, each set of RGEs is calculated separately and this information is then combined into one consistent version of SPheno . This code includes routines for calculating finite shifts in the gauge couplings and gaugino mass parameters. As an example or for instance$SO(10)$ motivated left-right models.

### Implementation

The implementation of models with a threshold scale that involved a gauge symmetry breaking is more involved. In general, the following steps are necessary:

1. For each regime, a separate model file for SARAH has to be created. These model file have to be saved in the subdirectories Regime-1, Regime-2, …of

[SARAH Directory]/Models/[Model]/

They are numbered beginning with the highest scale.

2. The SPheno input file for the higher scales must provide the following information:

• IntermediateScale = True

• RegimeNr = X

• A list of the heavy fields, which should be integrated out, the gauge sector below the threshold as well as the corresponding quantum numbers of the fields integrated. That’s needed to calculate the finite shifts at the threshold scale. For instance, the entries might read

HeavyFields = {Field_1, Field_2,..};
NextGauge =   {U[1],   SU[2],   SU[2],    SU[3]};
NextQN = {
{Field_1, 0,       2,       1,         1},
{Field_2, 1/3,     1,       2,         4},
...
};
3. All necessary information for combining the regimes to one SPheno is given in SPheno.m of the lowest scale.

• IntermediateScale = False

• RegimeNr = X

• The threshold scales: ThresholdScales = ...

• The boundary conditions for running up and down at each threshold scale:

 BoundaryConditionsUp[[x]] = {...};
BoundaryConditionsDown[[x]] = {...};

In the boundary conditions index1, index2, …can be used for defining sums over indices.

• The usual information for SPheno , defined in the sec. [sec:sphenoinputfile].

When starting the SPheno output of the lowest scale, automatically all other scales are evolved. Note, to calculate the RGEs of the different regimes requires SARAH to start one additional Mathematica kernel. For passing the information between the different Mathematica kernels a directory Meta in the model directory is created by SARAH . Also the screen output of Mathematica during the evaluation of the higher regimes is written to that directory (Output-Regime-X.m). So, the user can supervise the progress and see potential error messages. The necessary information of each regime for writing the combined source code for SPheno at the end is saved by SARAH in the files Regime-X.m.

### Example

The entire discussion in the following is based on the model Omega-Short delivered with SARAH. This model has a threshold scale [/itex]SU(2)_L \times SU(2)_R \times U(1)_{B-L} \rightarrow SU(2)_L \times U(1)_Y[/itex]

1. Boundary conditions: We don't further specify here all incredients of the model but just listen some boundary conditions which are important when matching the model to the MSSM.
[/itex]\begin{aligned} \label{eq:boundary_LR_1} Y_d = Y_Q^1 \cos \theta_1 - Y_Q^2 \sin \theta_1 \thickspace, &\qquad& Y_u = - Y_Q^1 \cos \theta_2 + Y_Q^2 \sin \theta_2 \thickspace, \\ Y_e = Y_L^1 \cos \theta_1 - Y_L^2 \sin \theta_1 \thickspace, &\qquad& Y_\nu = - Y_L^1 \cos \theta_2 + Y_L^2 \sin \theta_2 \thickspace,\end{aligned}[/itex] <br: where$R = \sin (\theta_1 - \theta_2)$. For the soft-trilinear couplings,$Y$ is replaced with$T$ in the expressions. For the sfermionic soft masses, we have
[/itex]\begin{aligned} m_{q}^2 =m_{u^c}^2 = m_{d^c}^2 &=& m_{Q^c}^2 \thickspace,\\ m_{l}^2=m_{e^c}^2 &=& m_{L^c}^2 \thickspace, \\ M_L = M_R &=& M_2 \thickspace.\end{aligned}[/itex]
while in the Higgs sector we need the relations
[/itex]\begin{aligned} m_{H_d}^2 &=& \cos^2 \theta_1 (m_\Phi^2)_{11} + \sin^2 \theta_1 (m_\Phi^2)_{22} - \sin \theta_1 \cos \theta_1 \left[ (m_\Phi^2)_{12} + (m_\Phi^2)_{21} \right] \thickspace, \\ m_{H_u}^2 &=& \cos^2 \theta_2 (m_\Phi^2)_{11} + \sin^2 \theta_2 (m_\Phi^2)_{22} - \sin \theta_2 \cos \theta_2 \left[ (m_\Phi^2)_{12} + (m_\Phi^2)_{21} \right] \thickspace, \end{aligned}[/itex]
In the gauge sector, we have to express the hypercharge coupling and the corresponding gaugino through
[/itex]\begin{aligned} g_1 & = & \frac{\sqrt{5} g_2 g_{BL}}{\sqrt{2 g_2^2 + 3 g_{BL}^2}} \thickspace , \\ \label{eq:boundary_LR_2} M_1 & = & \frac{2 g_2^2 M_{BL} + 3 g_{BL}^2 M_R}{2 g_2^2 + 3 g_{BL}^2} \thickspace.\end{aligned}[/itex]
2. Defining the SARAH models: one needs two model files
• Above the breaking scale the full model
• Below the breaking scale, only the particles and gauge groups of the MSSM survive as dynamical degrees of freedom. Therefore, the model file is almost identical that for the MSSM with only tiny differences:
• The superpotential contains the Weinberg operator

SuperPotential =  ... +  WOp/2 l.Hu.l.Hu;
• The neutrinos are massive and mix among each other

DEFINITION[EWSB][MatterSector]= { ...,{{FvL}, {FV, UV}}};
3. Input file for SPheno
• Above$SU(2)_R \times U(1)_{B-L}$ breaking scale: First, the regime must be flagged as an intermediate scale, and SARAHmust be told its position in the scale hierarchy (counted from GUT to low scale).
RegimeNr = 1;
IntermediateScale = True; 

Afterwards, we give a list with all particles which are integrated out at the threshold scale after gauge symmetry breaking.

HeavyFields = {Hpm1R1, ChiR1, Cha1r1, hhR1, AhR1,
FvR1, SVRr1, SH0r1[3], SHCr1[3],
SO1r1, SDLpR1, SDLppR1,
SDL0r1, SDRmmR1, DR3r1,
DL1r1, DL2r1, DL3r1, H0r1, HCr1};

The numbers in square brackets indicate that only the third generation and above is integrated out.

In order to calculate the finite shifts of the gauge couplings and gaugino masses, it is necessary to define the gauge sector of the next scale NextGauge as well as the quantum number of the fields which are integrated out with respect to those gauge groups.

NextGauge=           {U[1], SU[2], SU[3]};
NextQN = { {Hpm1R1,     -1,     1,    1},
{ChiR1,       0,     1,    1},
{Cha1r1,     -1,     1,    1},
{hhR1,        0,     1,    1},
{AhR1,        0,     1,    1},
{FvR1,        0,     1,    1},
{SVRr1,       0,     1,    1},
{SH0r1,    -1/2,     1,    1},
{SHCr1,     1/2,     2,    1},
{SO1r1,       0,     1,    1},
{SDLpR1,      1,     1,    1},
{SDLppR1,     2,     1,    1},
{SDL0r1,      1,     3,    1},
{SDRmmR1,    -2,     1,    1},
{DR3r1,      -2,     1,    1},
{DL1r1,       1,     1,    1},
{DL2r1,       1,     1,    1},
{DL3r1,       1,     3,    1},
{H0r1,     -1/2,     1,    1},
{HCr1,      1/2,     2,    1}
};

Finally, SARAHneeds information on the vacuum conditions. There are two different VEVs, and therefore, we need to choose two parameters which are fixed by the tadpole equations. As there is no closed analytical solution for them, we give an approximation obtained by neglecting the soft-breaking terms.

ParametersToSolveTadpoles = {Mdelta, Momega};

Mdelta -> - SignumMdelta ac1 vR/Sqrt[2],
Momega -> - SignumMomega ac1 vBL^2/(2 Sqrt[2] vR)
};
SubSolutionsTadpolesLoop={}; 
• Below$SU(2)_R \times U(1)_{B-L}$ breaking scale: The second scale is not an intermediate scale, and hence
RegimeNr = 2;
IntermediateScale = False; 

We make the following choice of free parameters of that model: to the set of standard mSugra parameters ([/itex]m_0,M_{1/2},A_0,\tan\beta,\text{sign}\mu[/itex]), we add$B_0$, the superpotential parameter$a$, the signs of$M_\Omega$ and$M_\Delta$, the two VEVs$v_R$ and$v_{BL}$, and the threshold scale. These are defined in the blocks MINPAR and EXTPAR.

MINPAR= {
{1,  m0},
{2,  m12},
{3,  TanBeta},
{4,  SignumMu},
{5,  Azero},
{6,  Bzero},
{7,  SignumMomega},
{8,  SignumMdelta},
{9,  aInput}};

EXTPAR = {
{100, vRinput},
{101, vBLinput},
{200, TScale}};

As in the, MSSM we fix the numerical values of$\mu$ and$B_\mu$ by the solving the two tadpole equations,

ParametersToSolveTadpoles = {\[Mu],B[\[Mu]]};

Furthermore, we use also the common definitions for the SUSY scale, see also sec. [sec:ExampleSPhenoMSSM].

RenormalizationScaleFirstGuess = m0^2 + 4 m12^2;
RenormalizationScale = MSu[1]*MSu[6];

We use as condition for the GUT scale

ConditionGUTscale = {gBL==g2, g1 == g2};

The second entry is necessary because it might be that the masses of the fields breaking the left-right symmetry are above the unification scale which has to be handled seperatly. At the GUT scale we use boundary conditions motivated by minimal supergravity: all scalar soft-breaking masses are proportional to$m_0$, the gaugino masses are proportional to$M_{1/2}$, the trilinear soft-breaking couplings are given by the corresponding superpotential parameter times$A_0$, and the bilinear soft-breaking couplings are set to$B_0$ times the superpotential parameter. The values for the coupling matrices$f$,$\alpha$ as well as for$\mu$ are read from the LesHouches input file.

BoundaryHighScale={
{T[YQ],       Azero*YQ},
{T[YL],       Azero*YL},
{f,           LHInput[f]},
{T[f],        Azero*LHInput[fm]},
{AlphaOm,     LHInput[AlphaOm]},
{T[AlphaOm],  Azero*LHInput[AlphaOm]},
{T[a],        Azero*aInput},
{B[Mdelta],   Bzero*Mdelta},
{B[Momega],   Bzero*Momega},
{B[Mu3],      Bzero*LHInput[Mu3]},
{mqL2,        DIAGONAL m0^2},
{mqR2,        DIAGONAL m0^2},
{mlL2,        DIAGONAL m0^2},
{mlR2,        DIAGONAL m0^2},
{mPhi2,       DIAGONAL m0^2},
{mdeltaL2,    m0^2},
{mdeltaLbar2, m0^2},
{mdeltaR2,    m0^2},
{mdeltaRbar2, m0^2},
{momegaL2,    m0^2},
{momegaR2,    m0^2},
{MassB,       m12},
{MassWL,      m12},
{MassG,       m12}
}; 

To glue the both regimes, we need to define the appropriate boundary conditions. First, we initialize the arrays

ThresholdScales = {TSCALE};

BoundaryConditionsUp = Table[{},{Length[ThresholdScales]}];
BoundaryConditionsDown = Table[{},{Length[ThresholdScales]}]; 

and then encode the equations eqs. ([eq:boundaryLR1])-([eq:boundaryLR2]). In order to keep the code short, we define

  ST1  = Sin[Theta1];
CT1  = Cos[Theta1];
ST2  = Sin[Theta2];
CT2  = Cos[Theta2];
ST21 = Sin[Theta2-Theta1];
CT21 = Cos[Theta2-Theta1];

Using these abbreviations, the boundary conditions can be written as

BoundaryConditionsUp[[1]] = {
{YQ[index1,index2,1],   (Yu[index1,index2] ST1 + Yd[index1,index2]ST2)/ST21 },
{YQ[index1,index2,2],   (Yu[index1,index2] CT1 + Yd[index1,index2]CT2)/ST21 },
{YL[index1,index2,1],   (Yv[index1,index2] ST1 + Ye[index1,index2]ST2)/ST21 },
{YL[index1,index2,2],   (Yv[index1,index2] CT1 + Ye[index1,index2]CT2)/ST21 },
{gBL,                    Sqrt[2] g1 g2 /Sqrt[5 g2^2 -3 g1^2]},
{Yv,                      LHInput[Yv]}
};

BoundaryConditionsDown[[1]] = {
{vR,           vRinput},
{vBL,          vBLinput},
{a,            aInput},
{Theta1,       ArcTan[RealPart[((vR*AlphaOm[1,2])/2 + Mu3[1,2])/Mu3[2,2]]]},
{Theta2,       ArcTan[RealPart[(-(vR*AlphaOm[1,2])/2 + Mu3[1,2])/Mu3[2,2]]]},
{g1,           Sqrt[5] g2 gBL/Sqrt[2 g2^2 + 3 gBL^2]},
{MassB,        (2 g2^2 MassB + 3 gBL^2 MassWL)/(2 g2^2 + 3 gBL^2)},
{MassWB,       MassWL},
{Yd[index1,index2], YQ[index1,index2,1] CT1 - YQ[index1,index2,2] ST1},
{Yu[index1,index2], - YQ[index1,index2,1] CT2 + YQ[index1,index2,2] ST2},
{Ye[index1,index2], YL[index1,index2,1] CT1 - YL[index1,index2,2] ST1},
{Yv[index1,index2], - YL[index1,index2,1] CT2 + YL[index1,index2,2] ST2},
{T[Yd][index1,index2], T[YQ][index1,index2,1] CT1 -T[YQ][index1,index2,2]ST1},
{T[Yu][index1,index2], -T[YQ][index1,index2,1] CT2+T[YQ][index1,index2,2] ST2},
{T[Ye][index1,index2], T[YL][index1,index2,1] CT1 -T[YL][index1,index2,2] ST1},
{mu2,          mqR2},
{md2,          mqR2},
{mq2,          mqR2},
{me2,          mlR2},
{ml2,          mlR2},
{mHd2, CT1^2 mPhi2[1,1] + ST1^2 mPhi2[2,2] - ST1 CT1(mPhi2[1,2] + mPhi2[2,1])},
{mHu2, CT2^2 mPhi2[1,1] + ST2^2 mPhi2[2,2] - ST2 CT2(mPhi2[1,2] + mPhi2[2,1])},
{WOp,          MatMul2[MatMul2[Yv,InverseMatrix[f],FortranFalse],Transpose[Yv],
FortranFalse]/vR}
};

Note that$Y_\nu$ does not appear in any of the model files and it is therefore necessary to fix the dimension of that matrix by hand

AdditionalVariablesSPheno={Yv[3,3]};

Several parameters are restricted to be real. In addition, it is helpful to choose initialization values for some parameters to stabilize the numerics in the first iteration

RealParameters = {TanBeta, vRinput,vBLinput,Theta1,Theta2,TScale};

InitializationValues = {
{Mu3IN[1,1], (Mu3IN[1,2]^2 - AlphaOmIN[1,2]^2 vRInput^2/4)/Mu3IN[2,2]},
{Theta1,   ArcTan[RealPart[-(Mu3IN[1,2]+AlphaOmIN[1,2]vRInput/2)/Mu3IN[2,2]]]},
{Theta2,   ArcTan[RealPart[(Mu3IN[1,2]-AlphaOmIN[1,2]vRInput/2)/Mu3IN[2,2]]]},
{Mdelta,   aInput*SignumMdelta*vRinput/2 },
{Momega,   SignumMomega*(aInput^2*vBLinput^2)/(8 Mdelta)}
}; 
4. Getting the SPheno code: One needs to start the MSSM part of the model and run MakeSPheno. SARAH will start automatically a second Mathematica kernel to generate all necessary information for the left-right symmetric phase.