---
abstract: 'The program `N2HDECAY` calculates the branching ratios and decay widths of the Higgs bosons of the Next-to-Two-Higgs-Doublet Model (N2HDM).  The code incorporates the dominant higher-order effects by including QCD corrections and off-shell decay modes.  The N2HDM is an extension of the Standard Model by a Higgs doublet and a real Higgs singlet.  Its phenomenology can change dramatically depending on which global symmetries are broken by electroweak symmetry breaking.  It can feature a large visible Higgs sector in the broken phase, behave like a 2HDM with scalar singlet dark matter in the dark singlet phase, or, in the inert doublet phase, extend an inert doublet model by mixing a singlet with the SM Higgs boson.  `N2HDECAY` provides precise predictions for the decays of the Higgs bosons in all of these phases.' 
address:
- 'Institute for Theoretical Physics, Karlsruhe Institute of Technology, Wolfgang-Gaede-Str.  1, 76131 Karlsruhe, Germany'
- 'Theory Group, Deutsches Elektronen-Synchrotron DESY, Notkestr.  85, 22607 Hamburg, Germany'
author:
- Isabell Engeln
- Margarete Mühlleitner
- Jonas Wittbrodt
bibliography:
- 'Bibliography.bib'
title: '`N2HDECAY`: Higgs Boson Decays in the Different Phases of the N2HDM'
---

\

[*Program Title:*]{} `N2HDECAY`\
[*Licensing provisions:*]{} MIT\
[*Programming language:*]{} Fortran\
[*A code to calculate the branching ratios and total widths of the Higgs bosons in the N2HDM including state-of-the-art QCD and off-shell corrections.  The code is available at* ]{}<https://gitlab.com/jonaswittbrodt/N2HDECAY>. 

Introduction
============

Experimental studies conducted since the discovery of the 125GeV Higgs boson [@Aad:2012tfa; @Chatrchyan:2012xdj] have established it to behave Standard Model (SM)-like.  The resulting constraints on models beyond the SM (BSM) force minimal extensions of the Higgs sector into specific corners of the parameter space where new-physics effects in the 125GeV Higgs boson are small.  This makes minimal extensions of the SM scalar sector, such as the Two-Higgs-Doublet Model (2HDM) [@Lee:1973iz; @Branco:2011iw], less useful for experimental analysis since potentially interesting collider signatures are likely to be incompatible with other observations.  Therefore, non-minimal extensions of the Higgs sector become more appealing since they can account for interesting novel collider signatures while remaining compatible with existing measurements. 

The Next-to-Two-Higgs-Doublet Model (N2HDM) [@Chen:2013jvg; @Drozd:2014yla; @Muhlleitner:2016mzt] adds a real, gauge-singlet scalar field to the Higgs sector of the 2HDM.  If this singlet acquires a vacuum expectation value (VEV) it mixes with the two CP-even Higgs bosons of the 2HDM.  The resulting enlargement of the parameter space with interesting collider signatures has been studied in [@Muhlleitner:2016mzt] using the first published version of the program `N2HDECAY`.  The N2HDM is similar to the 2HDM plus singlet (2HDM+S) framework [@Curtin:2013fra] where a complex singlet scalar is added to an aligned 2HDM.  However, in the N2HDM we do not assume the 125GeV Higgs boson to be exactly aligned. 

An additional feature of non-minimal extensions of the Higgs sector is that they may both modify the visible phenomenology of the model and add dark matter (DM) candidates.  Both dark singlet DM (for a recent study see [@Athron:2017kgt]) and inert doublet DM [@Deshpande:1977rw; @Ma:2006km; @LopezHonorez:2006gr; @Belyaev:2016lok; @Ilnicka:2015jba] can be realized within the N2HDM as different phases of electroweak symmetry breaking (EWSB).  In the dark singlet phase (DSP) the singlet field does not acquire a VEV and becomes a DM candidate.  In the inert doublet phase (IDP) one of the doublets has no VEV and forms a dark sector.  A phenomenological study of the dark phases of the N2HDM using `N2HDECAY` has been performed recently  [@Engeln:2018ywp]. 

The code `N2HDECAY` calculates branching ratios and total widths of the Higgs bosons in these three phases of the N2HDM.  It includes state-of-the-art higher-order QCD corrections and off-shell effects taken over from its parent code `HDECAY` [@Djouadi:1997yw; @Butterworth:2010ym; @Djouadi:2018xqq].  `N2HDECAY` therefore provides precise predictions for the phenomenology of the N2HDM Higgs bosons.  In this manual we first define the N2HDM with its three phases in \[sec:model\].  Afterwards in \[sec:corrections\], we discuss the higher-order corrections included in the code.  In \[sec:usage\] we give user instructions for the code, describe the input and output formats as well as additional tools distributed with the code.  Finally, \[app:output\] lists all decay modes calculated in `N2HDECAY` and their position in the output files while \[app:DSP\] and \[app:IDP\] contain the expressions for the parameter transformations and triple-Higgs couplings in the DSP and IDP, respectively. 

The Next-to Two Higgs Doublet Model {#sec:model}
===================================

The N2HDM [@Chen:2013jvg; @Drozd:2014yla; @Muhlleitner:2016mzt] extends the SM by a second Higgs doublet and a real Higgs singlet.  The resulting scalar potential for the doublets $\Phi_{1,2}$ and the singlet $\Phi_S$ reads $$\begin{aligned}
  V &= m_{11}^2 |\Phi_1|^2 + m_{22}^2 |\Phi_2|^2 - m_{12}^2 (\Phi_1^\dagger
  \Phi_2 + \text{h.c.})  + \frac{\lambda_1}{2} (\Phi_1^\dagger \Phi_1)^2 +
  \frac{\lambda_2}{2} (\Phi_2^\dagger \Phi_2)^2  \\
  &\hphantom{=} + \lambda_3
  (\Phi_1^\dagger \Phi_1) (\Phi_2^\dagger \Phi_2) + \lambda_4
  (\Phi_1^\dagger \Phi_2) (\Phi_2^\dagger \Phi_1) + \frac{\lambda_5}{2}
  [(\Phi_1^\dagger \Phi_2)^2 + \text{h.c.}]   \\
  &\hphantom{=} + \frac{m_S^2}{2}  \Phi_S^2 + \frac{\lambda_6}{8} \Phi_S^4 +
  \frac{\lambda_7}{2} (\Phi_1^\dagger \Phi_1) \Phi_S^2 +
  \frac{\lambda_8}{2} (\Phi_2^\dagger \Phi_2) \Phi_S^2 \,,
\end{aligned}
  \label{eq:n2hdmpot}$$ where all quartic couplings $\lambda_{1,\ldots,8}$ and mass parameters $m_{11,22,12,S}^2$ are real.  The first two lines of \[eq:n2hdmpot\], containing only the doublet fields, correspond to the standard 2HDM [@Branco:2011iw].  The last line contains the contributions of the singlet field. 

We impose two $\mathbb{Z}_2$ parities on the scalar potential \[eq:n2hdmpot\].  The symmetry called $\mathbb{Z}_2$, $$\begin{aligned}
  \Phi_1\to\Phi_1\,,\ \Phi_2\to-\Phi_2\,,\ \Phi_S\to\Phi_S\,, \label{eq:Z2}\end{aligned}$$ is the trivial generalization of the usual 2HDM $\mathbb{Z}_2$ symmetry to the N2HDM.  It can be extended to the Yukawa sector to prevent flavor changing neutral currents.  It is softly broken by $m_{12}^2\neq0$ and spontaneously broken if $\Phi_2$ develops a non-zero VEV.  The second symmetry, $\mathbb{Z}'_2$, $$\begin{aligned}
  \Phi_1\to\Phi_1\,,\ \Phi_2\to\Phi_2\,,\ \Phi_S\to-\Phi_S\,, \label{eq:Z2p}\end{aligned}$$ is spontaneously broken if the VEV of $\Phi_S$ is non-zero.  If any of these two parities remain unbroken after EWSB they lead to a dark sector containing a stable dark matter candidate. 

After EWSB we expand the Higgs fields around their VEVs $v_1$, $v_2$ and $v_s$ as $$\Phi_{i=\lbrace 1,2\rbrace} = \begin{pmatrix} \phi_i^+ \\ \frac{1}{\sqrt{2}} (v_i +
    \rho_i + i \eta_i) \end{pmatrix}\,,\quad
\Phi_S = v_S + \rho_S \,.  \label{eq:n2hdmfields}$$ The minimum conditions for the electroweak vacuum lead to the relations $$\begin{aligned}
v_2 m_{12}^2 - v_1 m_{11}^2 &= \frac{1}{2} (v_1^3 \lambda_1 +
v_2^2 v_1 \lambda_{345} + v_S^2 v_1 \lambda_7)\,, \label{eq:n2hdmmin1} \\
v_1 m_{12}^2 - v_2 m_{22}^2 &= \frac{1}{2} (v_1^2 v_2 \lambda_{345} +
v_2^3 \lambda_2 + v_S^2 v_2 \lambda_8)\,, \label{eq:n2hdmmin2} \\
- v_S m_S^2 &= \frac{1}{2} v_S (v_1^2 \lambda_7 + v_2^2 \lambda_8 + v_S^2
\lambda_6)\,,  \label{eq:n2hdmmin3}\end{aligned}$$ with $$\lambda_{345} \equiv \lambda_3 + \lambda_4 + \lambda_5 \,. 
\label{eq:l345}$$ We have not canceled any powers of VEVs in \[eq:n2hdmmin1,eq:n2hdmmin2,eq:n2hdmmin3\] to allow for the possibility of vanishing $v_{1,2,S}$. These conditions can be used to express the mass parameters $m_{11}^2$, $m_{22}^2$ and $m_{S}^2$ through the VEVs and allow for four different phases depending on which VEVs are non-zero. 

The Broken Phase
----------------

In the broken phase all CP- and charge-conserving VEVs ($v_1$, $v_2$ and $v_S$) are non-zero.  This breaks both $\mathbb{Z}_2$ and $\mathbb{Z}'_2$. Like in the 2HDM, the rotation into the mass basis with the matrix $$\begin{aligned}
  R_\beta = \begin{pmatrix}
  \cos\beta &  \sin\beta\\
  -\sin\beta & \cos\beta
\end{pmatrix}\end{aligned}$$ yields one physical CP-odd Higgs boson $A$ and the neutral Goldstone boson $G^0$ from $\eta_{1,2}$ as well as a pair of charged Higgs bosons $H^\pm$ and charged Goldstone bosons $G^\pm$ from $\phi^\pm_{1,2}$. However, $v_S\neq0$ leads to a mixing between the $\rho_{1,2}$ and $\rho_S$ resulting in three physical CP-even Higgs bosons $$\begin{pmatrix} H_1 \\ H_2 \\ H_3 \end{pmatrix} = R \begin{pmatrix}\rho_1 \\ \rho_2 \\ \rho_S\end{pmatrix}\,.$$ The $3\times3$ orthogonal mixing matrix $R$ is parametrized by three angles $\alpha_{1,2,3}$ as $$R = \begin{pmatrix}
  c_{\alpha_1} c_{\alpha_2} & s_{\alpha_1} c_{\alpha_2} & s_{\alpha_2}\\
  -(c_{\alpha_1} s_{\alpha_2} s_{\alpha_3} + s_{\alpha_1} c_{\alpha_3})
  & c_{\alpha_1} c_{\alpha_3} - s_{\alpha_1} s_{\alpha_2} s_{\alpha_3}
  & c_{\alpha_2} s_{\alpha_3} \\
  - c_{\alpha_1} s_{\alpha_2} c_{\alpha_3} + s_{\alpha_1} s_{\alpha_3} &
  -(c_{\alpha_1} s_{\alpha_3} + s_{\alpha_1} s_{\alpha_2} c_{\alpha_3})
  & c_{\alpha_2}  c_{\alpha_3}
  \end{pmatrix}\,,
  \label{eq:mixingmatrix}$$ where we introduced the short-hand notations $s_x\equiv\sin x$ and $c_x\equiv\cos x$. We assume that the three mass eigenstates $H_{1,2,3}$ are ordered as $m_{H_1}<m_{H_2}<m_{H_3}$.

                             $u$-type                    $d$-type                     leptons
  ----------------- --------------------------- --------------------------- ---------------------------
             type I  $\frac{R_{i2}}{s_\beta}$    $\frac{R_{i2}}{s_\beta}$    $\frac{R_{i2}}{s_\beta}$
            type II  $\frac{R_{i2}}{s_\beta} $   $\frac{R_{i1}}{c_\beta} $   $\frac{R_{i1}}{c_\beta} $
    lepton-specific  $\frac{R_{i2}}{s_\beta}$    $\frac{R_{i2}}{s_\beta}$    $\frac{R_{i1}}{c_\beta}$
            flipped  $\frac{R_{i2}}{s_\beta}$    $\frac{R_{i1}}{c_\beta}$    $\frac{R_{i2}}{s_\beta}$

  : Coupling coefficients $c(H_i ff)$ of the Yukawa couplings of the N2HDM Higgs bosons $H_i$. \[tab:yukcoup\]

The three CP-even Higgs bosons couple to the SM gauge bosons $V\in\{W^\pm,Z\}$ with an effective coupling normalized to the SM given by $$c(H_i VV)=R_{i 1}\cos\beta + R_{i2}\sin\beta\,.\label{eq:gaugecoupbroken}$$ To prevent tree-level flavor-changing neutral currents we extend the symmetry $\mathbb{Z}_2$, \[eq:Z2\], to the Yukawa sector.  As in the 2HDM [@Branco:2011iw], this results in the well-known four distinct Yukawa types.  The effective couplings for the different types normalized to their SM values are given in \[tab:yukcoup\].  The couplings of $A$ and $H^\pm$ to SM particles are the same as in the 2HDM and can be found in the corresponding literature [@Branco:2011iw].  All triple-Higgs couplings needed for the computation of the decay widths are given in the appendix of [@Muhlleitner:2016mzt]. 

In the broken phase we use the input parameters $$\alpha_{1,2,3}\,, \quad \tan\beta \,, \quad v \, ,
  \quad v_S \, , \quad m_{H_{1,2,3}} \,, \quad m_A \,, \quad m_{H^\pm}
  \,, \quad m_{12}^2 \,.  \label{eq:brokeninputpars}$$ The parameters $\tan\beta$ and $v$ are given by $$\begin{aligned}
  \tan\beta=\frac{v_2}{v_1}\,,\quad v=v_1^2+v_2^2\approx 246\, \text{GeV}\,,\label{eq:tbetadef}\end{aligned}$$ where the value of $v$ is the SM VEV.  The relations between this parameter set and the parameters of the scalar potential \[eq:n2hdmpot\] can be found in the appendix of [@Muhlleitner:2016mzt]. 

The Dark Singlet Phase
----------------------

The dark singlet phase (DSP) is characterized by $v_S=0$ such that $\mathbb{Z}'_2$, \[eq:Z2p\], remains exact after EWSB.  Any mixing between $\rho_{1,2}$ and $\rho_S$ is therefore forbidden and the CP-even mass eigenstates are $$\begin{pmatrix} H_1 \\ H_2 \\ H_D\\ \end{pmatrix} = R \begin{pmatrix}\rho_1 \\ \rho_2 \\ \rho_S\\\end{pmatrix}$$ with $$R=\begin{pmatrix}
    -\sin\alpha &\cos\alpha &0\\
    \cos\alpha &\sin\alpha &0\\
    0&0&1\\
  \end{pmatrix}\,.\label{eq:mixmatdsp}$$ The mixing matrix $R$ is chosen such that the mixing angle $\alpha$ follows the conventions for the 2HDM [@Branco:2011iw] if we identify $h\equiv H_1$ and $H\equiv -H_2$. For consistency with the other phases we will call the Higgs bosons $H_1$ and $H_2$ throughout this paper and assume them to be ordered as $m_{H_1}<m_{H_2}$. The singlet scalar $H_D$ is a DM candidate stabilized by $\mathbb{Z}'_2$. The charged and CP-odd sector are identical to the broken phase.  The DSP of the N2HDM is equivalent to adding scalar singlet DM to a 2HDM. 

The couplings of the Higgs bosons to SM particles in the DSP can be obtained from \[eq:gaugecoupbroken\] and \[tab:yukcoup\] by using the definition \[eq:mixmatdsp\].  The triple-Higgs couplings differ from the broken phase since all vertices involving an odd number of $H_D$ are forbidden by $\mathbb{Z}'_2$. The expressions for the relevant couplings can be found in \[app:tripleHDSP\]. 

In the DSP we cannot express all of the quartic couplings of the scalar potential \[eq:n2hdmpot\] through masses and mixing angles.  The two quartic Higgs-portal couplings $\lambda_7$ and $\lambda_8$ as well as the singlet self-coupling $\lambda_6$ remain free parameters.  However, $\lambda_6$ does not enter any of the three-particle vertices relevant for Higgs decays and can be ignored for the purposes of this code.  We therefore use the input parameters $$\alpha\,, \quad \tan\beta \,, \quad v \, ,
  \quad m_{H_{1,2}}\,, \quad m_{H_D}\,, \quad m_A \,, \quad m_{H^\pm}
  \,, \quad m_{12}^2\,,\quad \lambda_{7,8}  \,, \label{eq:dspinputpars}$$ with $\tan\beta$ and $v$ defined through \[eq:tbetadef\].  The transformations from these parameters to the parameters of \[eq:n2hdmpot\] are given in \[app:transDSP\]. 

The Inert Doublet Phase
-----------------------

The inert doublet phase (IDP) can only be realized if $m_{12}^2=0$ such that $\mathbb{Z}_2$, \[eq:Z2\], is not explicitly broken.  If additionally $v_2=0$ the doublet $\Phi_2$ is inert [@Deshpande:1977rw; @Ma:2006km; @LopezHonorez:2006gr; @Belyaev:2016lok; @Ilnicka:2015jba].  The component fields of $\Phi_2$ are all mass eigenstates and form a dark sector consisting of a CP-even scalar $H_D$, a CP-odd scalar $A_D$ and a charged scalar $H^\pm_D$. The only mixing in the IDP is between the singlet field $\rho_S$ and the doublet field $\rho_1$. This results in two visible CP-even neutral Higgs bosons $H_{1,2}$ through $$\begin{pmatrix} H_1 \\ H_2 \\ H_D\\ \end{pmatrix} = R \begin{pmatrix}\rho_1 \\ \rho_2 \\ \rho_S\\\end{pmatrix}$$ with $$R=\begin{pmatrix}
    \cos\alpha &  0 &\sin\alpha\\
    -\sin\alpha & 0 &\cos\alpha\\
    0&1&0\\
  \end{pmatrix}\,,\label{eq:mixmatIDP}$$ where we again take $m_{H_1}<m_{H_2}$.

The IDP is equivalent to extending an inert doublet model by a real singlet which mixes with the SM Higgs boson.  The couplings of the $H_{i=1,2}$ to any pair of SM particles $x$ are therefore simply rescaled by $$c(H_i xx) = R_{i1}\,.$$ Since this applies to the Yukawa couplings as well there are no different Yukawa types in the IDP.  The triple-Higgs couplings of the IDP can be found in \[app:tripleHIDP\]. 

In the IDP \[eq:n2hdmmin2\] is trivially satisfied so that we cannot express $m_{22}^2$ through the other parameters.  Additionally, the inert doublet self-coupling $\lambda_2$ and the portal coupling $\lambda_8$ cannot be expressed through masses and mixing angles.  However, $\lambda_2$ does not appear in any of the vertices relevant for Higgs boson decays and is not needed by the code.  The set of required input parameters for the IDP are therefore $$\alpha \,,\quad v \,,\quad v_S \,,\quad m_{H_{1,2}} \,,\quad m_{H_D} \,,\quad m_{A_D} \,,\quad m_{H^\pm_D} \,,\quad m_{22}^2 \,,\quad \lambda_8\,,\label{eq:idpinputpars}$$ where $v_1=v$ is the SM VEV.  The relations between the parameters of \[eq:idpinputpars\] and \[eq:n2hdmpot\] are given in \[app:transIDP\]. 

The SM+DM Phase
---------------

If $m_{12}^2=0$ and both $v_2$ and $v_S$ vanish the N2HDM becomes the SM extended by both an inert doublet and a scalar DM singlet.  This phase is currently not implemented in the code. 

Higgs Decays in N2HDECAY {#sec:corrections}
========================

The code `N2HDECAY` is an extension of the 2HDM implementation in `HDECAY` [@Djouadi:1997yw; @Butterworth:2010ym; @Harlander:2013qxa; @Djouadi:2018xqq] to the N2HDM.  The calculation of the decay widths includes the most relevant QCD corrections since they factorize and can be taken over from the SM or the MSSM.  The SM electroweak corrections, on the other hand, cannot be easily adapted to models with extended scalar sectors and are not taken into account.  They have been recently provided in [@Krause:2017mal].  The currently implemented higher-order corrections are

-   massive NLO corrections near threshold [@Braaten:1980yq; @Sakai:1980fa; @Inami:1980qp; @Drees:1989du; @Drees:1990dq] and massless $\mathcal{O}(\alpha_s^4)$ corrections far above threshold [@Gorishnii:1983cu; @Gorishnii:1990zu; @Gorishnii:1991zr; @Kataev:1993be; @Surguladze:1994gc; @Larin:1995sq; @Melnikov:1995yp; @Chetyrkin:1995pd; @Chetyrkin:1996sr; @Baikov:2005rw] to decays into quark pairs. 

-   $\text{N}^3\text{LO}$ corrections in the heavy-top-limit to decays into $gg$ [@Inami:1982xt; @Djouadi:1991tka; @Spira:1993bb; @Spira:1995rr; @Chetyrkin:1997iv; @Chetyrkin:1997un; @Kramer:1996iq; @Chetyrkin:1998mw; @Schroder:2005hy; @Chetyrkin:2005ia; @Baikov:2006ch]. 

-   QCD corrections to the quark loops in the loop-mediated ($W$, heavy-fermion and charged-Higgs-boson loops) decay into $\gamma\gamma$ including the full quark-mass dependence [@Zheng:1990qa; @Djouadi:1990aj; @Djouadi:1993ji; @Melnikov:1993tj; @Dawson:1992cy; @Inoue:1994jq; @Spira:1995rr; @Fleischer:2004vb; @Harlander:2005rq]. 

The QCD corrections to quark loops in the loop-induced $Z\gamma$ decay mode are numerically small [@Spira:1991tj; @Bonciani:2015eua; @Gehrmann:2015dua] and are not taken into account. 

The code additionally includes decays into off-shell gauge bosons and top-quarks where appropriate.  This includes

-   decays into a pair of off-shell massive gauge bosons [@Cahn:1988ru],

-   off-shell gauge bosons in all decays into a Higgs boson and a gauge boson [@Moretti:1994ds; @Djouadi:1995gv],

-   off-shell $t$-quarks in decays into $t$-quark pairs and the $t\bar{b}$, $t\bar{s}$ and $t\bar{d}$ decay modes of the charged Higgs boson $H^+$ [@Moretti:1994ds; @Djouadi:1995gv]. 

`N2HDECAY` does not include off-shell Higgs bosons in decays into pairs of Higgs bosons.  These off-shell corrections are only available in `HDECAY` assuming that the off-shell Higgs boson decays predominantly into $b\bar{b}$. This is justified for the 125GeV Higgs boson but may be a bad approximation for other Higgs bosons of the N2HDM.  Off-shell corrections are therefore disabled for all decays into Higgs boson pairs.  Finally, the code also calculates the branching ratios and total width of the top-quark to account for possibly allowed decays into $H^+b$. A complete list of all decay channels calculated by `N2HDECAY` including off-shell decays is given in \[tab:BP\_output,tab:DSP\_output,tab:IDP\_output\]. 

User Instructions {#sec:usage}
=================

The `N2HDECAY` source code is available at

<https://gitlab.com/jonaswittbrodt/N2HDECAY>. 

The compilation of the code requires only a Fortran compiler.  The code can be compiled either with the supplied `Makefile`, in which the correct compiler has to be selected manually, or using `cmake`.  The code compiles into a single executable called `n2hdecay`.  It reads input from a file called `n2hdecay.in` in the current working directory and writes its output files containing the branching ratios and total widths into the current directory.  Input file and output folder may be changed by providing arguments to the optional `-i` and `-o` command line flags, respectively. 

The input file is a plain text file containing the parameters of the model.  The first parameter to be specified is the `PHASE` which is 1 for the broken phase, 2 for the DSP and 3 for the IDP.  The parameter `YUKTYPE` specifies the Yukawa type.  The possible values are 1 (Type I), 2 (Type II), 3 (lepton-specific) and 4 (flipped).  In the next block the parameters $\tan\beta$, $m_{12}^2$, $m_A$, $m_{H^\pm}$ and $v_S$ can be specified.  This block is used for every phase.  If the value for one of the parameters is fixed by the phase the value given here is ignored.  After these initial blocks of common parameters there is a dedicated block containing the remaining parameters for each phase.  The final block of the input file contains the SM parameters defined as in `HDECAY` (see [@Djouadi:2018xqq]).  An overview of the required N2HDM parameters for each phase is shown in \[tab:inputpars\].  The SM VEV $v$, which is an input parameter in every phase (see \[eq:brokeninputpars,eq:dspinputpars,eq:idpinputpars\]), is calculated internally from the value of $G_F$ given in the SM parameter block. 

[0.32]{}

[l c]{} Parameter & Value\
`PHASE` & 1\
`YUKTYPE` & 1,2,3,4\
`tbeta` & $\tan\beta$\
`m12sq` & $m_{12}^2$\
`m_A` & $m_A$\
`m_Hc` & $m_{H^\pm}$\
`vs` & $v_S$\
\
`m_H1` & $m_{H_1}$\
`m_H2` & $m_{H_2}$\
`m_H3` & $m_{H_3}$\
`a1` & $\alpha_1$\
`a2` & $\alpha_2$\
`a3` & $\alpha_3$\

[0.32]{}

[l c]{} Parameter & Value\
`PHASE` & 2\
`YUKTYPE` & 1,2,3,4\
`tbeta` & $\tan\beta$\
`m12sq` & $m_{12}^2$\
`m_A` & $m_A$\
`m_Hc` & $m_{H^\pm}$\
`vs` & –\
\
`m_H1` & $m_{H_1}$\
`m_H2` & $m_{H_2}$\
`m_HD` & $m_{H_D}$\
`alpha` & $\alpha$\
`L7` & $\lambda_7$\
`L8` & $\lambda_8$\

[0.32]{}

[l c]{} Parameter & Value\
`PHASE` & 3\
`YUKTYPE` & –\
`tbeta` & –\
`m12sq` & –\
`m_A` & $m_{A_D}$\
`m_Hc` & $m_{H^\pm_D}$\
`vs` & $v_S$\
\
`m_H1` & $m_{H_1}$\
`m_H2` & $m_{H_2}$\
`m_HD` & $m_{H_D}$\
`alpha` & $\alpha$\
`L8` & $\lambda_8$\
`m22sq` & $m_{22}^2$\

Running the code produces a set of output files containing the branching ratios and decay widths of all Higgs bosons and the top-quark.  The files are named `br. Hx_N2HDM_y` where `Hx` denotes the decaying Higgs bosons.  For neutral `Hx`=`H1,H2,H3,A` the files with `y`=`a` contain the decays into SM fermions, `y`=`b` the decays into SM gauge bosons and `y`=`c,d` the decays into double-Higgs and gauge-Higgs final states as well as the total width.  The decays of the charged Higgs boson `Hx=H+` into fermions are given in the `y`=`a,b` files and `y`=`c` contains the decays into gauge-Higgs final states and the total width.  The gauge-Higgs decay modes of the dark-sector particles in the IDP `Hx=HD,AD,HD+` along with their total widths can be found in the files `br. Hx_N2HDM`.  The mass of the decaying particle is included in every file as the first entry.  The file `br.top` contains the mass of the charged Higgs boson, the top-quark branching ratios and its total width.  Each of the files includes a header specifying the content.  For a complete overview of the output files and their contents in the different phases see \[tab:BP\_output,tab:DSP\_output,tab:IDP\_output\]. 

For the user’s convenience we also provide a Python 3 script

`tools/IO_N2HDECAY.py`

with the utility functions:

-   `IO_N2HDECAY.write(filepath, phase, params)`: Generates an input file for `N2HDECAY` at `filepath` for the given `phase` and the parameters taken from the Python dictionary `params`. 

-   `IO_N2HDECAY.read(brfolder)`: Reads the results of an `N2HDECAY` run from the `brfolder` and returns them as a Python dictionary. 

In the subfolder `tests` we provide three example input files (one for each phase) with the corresponding results of the code.  The folder also contains a Python 3 script `testN2HDECAY.py` which can be used to compare the program output against these example results.  Invoking `make test` when compiling with `cmake` automatically calls this script for each phase. 

Finally, the source code also includes an up-to-date version of this manual and a `CHANGELOG` file detailing all changes in the code since the initial release. 

Decay Modes and Output Format {#app:output}
=============================

The \[tab:BP\_output,tab:DSP\_output,tab:IDP\_output\] in this section list all decay modes calculated by `N2HDECAY` and the content of the output files in all three phases. 

  file                                                                                                                                      
  ----------------- ------------------- ---------------------- ----------------------- ------------------------------ --------------------- ------------------------
  `br. H1_N2HDM_a`       $b\bar{b}$         $\tau\bar{\tau}$        $\mu\bar{\mu}$                $s\bar{s}$                $c\bar{c}$        $t^{(*)}\bar{t}^{(*)}$
  `br. H1_N2HDM_b`          $g g$            $\gamma\gamma$            $Z\gamma$         ${W^\pm}^{(*)}{W^\mp}^{(*)}$    $Z^{(*)}Z^{(*)}$    
  `br. H1_N2HDM_c`          $A A$             $Z^{(*)} A$        ${W^\pm}^{(*)} H^\mp$          $H^\pm H^\mp$           $\Gamma_\text{tot}$  
                                                                                                                                            
  `br. H2_N2HDM_a`       $b\bar{b}$         $\tau\bar{\tau}$        $\mu\bar{\mu}$                $s\bar{s}$                $c\bar{c}$        $t^{(*)}\bar{t}^{(*)}$
  `br. H2_N2HDM_b`          $g g$            $\gamma\gamma$            $Z\gamma$         ${W^\pm}^{(*)}{W^\mp}^{(*)}$    $Z^{(*)}Z^{(*)}$    
  `br. H2_N2HDM_c`        $H_1 H_1$              $A A$                $Z^{(*)} A$           ${W^\pm}^{(*)} H^\mp$          $H^\pm H^\mp$       $\Gamma_\text{tot}$
                                                                                                                                            
  `br. H3_N2HDM_a`       $b\bar{b}$         $\tau\bar{\tau}$        $\mu\bar{\mu}$                $s\bar{s}$                $c\bar{c}$        $t^{(*)}\bar{t}^{(*)}$
  `br. H3_N2HDM_b`          $g g$            $\gamma\gamma$            $Z\gamma$         ${W^\pm}^{(*)}{W^\mp}^{(*)}$    $Z^{(*)}Z^{(*)}$    
  `br. H3_N2HDM_c`        $H_1 H_1$            $H_1 H_2$               $H_2 H_2$                    $A A$                   $Z^{(*)} A$       ${W^\pm}^{(*)} H^\mp$
  `br. H3_N2HDM_d`      $H^\pm H^\mp$     $\Gamma_\text{tot}$                                                                                
                                                                                                                                            
  `br. A_N2HDM_a`        $b\bar{b}$         $\tau\bar{\tau}$        $\mu\bar{\mu}$                $s\bar{s}$                $c\bar{c}$        $t^{(*)}\bar{t}^{(*)}$
  `br. A_N2HDM_b`           $g g$            $\gamma\gamma$            $Z\gamma$                                                             
  `br. A_N2HDM_c`       $Z^{(*)}H_1$          $Z^{(*)}H_2$           $Z^{(*)}H_3$           ${W^\pm}^{(*)} H^\mp$       $\Gamma_\text{tot}$  
                                                                                                                                            
  `br. H+_N2HDM_a`       $c\bar{b}$       $\tau\bar{\nu_\tau}$    $\mu\bar{\nu_\mu}$              $u\bar{s}$                $c\bar{s}$           $t^{(*)}\bar{b}$
  `br. H+_N2HDM_b`       $c\bar{d}$            $u\bar{b}$          $t^{(*)}\bar{s}$            $t^{(*)}\bar{d}$                              
  `br. H+_N2HDM_c`    ${W^+}^{(*)} H_1$    ${W^+}^{(*)} H_2$       ${W^+}^{(*)} H_3$           ${W^+}^{(*)} A$          $\Gamma_\text{tot}$  
                                                                                                                                            
  `br.top`                $W^+ b$              $H^+ b$           $\Gamma_\text{tot}$                                                        

  : The branching ratios calculated in the broken phase and their location in the output files.  Particles marked with $^{(*)}$ are allowed to be off-shell.[]{data-label="tab:BP_output"}

  file                                                                                                                                          
  ----------------- --------------------- ---------------------- ----------------------- ------------------------------ ----------------------- ------------------------
  `br. H1_N2HDM_a`        $b\bar{b}$          $\tau\bar{\tau}$        $\mu\bar{\mu}$                $s\bar{s}$                 $c\bar{c}$         $t^{(*)}\bar{t}^{(*)}$
  `br. H1_N2HDM_b`           $g g$             $\gamma\gamma$            $Z\gamma$         ${W^\pm}^{(*)}{W^\mp}^{(*)}$     $Z^{(*)}Z^{(*)}$     
  `br. H1_N2HDM_c`         $H_D H_D$               $A A$                $Z^{(*)} A$           ${W^\pm}^{(*)} H^\mp$           $H^\pm H^\mp$        $\Gamma_\text{tot}$
                                                                                                                                                
  `br. H2_N2HDM_a`        $b\bar{b}$          $\tau\bar{\tau}$        $\mu\bar{\mu}$                $s\bar{s}$                 $c\bar{c}$         $t^{(*)}\bar{t}^{(*)}$
  `br. H2_N2HDM_b`           $g g$             $\gamma\gamma$            $Z\gamma$         ${W^\pm}^{(*)}{W^\mp}^{(*)}$     $Z^{(*)}Z^{(*)}$     
  `br. H2_N2HDM_c`         $H_D H_D$             $H_1 H_1$                 $A A$                   $Z^{(*)} A$            ${W^\pm}^{(*)} H^\mp$       $H^\pm H^\mp$
  `br. H2_N2HDM_d`    $\Gamma_\text{tot}$                                                                                                        
                                                                                                                                                
  `br. A_N2HDM_a`         $b\bar{b}$          $\tau\bar{\tau}$        $\mu\bar{\mu}$                $s\bar{s}$                 $c\bar{c}$         $t^{(*)}\bar{t}^{(*)}$
  `br. A_N2HDM_b`            $g g$             $\gamma\gamma$            $Z\gamma$                                                               
  `br. A_N2HDM_c`        $Z^{(*)}H_1$           $Z^{(*)}H_2$       ${W^\pm}^{(*)} H^\mp$       $\Gamma_\text{tot}$                               
                                                                                                                                                
  `br. H+_N2HDM_a`        $c\bar{b}$        $\tau\bar{\nu_\tau}$    $\mu\bar{\nu_\mu}$              $u\bar{s}$                 $c\bar{s}$            $t^{(*)}\bar{b}$
  `br. H+_N2HDM_b`        $c\bar{d}$             $u\bar{b}$          $t^{(*)}\bar{s}$            $t^{(*)}\bar{d}$                                
  `br. H+_N2HDM_c`     ${W^+}^{(*)} H_1$     ${W^+}^{(*)} H_2$        ${W^+}^{(*)} A$          $\Gamma_\text{tot}$                               
                                                                                                                                                
  `br.top`                 $W^+ b$               $H^+ b$           $\Gamma_\text{tot}$                                                          

  : The branching ratios calculated in the DSP and their location in the output files.  Particles marked with $^{(*)}$ are allowed to be off-shell.[]{data-label="tab:DSP_output"}

  file                                                                                                                                             
  ----------------- --------------------- ------------------------- --------------------- ------------------------------ ------------------------- ------------------------
  `br. H1_N2HDM_a`        $b\bar{b}$           $\tau\bar{\tau}$         $\mu\bar{\mu}$               $s\bar{s}$                  $c\bar{c}$          $t^{(*)}\bar{t}^{(*)}$
  `br. H1_N2HDM_b`           $g g$              $\gamma\gamma$             $Z\gamma$        ${W^\pm}^{(*)}{W^\mp}^{(*)}$      $Z^{(*)}Z^{(*)}$      
  `br. H1_N2HDM_c`         $H_D H_D$               $A_D A_D$             $Z^{(*)} A_D$        ${W^\pm}^{(*)} H_D^\mp$         $H_D^\pm H_D^\mp$       $\Gamma_\text{tot}$
                                                                                                                                                   
  `br. H2_N2HDM_a`        $b\bar{b}$           $\tau\bar{\tau}$         $\mu\bar{\mu}$               $s\bar{s}$                  $c\bar{c}$          $t^{(*)}\bar{t}^{(*)}$
  `br. H2_N2HDM_b`           $g g$              $\gamma\gamma$             $Z\gamma$        ${W^\pm}^{(*)}{W^\mp}^{(*)}$      $Z^{(*)}Z^{(*)}$      
  `br. H2_N2HDM_c`         $H_D H_D$               $H_1 H_1$               $A_D A_D$               $Z^{(*)} A_D$           ${W^\pm}^{(*)} H_D^\mp$     $H_D^\pm H_D^\mp$
  `br. H2_N2HDM_d`    $\Gamma_\text{tot}$                                                                                                           
                                                                                                                                                   
  `br. HD_N2HDM`         $Z^{(*)}A_D$       ${W^\pm}^{(*)} H_D^\mp$   $\Gamma_\text{tot}$                                                           
                                                                                                                                                   
  `br. AD_N2HDM`         $Z^{(*)}H_D$       ${W^\pm}^{(*)} H_D^\mp$   $\Gamma_\text{tot}$                                                           
                                                                                                                                                   
  `br. HD+_N2HDM`      ${W^+}^{(*)} H_D$       ${W^+}^{(*)} A_D$      $\Gamma_\text{tot}$                                                           
                                                                                                                                                   
  `br.top`                 $W^+ b$                 $H^+ b$           $\Gamma_\text{tot}$                                                           

  : The branching ratios calculated in the IDP and their location in the output files.  Particles marked with $^{(*)}$ are allowed to be off-shell.[]{data-label="tab:IDP_output"}

Parameter Transformations and Triple-Higgs Couplings in the DSP {#app:DSP}
===============================================================

In this section we present the transformations between the parameter set \[eq:dspinputpars\] and the parameters of \[eq:n2hdmpot\] in the DSP.  We also provide the formulas for the triple-Higgs couplings.  These expressions have first been derived in [@Engeln:2018ywp]. 

Parameter Transformations {#app:transDSP}
-------------------------

We first use the minimum conditions \[eq:n2hdmmin1,eq:n2hdmmin2\] to express the mass parameter $m_{11}^2$ and $m_{22}^2$ as $$\begin{aligned}
m_{11}^2 =&  \tan\beta m_{12}^2 -\frac{v^2}{2} \left(c_\beta^2 \lambda_1 +  s_\beta^2 \lambda_{345}\right)\,,\\
m_{22}^2  =& \cot\beta m_{12}^2-\frac{v^2}{2}\left(c_\beta^2 \lambda_{345} + s_\beta^2 \lambda_2 \right)\,,\end{aligned}$$ where we use the short-hand notations \[eq:l345\], $s_x\equiv\sin x$ and $c_x\equiv\cos x$. To obtain the remaining transformations we solve the linear system given by the rotations into the mass basis.  This yields $$\begin{aligned}
m_{S}^2 =&  -\frac{1}{2}\left(v_1^2 \lambda_7 + v_2^2 \lambda_8 - 2 m_{H_D}\right)\,,\\
\lambda_1 =& \frac{1}{v^2 c^2_\beta}\left[\left(\sum_{i=1}^2 m_{H_{i}}^2 R^2_{i1}\right) -m_{12}^2 \frac{s_\beta}{c_\beta}\right]\,,\\
\lambda_2 =& \frac{1}{v^2 s^2_\beta}\left[\left(\sum_{i=1}^2 m_{H_{i}}^2 R^2_{i2}\right) -m_{12}^2 \frac{c_\beta}{s_\beta}\right]\,,\\
\lambda_3 =& \frac{1}{v^2 c_\beta s_\beta}\left[\left(\sum_{i=1}^2 m_{H_{i}}^2 R_{i1} R_{i2}\right) - m_{12}^2\right] + \frac{2}{v^2} m_{H^\pm}^2\,,\\
\lambda_4 =& \frac{1}{v^2}\left(m_A^2 - 2 m_{H^\pm}^2\right) + \frac{1}{v^2 c_\beta s_\beta} m_{12}^2\,,\\
\lambda_5 =& - \frac{1}{v^2} m_A^2 + \frac{1}{v^2 c_\beta s_\beta}m_{12}^2\,,\end{aligned}$$ with $R_{ij}$ given by \[eq:mixmatdsp\]. 

Triple-Higgs Couplings {#app:tripleHDSP}
----------------------

The triple-Higgs couplings are defined through $$\begin{aligned}
  g(X_i X_j X_k)=\left.\frac{\partial^3\mathcal{L}}{\partial X_i \partial X_j \partial X_k}\right|_\text{VEV}\,,\label{eq:tripleHdef}\end{aligned}$$ where in the DSP $X_{i,j,k} \in \left\{H_1, H_2, A, H^\pm, H_D\right\}$. The couplings have first been derived in [@Engeln:2018ywp].  Using $H_{i,j}\in\{H_1,H_2\}$, $R_{ij}$ from \[eq:mixmatdsp\] and $\lambda_{34-5}\equiv\lambda_3+\lambda_4-\lambda_5$ they are given by $$\begin{aligned}
g(H_iH_iH_i) &=
\begin{aligned}[t]
3v \Big[
&c_\beta \left(
R^3_{i1} \lambda_1 + R_{i1}R^2_{i2} \lambda_{345} \right)\\
&+ s_\beta \left(
R^3_{i2} \lambda_2 + R_{i2}R^2_{i1} \lambda_{345}
\right)
\Big]\,,
\end{aligned}\\
g(H_iH_jH_j) &=
\begin{aligned}[t]
v \Big[
&c_\beta \left(
3 R_{i1} R^2_{j1} \lambda_1
+ (3 R_{i2}  R_{j1} R_{j2} + R_{i1}) \lambda_{345}
\right)\\
&+ v s_\beta \left(
3 R_{i2} R^2_{j2} \lambda_2
+ (3 R_{i1} R_{j1} R_{j2} + R_{i2}) \lambda_{345}
\right)\Big]\,,
\end{aligned}\\
g(H_iAA) &=
\begin{aligned}[t]
v \Big[&
c_\beta\left(
c_\beta s_\beta R_{i2} \left(\lambda_2 -2\lambda_5\right)
+ c^2_\beta R_{i1} \lambda_{34-5}
\right)\\
&+ s_\beta \left(
c_\beta s_\beta R_{i2} \left(\lambda_1 -2 \lambda_5\right)
+ s^2_\beta R_{i2} \lambda_{34-5}
\right)
\Big]\,,
\end{aligned}\\
g(H_iH^+H^-) &=
\begin{aligned}
v \Big[
& c_{\beta} \left(
s^2_\beta R_{i1} \lambda_1
+ c^2_\beta R_{i1} \lambda_3
- c_{\beta} s_{\beta} R_{i2} \left(\lambda_4 + \lambda_5\right)
\right)\\
&+ s_\beta \left(
 c^2_\beta R_{i2} \lambda_2
+ s^2_\beta R_{i2} \lambda_3
- c_{\beta} s_{\beta} R_{i1} \left(\lambda_4 + \lambda_5\right)
\right)
\Big]\,,
\end{aligned}\\
g(H_iH_DH_D) &= v\Big[c_\beta R_{i1} \lambda_7 + s_\beta R_{i2} \lambda_8\Big]\,.\end{aligned}$$ All other trilinear couplings among the $X_i$ vanish. 

Parameter Transformations and Triple-Higgs Couplings in the IDP {#app:IDP}
===============================================================

The relations between the parameters of the Lagrangian \[eq:n2hdmpot\] and the input parameter set \[eq:idpinputpars\] for the IDP as well as the triple-Higgs couplings were derived in [@Engeln:2018ywp].  In [@Engeln:2018ywp] the charges of $\Phi_1$ and $\Phi_2$ under the $\mathbb{Z}_2$ symmetry are switched.  Their convention can be restored from ours by replacing $$\begin{aligned}
  m_{11}^2\leftrightarrow m_{22}^2\,,\quad \lambda_1\leftrightarrow \lambda_2\,,\quad \lambda_7\leftrightarrow\lambda_8\,,\quad R_{i1}\leftrightarrow R_{i2}\end{aligned}$$ in all formulas given below. 

Parameter Transformations {#app:transIDP}
-------------------------

Using the minimum conditions \[eq:n2hdmmin1,eq:n2hdmmin3\] we obtain $$\begin{aligned}
m_{11}^2 &= -\frac{1}{2}\left(v^2\lambda_1+v_S^2\lambda_7\right),\\
m_{S}^2 &= -\frac{1}{2}\left(v_S^2\lambda_6+v^2\lambda_7\right).\end{aligned}$$ We can then solve for the parameters of the Lagrangian and obtain $$\begin{aligned}
\lambda_1 =& \frac{1}{v^2}\left(\sum_{i=1}^2 m_{H_{i}}^2  R^2_{i1}\right)\,,\\
\lambda_3 =& \frac{1}{v^2}\left(2\left(m_{H^{\pm}_D}^2 - m_{11}^2\right) - v_s^2\, \lambda_8\right)\,, \\
\lambda_4 =& \frac{1}{v^2}\left(m_{A_D}^2 + m_{H_{D}}^2 - 2  m_{H^{\pm}}^2\right)\,,\\
\lambda_5 =& \frac{1}{v^2}\left(m_{H_{D}}^{2}-m_{A_D}^{2}\right)\,,\\
\lambda_6 =& \frac{1}{v_s^2}\left(\sum_{i=1}^2 m_{H_{i}}^2  R^2_{i3}\right)\,,\\
\lambda_7 =& \frac{1}{v v_s}\left(\sum_{i=1}^2 m_{H_{i}}^2  R_{i1}  R_{i3}\right)\,,\end{aligned}$$ with $R_{ij}$ given by \[eq:mixmatIDP\]. 

Triple-Higgs Couplings {#app:tripleHIDP}
----------------------

The triple-Higgs couplings of the IDP are defined through \[eq:tripleHdef\] with $X_{i,j,k} \in \left\{H_1, H_2, H_D, A_D, H^\pm_D\right\}$. Due to the $\mathbb{Z}_2$ parity any couplings with an odd number of inert particles vanish.  The non-zero triple-Higgs couplings are given by $$\begin{aligned}
  g(H_i H_i H_i) &=
\begin{aligned}[t]
 &3 \lambda_1  v  R^3_{i1}
                  + 3 \lambda_6  v_s  R^3_{i3}\\
                & + 3 \lambda_7 \left( v R_{i1} R^2_{i3} + v_s   R^2_{i1}R_{i3}   \right)\,,
\end{aligned}\\
g(H_i H_j H_j) &=
\begin{aligned}[t]
 &3 \lambda_1  v   R_{i1}   R^2_{j1}
                  + 3 \lambda_6  v_s   R_{i3}   R^2_{j3} \\
                & + \lambda_7 \left[v   \left( R_{i1}   R^2_{j3} + 2 R_{i3}   R_{j1}   R_{j3}\right)\right.\\
                &  \left.\hphantom{+\lambda_7[} + v_s   \left(R_{i3}   R^2_{j1} + 2 R_{i1}   R_{j1}   R_{j3}\right)\right]\,,
\end{aligned}\\
g(H_i H_D H_D) &=
\frac{2}{v}\, \left(m_{H_{D}}^{2}-m_{22}^{2}\right) R_{i1} + \lambda_8 \, \frac{v_s}{v} \,\left(v    R_{i3} - v_s   R_{i1}\right)\,,\\
g(H_i H_D^{+} H_D^{-}) &= \dfrac{2}{v}\,\left(m_{H^{\pm}_D}^{2}-m_{22}^{2}\right) R_{i1} + \lambda_8\, \dfrac{v_s}{v}\,\left(v  R_{i3}-v_s  R_{i1}\right)\,,\\
g(H_i A_D A_D) &= \,\dfrac{2}{v}\,\left(m_{A_D}^{2}-m_{22}^{2}\right) R_{i1} + \lambda_8\, \dfrac{v_s}{v}\,\left(v  R_{i3}-v_s  R_{i1}\right).\end{aligned}$$ where $H_{i,j}\in\{H_1,H_2\}$ and $R_{ij}$ is given by \[eq:mixmatIDP\]. 

Acknowledgements {#acknowledgements .unnumbered}
================

We are grateful to Marco Sampaio and Rui Santos for helpful discussions and feedback.  We thank Michael Spira for providing and maintaining the code `HDECAY` that `N2HDECAY` is based upon.  JW acknowledges funding from the “PIER Helmholtz Graduate School”. 

Bibliography {#bibliography .unnumbered}
============
