FlexibleSUSY is hosted by Hepforge, IPPP Durham
HSSUSY

HSSUSY

HSSUSY (high scale supersymmetry) is an implementation of the Standard Model, matched to the MSSM at the SUSY scale, $M_\text{SUSY}$. The setup of HSSUSY is shown in the following figure.

HSSUSY_tower.svg

Boundary conditions

__High scale__

In HSSUSY, the `HighScale` variable is set to the SUSY scale, $M_{\text{SUSY}}$. At this scale the quartic Higgs coupling, $\lambda(M_\text{SUSY})$, is predicted from the matching to the MSSM using the full 1-loop and dominant 2-loop threshold corrections of $O((\alpha_t + \alpha_b)\alpha_s + (\alpha_t + \alpha_b)^2 + \alpha_b\alpha_\tau + \alpha_\tau^2)$ from [arXiv:1407.4081, arXiv:1504.05200, arXiv:1703.08166].

The 3- and partial 4-loop renormalization group equations of [arXiv:1303.4364, arXiv:1307.3536, arXiv:1508.00912, arXiv:1508.02680, arXiv:1604.00853] are used to run $\lambda(M_\text{SUSY})$ down to the electroweak scale $M_Z$ or $M_{\text{EWSB}}$.

If $M_{\text{SUSY}}$ is set to zero, $M_{\text{SUSY}} = \sqrt{m_{\tilde{t}_1}m_{\tilde{t}_2}}$ is used.

__Low scale__

The `LowScale` is set to $M_Z$. At this scale, the $\overline{\text{MS}}$ gauge and Yukawa couplings $g_{1,2,3}(M_Z)$, $Y_{u,d,e}(M_Z)$, as well as the SM vacuum expectation value (VEV), $v(M_Z)$, are calculated at the full 1-loop level from the known low-energy couplings $\alpha_{\text{em}}^{\text{SM(5)}}(M_Z)$, $\alpha_s^{\text{SM(5)}}(M_Z)$, from the pole masses $M_Z$, $M_e$, $M_\mu$, $M_\tau$, $M_t$ as well as from the $\overline{\text{MS}}$ masses $m_b^{\text{SM(5)}}(m_b)$, $m_c^{\text{SM(4)}}(m_c)$, $m_s(2\,\text{GeV})$, $m_d(2\,\text{GeV})$, $m_u(2\,\text{GeV})$. In addition to these 1-loop corrections, the known 2-loop and 3-loop QCD threshold corrections for $\alpha_s(M_Z)$ from [arXiv:hep-ph/9305305, arXiv:hep-ph/9708255, arXiv:hep-ph/9707474, arXiv:hep-ph/0004189] can be taken into account in addition by setting the threshold corrections flag appropriately. In the calculation of the Standard Model $\overline{\text{MS}}$ top Yukawa coupling, $y_t(M_Z)$, the known 2-loop [arXiv:hep-ph/9803493] and 3-loop [arXiv:hep-ph/9911434] QCD corrections can be taken into account. See Section FlexibleSUSY configuration block (FlexibleSUSY) for a description of the individual flags to enable/disable 2- and 3-loop threshold corrections in FlexibleSUSY.

__EWSB scale__

The Higgs and W boson pole masses, $M_h$ and $M_Z$ are calculated at the scale $M_{\text{EWSB}}$, which is an input parameter. We recommend to set $M_{\text{EWSB}} = M_t$.

Furthermore, the electroweak symmetry breaking condition is imposed at the scale $M_{\text{EWSB}}$ to fix the value of the bililear Higgs coupling $\mu^2(M_{\text{EWSB}})$.

Pole masses

The Higgs and W boson pole masses, $M_h$ and $M_Z$, are calculated at the full 1-loop level, including potential flavour mixing and momentum dependence. Depending on the given configuration flags, additional 2-, 3- and 4-loop corrections to the Higgs pole mass of $O(\alpha_t\alpha_s + \alpha_b\alpha_s)$ [arXiv:1407.4336] $O((\alpha_t + \alpha_b)^2)$ [arXiv:1205.6497] and $O(\alpha_\tau^2)$, as well as 3-loop corrections $O(\alpha_t^3 + \alpha_t^2\alpha_s + \alpha_t\alpha_s^2)$ [arXiv:1407.4336] and 4-loop corrections $O(\alpha_t\alpha_s^3)$ [arXiv:1508.00912] can be taken into account.

Note:
Note, that the 3-loop contributions $O(\alpha_t^3 + \alpha_t^2\alpha_s)$ are incomplete, because the corresponding 2-loop threshold corrections $O(\alpha_t^2 + \alpha_t\alpha_s)$ to the running top Yukawa coupling are not implemented yet.

Input parameters

HSSUSY takes the following physics parameters as input:

Parameter | Description | SLHA block/field | Mathematica symbol --------------------------------------------|----------------------------------------------------------------------|------------------|-------------------- $M_{\text{SUSY}}$ | SUSY scale | `EXTPAR[0]` | `MSUSY` $M_1(M_\text{SUSY})$ | Bino mass | `EXTPAR[1]` | `M1Input` $M_2(M_\text{SUSY})$ | Wino mass | `EXTPAR[2]` | `M2Input` $M_3(M_\text{SUSY})$ | Gluino mass | `EXTPAR[3]` | `M3Input` $\mu(M_\text{SUSY})$ | $\mu$-parameter | `EXTPAR[4]` | `MuInput` $m_A(M_\text{SUSY})$ | running CP-odd Higgs mass | `EXTPAR[5]` | `mAInput` $M_{\text{EWSB}}$ | scale at which the pole mass spectrum is calculated | `EXTPAR[5]` | `MEWSB` $A_t(M_\text{SUSY})$ | trililear stop coupling | `EXTPAR[7]` | `AtInput` $A_b(M_\text{SUSY})$ | trililear sbottom coupling | `EXTPAR[8]` | `AbInput` $A_\tau(M_\text{SUSY})$ | trililear stau coupling | `EXTPAR[9]` | `AtauInput` $\tan\beta(M_\text{SUSY})$ | $\tan\beta(M_\text{SUSY})=v_u(M_\text{SUSY})/v_d(M_\text{SUSY})$ | `EXTPAR[25]` | `TanBeta` $(m_{\tilde{q}}^2)_{ij}(M_\text{SUSY})$ | soft-breaking left-handed squark mass parameters | `MSQ2IN` | `msq2` $(m_{\tilde{u}}^2)_{ij}(M_\text{SUSY})$ | soft-breaking right-handed up-type squark mass parameters | `MSU2IN` | `msu2` $(m_{\tilde{d}}^2)_{ij}(M_\text{SUSY})$ | soft-breaking right-handed down-type squark mass parameters | `MSD2IN` | `msd2` $(m_{\tilde{l}}^2)_{ij}(M_\text{SUSY})$ | soft-breaking left-handed slepton mass parameters | `MSL2IN` | `msl2` $(m_{\tilde{e}}^2)_{ij}(M_\text{SUSY})$ | soft-breaking right-handed down-type slepton mass parameters | `MSE2IN` | `mse2`

The MSSM parameters are defined in the $\overline{\text{DR}}$ scheme at the scale $M_{\text{SUSY}}$.

In addition, HSSUSY defines further input parameters / flags to enable/disable higher order threshold corrections to the quartic Higgs coupling $\lambda(M_{\text{SUSY}})$ and to estimate the EFT and SUSY uncertainty:

Parameter | Description | Possible values | Recommended value | SLHA block/field | Mathematica symbol ---------------------------------------|-----------------------------------------------------------------------------------------------|-----------------|------------------------------|------------------|-------------------- $n$ | loop order for $\lambda^{(n)}(M_{\text{SUSY}})$ | 0, 1, 2 | 2 | `EXTPAR[100]` | `LambdaLoopOrder` $\Delta_{\alpha_t\alpha_s}$ | disable/enable 2-loop corrections to $\lambda(M_{\text{SUSY}})$ $O(\alpha_t\alpha_s)$ | 0, 1 | 1 | `EXTPAR[101]` | `TwoLoopAtAs` $\Delta_{\alpha_b\alpha_s}$ | disable/enable 2-loop corrections to $\lambda(M_{\text{SUSY}})$ $O(\alpha_b\alpha_s)$ | 0, 1 | 1 | `EXTPAR[102]` | `TwoLoopAbAs` $\Delta_{\alpha_t\alpha_b}$ | disable/enable 2-loop corrections to $\lambda(M_{\text{SUSY}})$ $O(\alpha_t\alpha_b)$ | 0, 1 | 1 | `EXTPAR[103]` | `TwoLoopAtAb` $\Delta_{\alpha_\tau\alpha_\tau}$ | disable/enable 2-loop corrections to $\lambda(M_{\text{SUSY}})$ $O(\alpha_\tau^2)$ | 0, 1 | 1 | `EXTPAR[104]` | `TwoLoopAtauAtau` $\Delta_{\alpha_t\alpha_t}$ | disable/enable 2-loop corrections to $\lambda(M_{\text{SUSY}})$ $O(\alpha_t^2)$ | 0, 1 | 1 | `EXTPAR[105]` | `TwoLoopAtAt` $\Delta_{\text{EFT}}$ | disable/enable corrections to $\lambda(M_{\text{SUSY}})$ $O(v^2/M_{\text{SUSY}}^2)$ | 0, 1 | 0 | `EXTPAR[200]` | `DeltaEFT` $\Delta_{y_t,g_3}$ | disable/enable 3-loop corrections from re-parametrization of $\lambda(M_{\text{SUSY}})$ in terms of $y_t^{\text{MSSM}}$, $g_3^{\text{MSSM}}$ | 0, 1 | 0 | `EXTPAR[201]` | `DeltaYt` $\Delta_{\text{OS}}$ | disable/enable conversion of stop masses to on-shell scheme | 0, 1 | 0 (= $\overline{\text{DR}}$ ) | `EXTPAR[202]` | `DeltaOS` $Q_\text{match}$ | scale at which $\lambda(Q_\text{match})$ is calculated | any real value | 0 (= $M_{\text{SUSY}})$ | `EXTPAR[203]` | `Qmatch`

Running HSSUSY

We recommend to run HSSUSY with the following configuration flags: In an SLHA input file we recommend to use:

~~~~~~~~~~~~~~~~~~~~~~~{.txt} Block FlexibleSUSY 0 1.0e-05 # precision goal 1 0 # max. iterations (0 = automatic) 2 0 # algorithm (0 = all, 1 = two_scale, 2 = semi_analytic) 3 1 # calculate SM pole masses 4 4 # pole mass loop order 5 4 # EWSB loop order 6 4 # beta-functions loop order 7 3 # threshold corrections loop order 8 1 # Higgs 2-loop corrections O(alpha_t alpha_s) 9 1 # Higgs 2-loop corrections O(alpha_b alpha_s) 10 1 # Higgs 2-loop corrections O((alpha_t + alpha_b)^2) 11 1 # Higgs 2-loop corrections O(alpha_tau^2) 12 0 # force output 13 1 # Top pole mass QCD corrections (0 = 1L, 1 = 2L, 2 = 3L) 14 1.0e-11 # beta-function zero threshold 15 0 # calculate observables (a_muon, ...) 16 0 # force positive majorana masses 17 0 # pole mass renormalization scale (0 = SUSY scale) 18 0 # pole mass renormalization scale in the EFT (0 = min(SUSY scale, Mt)) 19 0 # EFT matching scale (0 = SUSY scale) 20 2 # EFT loop order for upwards matching 21 1 # EFT loop order for downwards matching 22 0 # EFT index of SM-like Higgs in the BSM model 23 1 # calculate BSM pole masses 24 123111321 # individual threshold correction loop orders 25 0 # ren. scheme for Higgs 3L corrections (0 = DR, 1 = MDR) 26 1 # Higgs 3-loop corrections O(alpha_t alpha_s^2) 27 1 # Higgs 3-loop corrections O(alpha_b alpha_s^2) 28 1 # Higgs 3-loop corrections O(alpha_t^2 alpha_s) 29 1 # Higgs 3-loop corrections O(alpha_t^3) 30 1 # Higgs 4-loop corrections O(alpha_t alpha_s^3) ~~~~~~~~~~~~~~~~~~~~~~~

In the Mathematica interface we recommend to use:

~~~~~~~~~~~~~~~~~~~~~~~{.m} handle = FSHSSUSYOpenHandle[ fsSettings -> { precisionGoal -> 1.*^-5, (* FlexibleSUSY[0] *) maxIterations -> 0, (* FlexibleSUSY[1] *) solver -> 0, (* FlexibleSUSY[2] *) calculateStandardModelMasses -> 1, (* FlexibleSUSY[3] *) poleMassLoopOrder -> 4, (* FlexibleSUSY[4] *) ewsbLoopOrder -> 4, (* FlexibleSUSY[5] *) betaFunctionLoopOrder -> 4, (* FlexibleSUSY[6] *) thresholdCorrectionsLoopOrder -> 3,(* FlexibleSUSY[7] *) higgs2loopCorrectionAtAs -> 1, (* FlexibleSUSY[8] *) higgs2loopCorrectionAbAs -> 1, (* FlexibleSUSY[9] *) higgs2loopCorrectionAtAt -> 1, (* FlexibleSUSY[10] *) higgs2loopCorrectionAtauAtau -> 1, (* FlexibleSUSY[11] *) forceOutput -> 0, (* FlexibleSUSY[12] *) topPoleQCDCorrections -> 1, (* FlexibleSUSY[13] *) betaZeroThreshold -> 1.*^-11, (* FlexibleSUSY[14] *) forcePositiveMasses -> 0, (* FlexibleSUSY[16] *) poleMassScale -> 0, (* FlexibleSUSY[17] *) eftPoleMassScale -> 0, (* FlexibleSUSY[18] *) eftMatchingScale -> 0, (* FlexibleSUSY[19] *) eftMatchingLoopOrderUp -> 2, (* FlexibleSUSY[20] *) eftMatchingLoopOrderDown -> 1, (* FlexibleSUSY[21] *) eftHiggsIndex -> 0, (* FlexibleSUSY[22] *) calculateBSMMasses -> 1, (* FlexibleSUSY[23] *) thresholdCorrections -> 123111321, (* FlexibleSUSY[24] *) higgs3loopCorrectionRenScheme -> 0,(* FlexibleSUSY[25] *) higgs3loopCorrectionAtAsAs -> 1, (* FlexibleSUSY[26] *) higgs3loopCorrectionAbAsAs -> 1, (* FlexibleSUSY[27] *) higgs3loopCorrectionAtAtAs -> 1, (* FlexibleSUSY[28] *) higgs3loopCorrectionAtAtAt -> 1, (* FlexibleSUSY[29] *) higgs4loopCorrectionAtAsAsAs -> 1, (* FlexibleSUSY[30] *) parameterOutputScale -> 0 (* MODSEL[12] *) }, ... ]; ~~~~~~~~~~~~~~~~~~~~~~~

In the Section Mathematica interface an example Mathematica script can be found, which illustrates how to perform a parameter scan using the HSSUSY model.

Uncertainty estimate of the predicted Higgs pole mass

In the file model_files/HSSUSY/HSSUSY_uncertainty_estimate.m FlexibleSUSY provides the Mathematica function `CalcHSSUSYDMh[]`, which calculates the Higgs pole mass at the 3-loop level with HSSUSY and performs an uncertainty estimate of missing higher order corrections. Three main sources of the theory uncertainty are taken into account:

  • _SM uncertainty_: Missing higher order corrections in the calculation of the running Standard Model top Yukawa coupling and in the calculation of the Higgs pole mass. The uncertainty from this source is estimated by (i) switching on/off the 3-loop QCD contributions in the calculation of the running top Yukawa coupling $y_t^{\text{SM}}(M_Z)$ from the top pole mass and by (ii) varying the renormalization scale at which the Higgs pole mass is calculated within the interval $[M_{\text{EWSB}}/2, 2 M_{\text{EWSB}}]$.
  • _EFT uncertainty_: Missing terms of $O(v^2/M_{\text{SUSY}}^2)$. These missing terms are estimated by adding 1-loop terms of the form $v^2/M_{\text{SUSY}}^2$ to the quartic Higgs coupling $\lambda(M_\text{SUSY})$.
  • _SUSY uncertainty_: Missing higher order corrections in the calculation of the quartic Higgs coupling $\lambda(M_\text{SUSY})$. This uncertainty is estimated by (i) varying the matching scale within the interval $[M_{\text{SUSY}}/2, 2 M_{\text{SUSY}}]$ and by (ii) re-parametrization of $\lambda(M_\text{SUSY})$ in terms of $y_t^{\text{MSSM}}(M_\text{SUSY})$ and $g_3^{\text{MSSM}}(M_\text{SUSY})$.

The following code snippet illustrates the calculation of the Higgs pole mass calculated at the 3-loop level with HSSUSY as a function of the SUSY scale (red solid line), together with the estimated uncertainty (grey band).

Get["models/HSSUSY/HSSUSY_librarylink.m"];
Get["model_files/HSSUSY/HSSUSY_uncertainty_estimate.m"];

settings = {
    precisionGoal -> 1.*^-5,
    poleMassLoopOrder -> 3,
    ewsbLoopOrder -> 3,
    betaFunctionLoopOrder -> 3,
    thresholdCorrectionsLoopOrder -> 3,
    thresholdCorrections -> 123111321
};

smpars = {
    alphaEmMZ -> 1/127.916, (* SMINPUTS[1] *)
    GF -> 1.166378700*^-5,  (* SMINPUTS[2] *)
    alphaSMZ -> 0.1184,     (* SMINPUTS[3] *)
    MZ -> 91.1876,          (* SMINPUTS[4] *)
    mbmb -> 4.18,           (* SMINPUTS[5] *)
    Mt -> 173.34,           (* SMINPUTS[6] *)
    Mtau -> 1.77699,        (* SMINPUTS[7] *)
    Mv3 -> 0,               (* SMINPUTS[8] *)
    MW -> 80.385,           (* SMINPUTS[9] *)
    Me -> 0.000510998902,   (* SMINPUTS[11] *)
    Mv1 -> 0,               (* SMINPUTS[12] *)
    Mm -> 0.1056583715,     (* SMINPUTS[13] *)
    Mv2 -> 0,               (* SMINPUTS[14] *)
    md2GeV -> 0.00475,      (* SMINPUTS[21] *)
    mu2GeV -> 0.0024,       (* SMINPUTS[22] *)
    ms2GeV -> 0.104,        (* SMINPUTS[23] *)
    mcmc -> 1.27,           (* SMINPUTS[24] *)
    CKMTheta12 -> 0,
    CKMTheta13 -> 0,
    CKMTheta23 -> 0,
    CKMDelta -> 0,
    PMNSTheta12 -> 0,
    PMNSTheta13 -> 0,
    PMNSTheta23 -> 0,
    PMNSDelta -> 0,
    PMNSAlpha1 -> 0,
    PMNSAlpha2 -> 0,
    alphaEm0 -> 1/137.035999074,
    Mh -> 125.09
};

HSSUSYCalcMh[MS_, TB_, Xtt_] :=
    CalcHSSUSYDMh[
        fsSettings -> settings,
        fsSMParameters -> smpars,
        fsModelParameters -> {
            TanBeta -> TB,
            MEWSB -> 173.34,
            MSUSY -> MS,
            M1Input -> MS,
            M2Input -> MS,
            M3Input -> MS,
            MuInput -> MS,
            mAInput -> MS,
            AtInput -> (Xtt + 1/TB) MS,
            AbInput -> 0,
            AtauInput -> 0,
            msq2 -> MS^2 IdentityMatrix[3],
            msu2 -> MS^2 IdentityMatrix[3],
            msd2 -> MS^2 IdentityMatrix[3],
            msl2 -> MS^2 IdentityMatrix[3],
            mse2 -> MS^2 IdentityMatrix[3],
            LambdaLoopOrder -> 2,
            TwoLoopAtAs -> 1,
            TwoLoopAbAs -> 1,
            TwoLoopAtAb -> 1,
            TwoLoopAtauAtau -> 1,
            TwoLoopAtAt -> 1
        }
   ];

LinearRange[start_, stop_, steps_] :=
    Range[start, stop, (stop - start)/steps];

Xtt = Sqrt[6];
TB  = 5;

data = ParallelMap[
    { N[#], Sequence @@ HSSUSYCalcMh[#, TB, Xtt] }&,
    LinearRange[500, 10^4, 100]
];

MhMin[{MS_, Mh_, DMh_}]  := {MS, Mh - DMh};
MhMax[{MS_, Mh_, DMh_}]  := {MS, Mh + DMh};
MhBest[{MS_, Mh_, DMh_}] := {MS, Mh};

dataMhMin  = MhMin  /@ data;
dataMhMax  = MhMax  /@ data;
dataMhBest = MhBest /@ data;

plot2 = ListLinePlot[dataMhBest,
                     PlotStyle -> {Red, Thick}];

plot1 = ListLinePlot[{dataMhMax, dataMhMin},
                     PlotStyle -> LightGray,
                     Filling -> {1 -> {{2}, LightGray}},
                     PlotRange -> All];

plot = Show[{plot1, plot2},
            BaseStyle -> {FontSize -> 16, FontFamily -> "Helvetica"},
            PlotLabel -> Style["\*SubscriptBox[X, t] = 2.44949 \*SubscriptBox[M, S], tan\[Beta] = 5"],
            PlotRange -> Automatic,
            Axes -> False, Frame -> True,
            FrameLabel -> {Style["\*SubscriptBox[M, S] / GeV"],
                           Style["\*SubscriptBox[M, h] / GeV"]}];

Export["HSSUSY_Mh_MS.png", plot, ImageSize -> 600];

When this script is executed, the following figure is produced:

HSSUSY_Mh_MS.png