Skip to content

GitLab

  • Projects
  • Groups
  • Snippets
  • Help
    • Loading...
  • Help
    • Help
    • Support
    • Community forum
    • Submit feedback
    • Contribute to GitLab
  • Sign in / Register
SARAH SARAH
  • Project overview
    • Project overview
    • Details
    • Activity
  • Packages & Registries
    • Packages & Registries
    • Container Registry
  • Analytics
    • Analytics
    • Repository
    • Value Stream
  • Wiki
    • Wiki
  • Members
    • Members
  • Activity
Collapse sidebar
  • GOODSELL Mark
  • SARAHSARAH
  • Wiki
  • Vertices

Last edited by Martin Gabelmann Jun 28, 2019
Page history

Vertices

Vertices

Calculation of vertices

Vertices are not automatically calculated during the initialization of a model like this is done for mass matrices and tadpole equations. However, the calculation can be started very easily. In general, SARAH is optimized for the extraction of three- and four-point interactions with renormalizable operators. That means, usually only the following generic interactions are taken into account in the calculations: interactions of two fermions or two ghosts with one scalar or vector bosons (FFS, FFV, GGS, GGV ), interactions of three or four scalars or vector bosons (SSS, SSSS, VVV, VVVV ), as well as interactions of two scalars with one or two vector bosons (SSV, SSVV ) or two vector bosons with one scalar (SVV ). In this context, vertices not involving fermions are calculated by

V(\eta_a,\eta_b,\eta_c) \equiv i \frac{\partial^3 \mathfrak{L}}{\partial \eta_a \partial \eta_b \partial \eta_c} = C \Gamma \\ V(\eta_a,\eta_b,\eta_c, \eta_d) \equiv= i \frac{\partial^4 \mathfrak{L}}{\partial \eta_a \partial \eta_b \partial \eta_c \partial \eta_d} = C \Gamma

Here,ηi are either scalars, vector bosons, or ghosts. The results are expressed by a coefficientC which is a Lorentz invariant and a Lorentz factorΓ which involvesγμ,pμ, orgμ**ν. Vertices for Dirac fermions are first expressed in terms of Weyl fermions. The two vertices are then calculated separately. Taking two Dirac fermionsFa = (Ψa1, Ψa2*),Fb = (Ψb1, Ψb2*) and distinguishing the two cases for fermion–vector and fermion–scalar couplings, the vertices are calculated and expressed by

V(\bar F_a,F_b,V_c) = \{V(\Psi_a^{1 \*},\Psi_b^1,V_c), V(\Psi_a^2,\Psi_b^{2\*},V_c)\} \equiv \{C^L \gamma_\mu P_L, C^R \gamma_\mu P_R\} \\ V(\bar F_a,F_b,S_c) = \{V(\Psi_a^2,\Psi_b^1,S_c), V(\Psi_a^{1\*},\Psi_b^{2\*},V_c)\} \equiv \{C^L P_L, C^R P_R\}

Here, the polarization operatorsPL, R are used. The user can either calculate specific vertices for a particular set of external states or call functions that SARAH derives all existing interactions from the Lagrangian. The first option might be useful to check the exact structure of single vertices, while the second one is needed to get all vertices to write model files for other tools.

Calculating specific vertices with SARAH

Vertices are calculated by

Vertex[ParticleList,Options]

ParticleList is a list containing the involved fields. The following Options are supported by the Vertex command:

  1. Eigenstates, Value: $EIGENSTATES, Default: Last entry in NameOfState Fixes the considered eigenstates
  2. UseDependences: Value True or False, Default: False Optional relations between the parameters (see section [parameterFile]) will be used, if UseDependences is set to True.

The output of Vertex is an array:

{{ParticleList},{{Coefficient 1, Lorentz 1},{Coefficient 2, Lorentz 2},...}

First, the list of the involved particles is given and the indices are inserted. The second part consists of the value of the vertex and can be also a list, if different Lorentz structures are possible. In the part independent of any Lorentz index the following symbols can appear

  1. Delta[a,b]: Kronecker deltaδα**β
  2. ThetaStep[i,j]: Step functionΘi**j
  3. Lam[t,a,b]: Gell-Mann matrixλα**βt
  4. LambdaProd[x,y][a,b]: Matrix product of two Gell-Mann matrices(λxλy)α**β
  5. Sig[t,a,b]: Pauli matrixσα**βt
  6. SigmaProd[x,y][a,b]: Matrix product of two Pauli matrices(σxσy)α**β
  7. fSU3[i,j,k]: Structure constants ofS**U(3) fijk
  8. fSU2[i,j,k]: Structure constants ofS**U(2) ϵijk
  9. FST[SU[N]][i,j,k]: Structure constants ofS**U(N) fNijk
  10. TA[SU[N]][a,i,j]: Generator ofS**U(N) Ti**ja
  11. Couplings, e.g. g1, g2, g3, Ye[a,b], Yd[a,b], Yu[a,b], …
  12. Mixing matrices, e.g. ZD[a,b], …

The part transforming under the Lorentz group can consist of

  1. gamma[lor]: Gamma matrixγμ
  2. g[lor1,lor2]: Metric tensorgμ**ν
  3. Mom[particle,lor]: MomentumpPμ of particleP
  4. PL, PR: Polarization operatorsP_L = \frac{1 - \gamma_5}{2},P_R = \frac{1+\gamma_5}{2}
  5. 1: If the vertex is a Lorentz scalar.
  6. LorentzProduct[_,_]: A non commutative product of terms transforming under the Lorentz group

Examples

Some examples to clarify the usage and output of Vertex:

  1. One possible Lorentz structure: Vertex[{hh,Ah,Z}] leads to the vertex of scalar and a pseudo scalar Higgs with aZ-boson {{hh[{gt1}], Ah[{gt2}], VZ[{lt3}]}, {((MA[gt2,1]MH[gt1,1] - MA[gt2,2]MH[gt1,2])(g2Cos[ThetaW]+g1*Sin[ThetaW]))/2, Mom[Ah[{gt2}], lt3] - Mom[hh[{gt1}],lt3]}}

    The output is divided in two parts. First, the involved particles are given, second, the value of the vertex is given. This second part is again split in two parts: one is the Lorentz independent part and the second part defines the transformation under the Lorentz group.

  2. Several possible Lorentz structures Vertex[{bar[Fd],Fd,hh}] is the interaction between d-quarks and a Higgs: {{bar[Fd[{gt1, ct1}]], Fd[{gt2, ct2}], hh[{gt3}]}, {((-I)*Delta[ct1,ct2]*Delta[gt1,gt2]*MH[gt3,2]*Yd[gt2,gt1])/Sqrt[2],PL}, {((-I)*Delta[ct1,ct2]*Delta[gt1,gt2]*MH[gt3,2]*Yd[gt1,gt2])/Sqrt[2],PR}}

    Obviously, there are three parts: one for the involved particles and two for the different Lorentz structures. PL and PR are the polarization projectors$P_L = \frac{1}{2} (1 - \gamma_5), P_R = \frac{1}{2} (1 + \gamma_5)$.

  3. Changing the considered eigenstates and using Weyl fermions It is also possible to calculate the vertices for Weyl fermions and/or to consider the gauge eigenstates. For instance, Vertex[{fB, FdL, conj[SdL]}, Eigenstates -> GaugeES]

    returns

    {{fB, FdL[{gt2, ct2}], conj[SdL[{gt3, ct3}]]},
     {((-I/3)*g1*Delta[ct2, ct3]*Delta[gt2, gt3])/Sqrt[2],1}}
  4. Using dependences With Vertex[{conj[Se], Se, VP}, UseDependences -> True]g1 andg2 are replaced by the electric chargee. This and similar relations can be defined in parameters.m). {{conj[Se[{gt1}]], Se[{gt2}], VP[{lt3}]}, {(-I)eDelta[gt1,gt2],-Mom[conj[Se[{gt1}]],lt3]+Mom[Se[{gt2}],lt3]}}

  5. Fixing the generations It is possible to give the indices of the involved particles already as input Vertex[{hh[{1}], hh[{1}], Ah[{2}], Ah[{2}]}]

    leads to

    {{hh[{1}], hh[{1}], Ah[{2}], Ah[{2}]},
     {(-I/4)*(g1^2 + g2^2)*Cos[2*\[Alpha]]*Cos[2*\[Beta]], 1}}

    Obviously, the given definition of the mixing matrices for the Higgs fields were automatically inserted. If the indices are fixed by a replacement, the definition of the mixing matrix wouldn’t be used

    Vertex[{hh, hh, Ah, Ah}] /. {gt1->1, gt2->1,gt3->2, gt3->2}

    returns

    {{hh[{1}], hh[{1}], Ah[{2}], Ah[{gt4}]},
     {(-I/4)*(g1^2 + g2^2)*(conj[ZA[2, 1]]*conj[ZA[gt4, 1]] -
       conj[ZA[2, 2]]*conj[ZA[gt4, 2]])*(conj[ZH[1, 1]]^2 - conj[ZH[1, 2]]^2), 1}}

    However,

    Vertex[{hh, hh, Ah, Ah}] /. {gt1->1, gt2->1,gt3->2, gt3->2} /.subAlways

    leads to the former expression using the mixing angleβ.

  6. Effective operators In effective theories also interactions between two fermions and two scalars are possible. As example an effective vertex for a model in which the gluino was integrated out: Vertex[{Fd, Fd, conj[Sd], conj[Sd]}]

    Returns

    {{Fd[{gt1, ct1}], Fd[{gt2, ct2}], conj[Sd[{gt3, ct3}]], conj[Sd[{gt4, ct4}]]},
    {-(g3^2*(sum[j1, 1, 8, (Lam[j1, ct3, ct2]*Lam[j1, ct4, ct1])/Mass[fG][j1]]*
               ZD[gt3, gt2]*ZD[gt4, gt1] +
             sum[j1, 1, 8, (Lam[j1, ct3, ct1]*Lam[j1, ct4, ct2])/Mass[fG][j1]]*
               ZD[gt3, gt1]*ZD[gt4, gt2])),
        LorentzProduct[PL, PL]}, {0, LorentzProduct[PR, PL]},
    {g3^2*(sum[j1, 1, 8, (Lam[j1, ct2, ct3]*Lam[j1, ct4, ct1])/Mass[fG][j1]]*
           ZD[gt3, 3 + gt2]*ZD[gt4, gt1] +
           sum[j1, 1, 8, (Lam[j1, ct2, ct4]*Lam[j1, ct3, ct1])/Mass[fG][j1]]*
           ZD[gt3, gt1]*ZD[gt4, 3 + gt2]),
        LorentzProduct[PL, PR]}, {0, LorentzProduct[PR, PR]},
    {0, LorentzProduct[gamma,PL, PL]}, {0, LorentzProduct[gamma, PR, PL]},
    {0, LorentzProduct[gamma, PL, PR]}, {0, LorentzProduct[gamma,PR, PR]}}

    Obviously, SARAH checks the eight possible operators (4 different combination of polarization operators with and without aγ matrix) and returns the result for each operator.

Calculating all vertices at once

To calculate all vertices at once for a given model, use

MakeVertexList[Eigenstates,Options]

First, the name of the eigenstates has to be given. The possible options are:

  1. effectiveOperators, Values: True or False, Default: False If also higher dimensional operators should be calculated. By default, this concerns only four point interactions.
  2. SixParticleInteractions, Values: True or False, Default: False If also the six-point interactions should be calculated.
  3. GenericClasses, Values: All or a list of generic types, Default: All Calculates the vertices only for the given types of interaction

The results are saved in list named

 SA`VertexList[Type]

with Type = { SSS, SSSS, SSVV, SSV, SVV, FFS, FFV, VVV, VVVV, GGS, GGV, ASS}.

Output

The vertices are exported into many different formats. They are saved in the SARAH internal format and they can be written to LaTeX files. The main purpose is the export into formats which can be used with other tools. SARAH writes model files for

  • FeynArts
  • WHIZARD
  • CalcHep/CompHep
  • in the UFO format.

The UFO format is supported by MadGraph, Herwigg+ and Sherpa. Thus, by the output of the vertices into these different format, SARAH provides an implementation of a given model in a wide range of HEP tools. In addition, SARAH generates also Fortran-code to implement all vertices in SPheno.

See also

Clone repository

Home

Index

  • Additional terms in Lagrangian
  • Advanced usage of FlavorKit
  • Advanced usage of FlavorKit to calculate new Wilson coefficients
  • Advanced usage of FlavorKit to define new observables
  • Already defined Operators in FlavorKit
  • Already defined observables in FlavorKit
  • Auto-generated templates for particles.m and parameters.m
  • Automatic index contraction
  • Basic definitions for a non-supersymmetric model
  • Basic definitions for a supersymmetric model
  • Basic usage of FlavorKit
  • Boundary conditions in SPheno
  • CalcHep CompHep
  • Calculation of flavour and precision observables with SPheno
  • Checking the particles and parameters within Mathematica
  • Checks of implemented models
  • Conventions
  • Decay calculation with SPheno
  • Defined FlavorKit parameters
  • Definition of the properties of different eigenstates
  • Delete Particles
  • Different sets of eigenstates
  • Diphoton and digluon vertices with SPheno
  • Dirac Spinors
  • FeynArts
  • Fine-Tuning calculations with SPheno
  • Flags for SPheno Output
  • Flags in SPheno LesHouches file
  • FlavorKit
  • FlavorKit Download and Installation
  • Flavour Decomposition
  • GUT scale condition in SPheno
  • Gauge Symmetries SUSY
  • Gauge Symmetries non-SUSY
  • Gauge fixing
  • Gauge group constants
  • General information about Field Properties
  • General information about model implementations
  • Generating files with particle properties
  • Generic RGE calculation
  • Global Symmetries SUSY
  • Global Symmetries non-SUSY
  • Handling of Tadpoles with SPheno
  • Handling of non-fundamental representations
  • HiggsBounds
  • Higher dimensionsal terms in superpotential
  • Input parameters of SPheno
  • Installation
  • Installing Vevacious
  • LHCP
  • LHPC
  • LaTeX
  • Lagrangian
  • Loop Masses
  • Loop calculations
  • Loop functions
  • Low or High scale SPheno version
  • Main Commands
  • Main Model File
  • Matching to the SM in SPheno
  • MicrOmegas
  • ModelOutput
  • Model files for Monte-Carlo tools
  • Model files for other tools
  • Models with Thresholds in SPheno
  • Models with another gauge group at the SUSY scale
  • Models with several generations of Higgs doublets
  • More precise mass spectrum calculation
  • No SPheno output possible
  • Nomenclature for fields in non-supersymmetric models
  • Nomenclature for fields in supersymmetric models
  • One-Loop Self-Energies and Tadpoles
  • One-Loop Threshold Corrections in Scalar Sectors
  • Options SUSY Models
  • Options non-SUSY Models
  • Parameters.m
  • Particle Content SUSY
  • Particle Content non-SUSY
  • Particles.m
  • Phases
  • Potential
  • Presence of super-heavy particles
  • RGE Running with Mathematica
  • RGEs
  • Renormalisation procedure of SPheno
  • Rotations angles in SPheno
  • Rotations in gauge sector
  • Rotations in matter sector
  • SARAH in a Nutshell
  • SARAH wiki
  • SLHA input for Vevacious
  • SPheno
  • SPheno Higgs production
  • SPheno Output
  • SPheno and Monte-Carlo tools
  • SPheno files
  • SPheno mass calculation
  • SPheno threshold corrections
  • Setting up SPheno.m
  • Setting up Vevacious
  • Setting up the SPheno properties
  • Special fields and parameters in SARAH
  • Superpotential
  • Support of Dirac Gauginos
  • Supported Models
  • Supported gauge sectors
  • Supported global symmetries
  • Supported matter sector
  • Supported options for symmetry breaking
  • Supported particle mixing
  • Tadpole Equations
  • The renormalisation scale in SPheno
  • Tree-level calculations
  • Tree Masses
  • Two-Loop Self-Energies and Tadpoles
  • UFO
  • Usage of tadpoles equations
  • Using SPheno for two-loop masses
  • Using auxiliary parameters in SPheno
  • VEVs
  • Vertices
  • Vevacious
  • WHIZARD