Higgs to two photon

Initialization

Quit[];

$Language="English";
SetDirectory[NotebookDirectory[]];
<<FeynArts`
<<FormCalc`

h -> γ γ Matrix Element Calculation

topo1=CreateTopologies[1,1->2,ExcludeTopologies->{Internal}];
diag1=InsertFields[topo1,{S[1]}->{V[1],V[1]},InsertionLevel->{Particles},Restrictions->NoLightFHCoupling];
Paint[diag1];

ClearProcess[];
amp=CalcFeynAmp[CreateFeynAmp[diag1],FermionChains->Chiral]

To confirm that this diagram converges...

UVDivergentPart[amp]
divergence=%[[1]]

PolarizationSum directory to the amplitude invokes SquaredME, HelicityME, and ColourME (with _Hel=0 option).
[but you should not do so because you will lose traceability; instead I recommend you to do step by step as before.]

mat=PolarizationSum[amp,GaugeTerms->False]//.Subexpr[]//.Abbr[]//Simplify

Numerical Evaluation

Install["LoopTools"]

VALUES[higgs_]:={
Den[p2_,m2_]:>1/(p2-m2),
MT->173,MT2->(173^2),
MW->80.4,MW2->(80.4^2),
Alfa->1/128,Alfa2->1/(128^2),
SW2->0.23,
MH->higgs,MH2->(higgs^2)}
MyMatrixElement[higgs_]:=mat//.VALUES[higgs]//Simplify

MyMatrixElement[125]

Compare to a literature

f[t_]:=If[t<1,-(1/4)(((Log[((1+(Sqrt[1-t]))/(1-(Sqrt[1-t])))]-I π))^2),(((ArcSin[(Sqrt[1/t])]))^2)]
Fb[t_]:=2+3t+3t(2-t)f[t]
Ff[t_]:=-2t(1+(1-t)f[t])
Fs[t_]:=t(1-t f[t])
InSomePaper[higgs_]:=((Alfa (MH2^2))/(8π SW2 MW2))((Abs[Alfa*3*((((2/3)))^2)Ff[(((((2MT)/MH)))^2)]+Alfa*1*(1^2)Fb[(((((2MW)/MH)))^2)]])^2)//.VALUES[higgs]

Plot[{MyMatrixElement[mh],InSomePaper[mh]},{mh,0,200}]

Plot[{MyMatrixElement[mh]/InSomePaper[mh]},{mh,0,200},PlotRange->{{0,200},{0,2}}]

The partial decay rate will be:

MyDecayRate[higgs_]:=(1/(32π higgs))MyMatrixElement[higgs]

MyDecayRate[125]

Compare this value with literature!

Of course we can (and should) ascertain the gauge independence. Here we utilize PaVeReduce.

ClearProcess[];
amp=CalcFeynAmp[CreateFeynAmp[diag1],FermionChains->Chiral,PaVeReduce->True]//.Subexpr[]//.Abbr[]
matWithGaugeTerms=PolarizationSum[amp,GaugeTerms->True]//.Subexpr[]//.Abbr[]//.Finite->1//FullSimplify