a Code for the Combination of Indirect and Direct Constraints on High Energy Physics Models Logo
Public Member Functions | Private Attributes | List of all members
GeneralTHDMMatching Class Reference

#include <GeneralTHDMMatching.h>

+ Inheritance diagram for GeneralTHDMMatching:

Detailed Description

Author
HEPfit Collaboration
Copyright
GNU General Public License

Definition at line 25 of file GeneralTHDMMatching.h.

Public Member Functions

virtual double C10Bll (double xt, double xHp, gslpp::complex su)
 
virtual std::vector< WilsonCoefficient > & CMBMll (QCD::lepton lepton)
 
virtual std::vector< WilsonCoefficient > & CMbsg ()
 operator basis: current current; qcd penguins; magnetic and chromomagnetic penguins; semileptonic More...
 
virtual std::vector< WilsonCoefficient > & CMdbs2 ()
 Calculates the function of Eq. (68) of 1607.06292. More...
 
virtual std::vector< WilsonCoefficient > & CMdiujleptonknu (int i, int j, int k)
 
virtual std::vector< WilsonCoefficient > & CMgminus2mu ()
 Wilson coefficient for \( (g-2)_{\mu} \). More...
 
virtual gslpp::complex CPboxBll (double xt, double xHp, gslpp::complex su, gslpp::complex sd, gslpp::complex sl)
 
virtual gslpp::complex CphiU (double xHp, double xt, double vev, double xphi, double mu, double Ri1, double Ri2, double Ri3, double mHi_2, double lambda3, double Relambda7, double Imlambda7, gslpp::complex su, gslpp::complex sd)
 
virtual gslpp::complex CPZUBll (double xt, double xHp, double sW2, gslpp::complex su, gslpp::complex sd)
 
virtual gslpp::complex CSboxBll (double xt, double xHp, gslpp::complex su, gslpp::complex sd, gslpp::complex sl)
 
virtual double f1 (double xHp, double xt)
 
virtual double f10 (double xHp, double xt)
 
double F1oneloopgm2 (const double ratio_sq)
 Loop \( (g-2)_{\mu} \) at one loop from 1502.04199. The complete loop function is \int_0^1 dx \frac{x^2(2-x)}{ratio_sq*x^2-x+1} but for small values of ratio_sq we can expand and take only the first term for which the functions becomes (-7.0/6.0-log(ratio_sq)) More...
 
double F1tildetwoloopgm2 (const double ratio_sq)
 Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{ratio_sq}{2}*\int_0^1 dx \frac{1}{ratio_sq-x(1-x)}\log(\frac{ratio_sq}{x(1-x)}) More...
 
double F1twoloopgm2 (const double ratio_sq)
 Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{ratio_sq}{2}*\int_0^1 dx \frac{2x(1-x)-1}{ratio_sq-x(1-x)}\log(\frac{ratio_sq}{x(1-x)}) More...
 
virtual double f2 (double xHp, double xt)
 
double F2oneloopgm2 (const double ratio_sq)
 Loop \( (g-2)_{\mu} \) at one loop from 1502.04199. The complete loop function is \int_0^1 dx \frac{-x^3}{ratio_sq*x^2-x+1} but for small values of ratio_sq we can expand and take only the first term for which the functions becomes (11.0/6.0+log(ratio_sq)) More...
 
double F2twoloopgm2 (const double ratio_sq)
 Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{1}{2}*\int_0^1 dx \frac{x(x-1)}{ratio_sq-x(1-x)}\log(\frac{ratio_sq}{x(1-x)}) More...
 
virtual double f3 (double xHp, double xt)
 
double F3oneloopgm2 (const double ratio_sq)
 Loop \( (g-2)_{\mu} \) at one loop from 1502.04199. The complete loop function is \int_0^1 dx \frac{x^2(1-x)}{ratio_sq*x*(1-x)-x} but for small values of ratio_sq we can expand and take only the first term for which the functions becomes (-1/12. - ratio_sq/60.) More...
 
double F3twoloopgm2 (const double ratio_sq)
 Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{1}{2}*\int_0^1 dx \frac{x(3x(4x-1)+10)ratio_sq-x(1-x)}{ratio_sq-x(1-x)}\log(\frac{ratio_sq}{x(1-x)}) More...
 
virtual double f4 (double xHp, double xt)
 
double F4twoloopgm2 (const double ratio_sq)
 Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. In the notation of [1502.04199] this corresponds to (Q_t*x+Q_b*(1-x))x(1+x)G(ratio_sq,0), i.e. we neglect mb/mHp and mb/mW which is a great approximation The complete loop function is \int_0^1 dx (2/3 x-1/3 *(1-x))x(1+x)\frac{\log{\frac{ratio_sq*x}{x(1-x)}}}{x(1-x-ratio_sq*x)}) More...
 
virtual double f5 (double xHp, double xt)
 
double F5twoloopgm2 (const double ratio_sq)
 Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. In the notation of [1502.04199] this corresponds to (Q_t*x+Q_b*(1-x))x(1-x)G(ratio_sq,0), i.e. we neglect mb/mHp and mb/mW which is a great approximation The complete loop function is \int_0^1 dx (2/3 x-1/3 *(1-x))x(1-x)\frac{\log{\frac{ratio_sq*x}{x(1-x)}}}{x(1-x-ratio_sq*x)}) More...
 
virtual double f6 (double xHp, double xt)
 
virtual double f7 (double xHp, double xt)
 
virtual double f8 (double xHp, double xt)
 
virtual double f9 (double xHp, double xt)
 
virtual gslpp::complex g0 (double xHp, double xt, gslpp::complex su, gslpp::complex sd)
 
virtual gslpp::complex g1a (double xHp, double xt, gslpp::complex su, gslpp::complex sd)
 
virtual gslpp::complex g2a (double xHp, double xt, gslpp::complex su, gslpp::complex sd)
 
virtual gslpp::complex g3a (double xHp, double xt, gslpp::complex su, gslpp::complex sd)
 
 GeneralTHDMMatching (const GeneralTHDM &GeneralTHDM_i)
 
const Polylogarithms getPolyLog () const
 
virtual double gminus2muLO ()
 Calculates amplitudes for \( (g-2)_{\mu} \) at one loop from 1502.04199, before [Broggio:2014mna] was used. More...
 
virtual double gminus2muNLO ()
 Calculates amplitudes for \( (g-2)_{\mu} \) at approximate two-loop from [Broggio:2014mna]. More...
 
virtual gslpp::complex lambdaHHphi (double lambda3, double Relambda7, double Imlambda7, double Ri1, double Ri2, double Ri3)
 
virtual gslpp::complex neglog (gslpp::complex argument)
 Calculates the log of a negative number. More...
 
virtual gslpp::complex negpow (double basis, double exp)
 Calculates the power root of a negative number. More...
 
virtual gslpp::complex negsquareroot (double x)
 Calculates amplitudes for \( (g-2)_{\mu} \) at approximate two-loop from [Cherchiglia:2016eui]. More...
 
gslpp::complex setWCbsg (int i, gslpp::complex sigmau, gslpp::complex sigmad, double mt, double mhp, double mu, orders order)
 
void updateGTHDMParameters ()
 
- Public Member Functions inherited from StandardModelMatching
virtual std::vector< WilsonCoefficient > & CMdbd2 ()
 \( \Delta B = 2 \), \( B_{d} \) More...
 
virtual std::vector< WilsonCoefficient > & CMdd2 ()
 \( \Delta C = 2 \), More...
 
 StandardModelMatching (const StandardModel &SM_i)
 
void updateSMParameters ()
 Updates to new Standard Model parameter sets. More...
 
virtual ~StandardModelMatching ()
 
- Public Member Functions inherited from ModelMatching
virtual std::vector< WilsonCoefficient > & CMbnlep (int a)=0
 
virtual std::vector< WilsonCoefficient > & CMbnlepCC (const int a)=0
 
virtual std::vector< WilsonCoefficient > & CMBXsnn (QCD::lepton lepton)=0
 
virtual std::vector< WilsonCoefficient > & CMd1 ()=0
 
virtual std::vector< WilsonCoefficient > & CMd1Buras ()=0
 
virtual std::vector< WilsonCoefficientNew > & CMDF1 (std::string blocks, unsigned int nops)=0
 
virtual std::vector< WilsonCoefficient > & CMprimeBMll (QCD::lepton lepton)=0
 
virtual std::vector< WilsonCoefficient > & CMprimebsg ()=0
 
virtual ~ModelMatching ()
 

Private Attributes

gslpp::complex CWbsgArrayLO [8]
 
gslpp::complex CWbsgArrayNLO [8]
 
gslpp::complex CWbsgArrayNNLO [8]
 
double GF
 
WilsonCoefficient mcBMll
 
WilsonCoefficient mcbsg
 
WilsonCoefficient mcbsmm
 
WilsonCoefficient mcdbs2
 
WilsonCoefficient mcgminus2mu
 
WilsonCoefficient mculeptonnu
 
double mhpbsg
 
double mMU
 
double mtbsg
 
double mubsg
 
gslpp::matrix< gslpp::complex > myCKM
 
const GeneralTHDM & myGTHDM
 
const Polylogarithms PolyLog
 
gslpp::complex sd
 
gslpp::complex sl
 
gslpp::complex su
 

Constructor & Destructor Documentation

◆ GeneralTHDMMatching()

GeneralTHDMMatching::GeneralTHDMMatching ( const GeneralTHDM &  GeneralTHDM_i)

Definition at line 15 of file GeneralTHDMMatching.cpp.

15 :
16
17StandardModelMatching(GeneralTHDM_i),
18myGTHDM(GeneralTHDM_i),
19myCKM(3, 3, 0.),
20mcdbs2(5, NDR, NLO),
21mculeptonnu(5, NDR, LO),
22mcBMll(13, NDR, NLO),
23mcbsg(8, NDR, NNLO),
24mcgminus2mu(2, NDR, NLO),
25mcbsmm(8, NDR, NNLO, NLO_QED22) {
26}
NNLO
@ NNLO
Definition: OrderScheme.h:36
LO
@ LO
Definition: OrderScheme.h:34
NLO
@ NLO
Definition: OrderScheme.h:35
NDR
@ NDR
Definition: OrderScheme.h:21
NLO_QED22
@ NLO_QED22
Definition: OrderScheme.h:63
GeneralTHDMMatching::mculeptonnu
WilsonCoefficient mculeptonnu
Definition: GeneralTHDMMatching.h:370
GeneralTHDMMatching::myGTHDM
const GeneralTHDM & myGTHDM
Definition: GeneralTHDMMatching.h:367
GeneralTHDMMatching::mcBMll
WilsonCoefficient mcBMll
Definition: GeneralTHDMMatching.h:370
GeneralTHDMMatching::mcgminus2mu
WilsonCoefficient mcgminus2mu
Definition: GeneralTHDMMatching.h:370
GeneralTHDMMatching::mcbsg
WilsonCoefficient mcbsg
Definition: GeneralTHDMMatching.h:370
GeneralTHDMMatching::mcbsmm
WilsonCoefficient mcbsmm
Definition: GeneralTHDMMatching.h:370
GeneralTHDMMatching::myCKM
gslpp::matrix< gslpp::complex > myCKM
Definition: GeneralTHDMMatching.h:369
GeneralTHDMMatching::mcdbs2
WilsonCoefficient mcdbs2
Definition: GeneralTHDMMatching.h:370
StandardModelMatching::StandardModelMatching
StandardModelMatching(const StandardModel &SM_i)
Definition: StandardModelMatching.cpp:16

Member Function Documentation

◆ C10Bll()

double GeneralTHDMMatching::C10Bll ( double  xt,
double  xHp,
gslpp::complex  su 
)
virtual
Returns
C10 Wilson coefficient of \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1067 of file GeneralTHDMMatching.cpp.

1067 {
1068
1069 double C10 = su.abs2() * xt * xt / 8 * (1 / (xHp - xt) + xHp / ((xHp - xt)*(xHp - xt))*(log(xt) - log(xHp)));
1070 return C10;
1071}

◆ CMBMll()

std::vector< WilsonCoefficient > & GeneralTHDMMatching::CMBMll ( QCD::lepton  lepton)
virtual

Implements ModelMatching.

Definition at line 1268 of file GeneralTHDMMatching.cpp.

1268 {
1269 updateGTHDMParameters();
1270
1271
1272 //From 1404.5865
1273
1274 //complex i
1275 gslpp::complex i = gslpp::complex::i();
1276
1277
1278 double Muw = myGTHDM.getMuw();
1279 // double Mut = myGTHDM.getMut();
1280 double mHp2 = myGTHDM.getmHp2();
1281 double MW = myGTHDM.Mw();
1282 double Mt_muw = myGTHDM.Mrun(Muw, myGTHDM.getQuarks(QCD::TOP).getMass_scale(),
1283 myGTHDM.getQuarks(QCD::TOP).getMass(), QCD::TOP, FULLNNLO);
1284 /*double mt_mt = myGTHDM.Mrun(Mut, myGTHDM.getQuarks(QCD::TOP).getMass_scale(),
1285 myGTHDM.getQuarks(QCD::TOP).getMass(), FULLNNLO);*/
1286 double mb = myGTHDM.getQuarks(QCD::BOTTOM).getMass();
1287
1288 double ml = myGTHDM.getLeptons(lepton).getMass();
1289
1290 // mlep = SM.getLeptons(lep).getMass();
1291
1292
1293 double xt = (Mt_muw * Mt_muw) / (MW * MW);
1294 double xHp = (mHp2) / (MW * MW);
1295 double vev = myGTHDM.v();
1296 double sW2 = myGTHDM.sW2();
1297 double mHl = myGTHDM.getMHl();
1298 double mH1_2 = mHl*mHl;
1299
1300
1301 //mu contains the missalignemtn dependece. It should be mu -> CR(mu0) - log(mu/mu0). Eq (22)
1302
1303
1304
1305 // mu relates the high alignment scale (\Lambda), not know (which can be until the Plank scale), to the electroweak scale
1306 //double mu = log(Lambda/MW);
1307 // Lambda is not a paraemeter of the model so now mu is set to 0 -> the model is aligned at the ew scale
1308 double mu = 0;
1309
1310 double Imlambda7 = myGTHDM.getMyGTHDMCache()->Imlambda7;
1311 double Relambda7 = myGTHDM.getMyGTHDMCache()->Relambda7;
1312 double lambda3 = myGTHDM.getMyGTHDMCache()->lambda3;
1313
1314 double R11 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(0,0);
1315 double R12 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(0,1);
1316 double R13 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(0,2);
1317 double R21 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(1,0);
1318 double R22 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(1,1);
1319 double R23 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(1,2);
1320 double R31 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(2,0);
1321 double R32 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(2,1);
1322 double R33 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(2,2);
1323
1324 gslpp::complex sl = myGTHDM.getNl_11();
1325 gslpp::complex su = myGTHDM.getNu_11();
1326 gslpp::complex sd = myGTHDM.getNd_11();
1327
1328
1329
1330 //Mass of the physical scalars
1331
1332 double mH2_2 = myGTHDM.getmH2sq();
1333 double mH3_2 = myGTHDM.getmH3sq();
1334
1335
1336 double xphi1 = mH1_2 / (MW * MW);
1337 double xphi2 = mH2_2 / (MW * MW);
1338 double xphi3 = mH3_2 / (MW * MW);
1339
1340 //Yukawa couplings. Eq. (19)
1341
1342 gslpp::complex yl1 = R11 + (R12 + i * R13) * sl;
1343 gslpp::complex yl2 = R21 + (R22 + i * R23) * sl;
1344 gslpp::complex yl3 = R31 + (R32 + i * R33) * sl;
1345
1346 gslpp::complex CSboxU = CSboxBll(xt, xHp, su, sd, sl);
1347 gslpp::complex CPboxU = CPboxBll(xt, xHp, su, sd, sl);
1348 gslpp::complex CPZU = CPZUBll(xt, xHp, sW2, su, sd);
1349
1350 gslpp::complex CSphi1U = yl1.real() * CphiU(xHp, xt, vev, xphi1, mu, R11, R12, R13, mH1_2, lambda3, Relambda7, Imlambda7, su, sd);
1351 gslpp::complex CSphi2U = yl2.real() * CphiU(xHp, xt, vev, xphi2, mu, R21, R22, R23, mH2_2, lambda3, Relambda7, Imlambda7, su, sd);
1352 gslpp::complex CSphi3U = yl3.real() * CphiU(xHp, xt, vev, xphi3, mu, R31, R32, R33, mH3_2, lambda3, Relambda7, Imlambda7, su, sd);
1353
1354 gslpp::complex CPphi1U = i * yl1.imag() * CphiU(xHp, xt, vev, xphi1, mu, R11, R12, R13, mH1_2, lambda3, Relambda7, Imlambda7, su, sd);
1355 gslpp::complex CPphi2U = i * yl2.imag() * CphiU(xHp, xt, vev, xphi2, mu, R21, R22, R23, mH2_2, lambda3, Relambda7, Imlambda7, su, sd);
1356 gslpp::complex CPphi3U = i * yl3.imag() * CphiU(xHp, xt, vev, xphi3, mu, R31, R32, R33, mH3_2, lambda3, Relambda7, Imlambda7, su, sd);
1357
1358
1359 //Total 2HDM Wilson coefficients CS and CP PART. Eq. (31)-(33) without SM part
1360
1361 gslpp::complex CSphiU = CSboxU + CSphi1U + CSphi2U + CSphi3U;
1362 gslpp::complex CPphiU = CPboxU + CPZU + CPphi1U + CPphi2U + CPphi3U;
1363
1364
1365 vmcBMll = StandardModelMatching::CMBMll(lepton);
1366 switch (mcbsg.getScheme()) {
1367 case NDR:
1368
1369 break;
1370 default:
1371 std::stringstream out;
1372 out << mcBMll.getScheme();
1373 throw std::runtime_error("GeneralTHDMMatching::CMBMll(): scheme " + out.str() + "not implemented");
1374 }
1375 mcBMll.setMu(Muw);
1376
1377 switch (mcBMll.getOrder()) {
1378 case NNLO:
1379 case NLO:
1380 case LO:
1381 mcBMll.setCoeff(9, C10Bll(xt, xHp, su) / (sW2), LO);
1382 mcBMll.setCoeff(10, (CSphiU * mb * ml) / (MW * MW * sW2), LO);
1383 mcBMll.setCoeff(11, (CPphiU * mb * ml) / (MW * MW * sW2), LO);
1384 break;
1385 default:
1386 std::stringstream out;
1387 out << mcBMll.getOrder();
1388 throw std::runtime_error("GeneralTHDMMatching::CMBMll(): order " + out.str() + "not implemeted");
1389 }
1390 vmcBMll.push_back(mcBMll);
1391 return (vmcBMll);
1392}
FULLNNLO
@ FULLNNLO
Definition: OrderScheme.h:39
GeneralTHDMMatching::CphiU
virtual gslpp::complex CphiU(double xHp, double xt, double vev, double xphi, double mu, double Ri1, double Ri2, double Ri3, double mHi_2, double lambda3, double Relambda7, double Imlambda7, gslpp::complex su, gslpp::complex sd)
Definition: GeneralTHDMMatching.cpp:1260
GeneralTHDMMatching::C10Bll
virtual double C10Bll(double xt, double xHp, gslpp::complex su)
Definition: GeneralTHDMMatching.cpp:1067
GeneralTHDMMatching::CPboxBll
virtual gslpp::complex CPboxBll(double xt, double xHp, gslpp::complex su, gslpp::complex sd, gslpp::complex sl)
Definition: GeneralTHDMMatching.cpp:1083
GeneralTHDMMatching::CSboxBll
virtual gslpp::complex CSboxBll(double xt, double xHp, gslpp::complex su, gslpp::complex sd, gslpp::complex sl)
Definition: GeneralTHDMMatching.cpp:1073
GeneralTHDMMatching::CPZUBll
virtual gslpp::complex CPZUBll(double xt, double xHp, double sW2, gslpp::complex su, gslpp::complex sd)
Definition: GeneralTHDMMatching.cpp:1093
ModelMatching::CMBMll
virtual std::vector< WilsonCoefficient > & CMBMll(QCD::lepton lepton)=0
QCD::BOTTOM
@ BOTTOM
Definition: QCD.h:329
QCD::TOP
@ TOP
Definition: QCD.h:328
WilsonCoefficient::setCoeff
void setCoeff(const gslpp::vector< gslpp::complex > &z, orders order_i)
Definition: WilsonCoefficient.h:34
WilsonTemplate::getScheme
schemes getScheme() const
Definition: WilsonTemplate.h:98
WilsonTemplate::setMu
virtual void setMu(double mu)
Definition: WilsonTemplate.h:92
WilsonTemplate::getOrder
orders getOrder() const
Definition: WilsonTemplate.h:65
lambda3
An observable class for the quartic Higgs potential coupling .
Definition: THDMquantities.h:428
mH1_2
tanbeta constructor.
Definition: GeneralTHDMquantities.h:46
mH2_2
The mass of the SM Higgs.
Definition: GeneralTHDMquantities.h:68
mH3_2
The mass of the SM Higgs.
Definition: GeneralTHDMquantities.h:91

◆ CMbsg()

std::vector< WilsonCoefficient > & GeneralTHDMMatching::CMbsg ( )
virtual

operator basis: current current; qcd penguins; magnetic and chromomagnetic penguins; semileptonic

Returns
GeneralTHDM Wilson coefficients, Misiak basis, for \( B \rightarrow X_{s} \gamma, l^{+} l^{-} \)

Implements ModelMatching.

Definition at line 1468 of file GeneralTHDMMatching.cpp.

1468 {
1469 vmcbsg = StandardModelMatching::CMbsg();
1470
1471 if (!myGTHDM.getATHDMflag()) {
1472 throw std::runtime_error("bsgamma is only available in the ATHDM at the moment.");
1473 return (vmcbsg);
1474 }
1475
1476 double Muw = myGTHDM.getMuw();
1477
1478 double Mut = myGTHDM.getMut();
1479 double mt = myGTHDM.Mrun(Mut, myGTHDM.getQuarks(QCD::TOP).getMass_scale(),
1480 myGTHDM.getQuarks(QCD::TOP).getMass(), QCD::TOP, FULLNNLO);
1481 double mHp = myGTHDM.getmHp();
1482
1483 gslpp::complex sigmau = myGTHDM.getNu_11();
1484 gslpp::complex sigmad = myGTHDM.getNd_11();
1485
1486 gslpp::complex co = 1.; // (- 4. * GF / sqrt(2)) * SM.computelamt_s(); THIS SHOULD ALREADY BE IMPLEMENTED IN THE OBSERVABLE
1487 mcbsg.setMu(Muw);
1488
1489 switch (mcbsg.getOrder()) {
1490 case NNLO:
1491 for (int j = 0; j < 8; j++) {
1492 mcbsg.setCoeff(j, co * myGTHDM.Alstilde5(Muw) * myGTHDM.Alstilde5(Muw) * setWCbsg(j, sigmau, sigmad, mt, mHp, Muw, NNLO), NNLO);
1493 }
1494 case NLO:
1495 for (int j = 0; j < 8; j++) {
1496 mcbsg.setCoeff(j, co * myGTHDM.Alstilde5(Muw) * setWCbsg(j, sigmau, sigmad, mt, mHp, Muw, NLO), NLO);
1497 }
1498 case LO:
1499 for (int j = 0; j < 8; j++) {
1500 mcbsg.setCoeff(j, co * setWCbsg(j, sigmau, sigmad, mt, mHp, Muw, LO), LO);
1501 }
1502 break;
1503 default:
1504 std::stringstream out;
1505 out << mcbsg.getOrder();
1506 throw std::runtime_error("THDMMatching::CMbsg(): order " + out.str() + "not implemented");
1507 }
1508
1509 vmcbsg.push_back(mcbsg);
1510 return (vmcbsg);
1511}
GeneralTHDMMatching::setWCbsg
gslpp::complex setWCbsg(int i, gslpp::complex sigmau, gslpp::complex sigmad, double mt, double mhp, double mu, orders order)
Definition: GeneralTHDMMatching.cpp:1517
ModelMatching::CMbsg
virtual std::vector< WilsonCoefficient > & CMbsg()=0

◆ CMdbs2()

std::vector< WilsonCoefficient > & GeneralTHDMMatching::CMdbs2 ( )
virtual

Calculates the function of Eq. (68) of 1607.06292.

Calculates the function of Eq. (68) of 1607.06292

Returns
function of Eq. (68) of 1607.06292
GeneralTHDM Wilson coefficients for \( B_s \to \bar{B_s}\) according to [Geng:1988bq], [Deschamps:2009rh]

Reimplemented from StandardModelMatching.

Definition at line 967 of file GeneralTHDMMatching.cpp.

967 {
968
969 double Mut = myGTHDM.getMut();
970 double xt = x_t(Mut);
971 double GF = myGTHDM.getGF();
972 double MW = myGTHDM.Mw();
973 gslpp::complex co = GF / 4. / M_PI * MW * myGTHDM.getCKM().computelamt_s();
974 //double tanb = myGTHDM.gettanb();
975 double mHp2 = myGTHDM.getmHp2();
976 double xHW = mHp2 / (MW * MW);
977 double xtH = xt / xHW;
978 double mb = myGTHDM.getQuarks(QCD::BOTTOM).getMass();
979 double sd = (myGTHDM.getNd_11()).real();
980 double su = (myGTHDM.getNu_11()).real();
981 double SWH = xtH * ((2.0 * xHW - 8.0) * log(xtH) / ((1.0 - xHW)*(1.0 - xtH)*(1.0 - xtH)) + 6.0 * xHW * log(xt) / ((1.0 - xHW)*(1.0 - xt)*(1.0 - xt))-(8.0 - 2.0 * xt) / ((1.0 - xt)*(1.0 - xtH))) * su*su; //su*su = sigu.abs2()
982 double SHH = xtH * ((1.0 + xtH) / ((1.0 - xtH)*(1.0 - xtH)) + 2.0 * xtH * log(xtH) / ((1.0 - xtH)*(1.0 - xtH)*(1.0 - xtH))) * su * su * su*su; //su*su*su*su = sigu.abs2()*sigu.abs2()
983 double C1bsSRR = 4.0 * mb * mb * xt * xtH * sd * su / (mHp2 * (xtH - 1.0)*(xtH - 1.0)*(xtH - 1.0)) //sd*su = sigd*sigu.conjugate()
984 * (sd * su * (2.0 * (xtH - 1.0)-(xtH + 1.0) * log(xtH)) //sd*su = sigd*sigu.conjugate()
985 + (2.0 * xt * xt * (xtH - 1.0)*(xtH - 1.0)*(xtH - 1.0) * log(xt) / (xt - 1.0)
986 + 2.0 * xt * (xtH - 1.0)*((xt - xtH)*(xtH - 1.0)+(xtH - xt * xtH) * log(xtH))) / ((xt - 1.0)*(xt - xtH) * xtH));
987
988 vmcds = StandardModelMatching::CMdbs2();
989 mcdbs2.setMu(Mut);
990
991 switch (mcdbs2.getOrder()) {
992 case NNLO:
993 case NLO:
994 case LO:
995 mcdbs2.setCoeff(0, co * co * xt * (SWH + SHH), LO);
996 break;
997 default:
998 std::stringstream out;
999 out << mcdbs2.getOrder();
1000 throw std::runtime_error("GeneralTHDMMatching::CMdbs2(): order " + out.str() + "not implemented");
1001 }
1002
1003 vmcds.push_back(mcdbs2);
1004 //The following are the primed coefficients.
1005 mcdbs2.setMu(Mut);
1006
1007 switch (mcdbs2.getOrder()) {
1008 case NNLO:
1009 case NLO:
1010 case LO:
1011 mcdbs2.setCoeff(1, co * co * C1bsSRR, LO);
1012 break;
1013 default:
1014 std::stringstream out;
1015 out << mcdbs2.getOrder();
1016 throw std::runtime_error("GeneralTHDMMatching::CMdbs2(): order " + out.str() + "not implemented");
1017 }
1018
1019 vmcds.push_back(mcdbs2);
1020
1021 return (vmcds);
1022}
StandardModelMatching::CMdbs2
virtual std::vector< WilsonCoefficient > & CMdbs2()
,
Definition: StandardModelMatching.cpp:1016

◆ CMdiujleptonknu()

std::vector< WilsonCoefficient > & GeneralTHDMMatching::CMdiujleptonknu ( int  i,
int  j,
int  k 
)
virtual
Returns
GeneralTHDM Wilson coefficients for \( \bar{d}_i u_j \bar{\nu} \ell_k \) operators in the JMS basis ordered as CnueduVLLkkij, CnueduVLRkkij, CnueduSRRkkij, CnueduSRLkkij, CnueduTRRkkij

Definition at line 1398 of file GeneralTHDMMatching.cpp.

1398 {
1399
1400
1401 if (!myGTHDM.getATHDMflag()) {
1402 throw std::runtime_error("CMdiujleptonknu is only available in the ATHDM at the moment.");
1403 return (vmcbsg);
1404 }
1405
1406
1407 double Muw = myGTHDM.getMuw();
1408 double GF = myGTHDM.getGF();
1409 myCKM = myGTHDM.getVCKM();
1410 QCD::meson mymeson;
1411 double mu;
1412 double md;
1413 if (i == 0 && j == 0) {
1414 mymeson = QCD::P_P;
1415 mu = myGTHDM.getQuarks(QCD::UP).getMass();
1416 md = myGTHDM.getQuarks(QCD::DOWN).getMass();
1417 } else if (i == 1 && j == 0) {
1418 mymeson = QCD::K_P;
1419 mu = myGTHDM.getQuarks(QCD::UP).getMass();
1420 md = myGTHDM.getQuarks(QCD::STRANGE).getMass();
1421 } else if (i == 2 && j == 0) {
1422 mymeson = QCD::B_P;
1423 mu = myGTHDM.getQuarks(QCD::UP).getMass();
1424 md = myGTHDM.getQuarks(QCD::BOTTOM).getMass();
1425 } else if (i == 2 && j == 1) {
1426 mymeson = QCD::B_C;
1427 mu = myGTHDM.getQuarks(QCD::CHARM).getMass();
1428 md = myGTHDM.getQuarks(QCD::BOTTOM).getMass();
1429 } else if (i == 0 && j == 1) {
1430 mymeson = QCD::D_P;
1431 mu = myGTHDM.getQuarks(QCD::CHARM).getMass();
1432 md = myGTHDM.getQuarks(QCD::DOWN).getMass();
1433 } else if (i == 1 && j == 1) {
1434 mymeson = QCD::D_S;
1435 mu = myGTHDM.getQuarks(QCD::CHARM).getMass();
1436 md = myGTHDM.getQuarks(QCD::STRANGE).getMass();
1437 } else
1438 throw std::runtime_error("GeneralTHDMMatching::CMuleptonnu(): doesn't include that meson");
1439
1440 double mP = myGTHDM.getMesons(mymeson).getMass();
1441 gslpp::complex zetad = myGTHDM.getNd_11();
1442 gslpp::complex zetal = myGTHDM.getNl_11();
1443 gslpp::complex zetau = myGTHDM.getNu_11();
1444 double mHp2 = myGTHDM.getmHp2();
1445
1446
1447 vmculeptonnu = StandardModelMatching::CMdiujleptonknu(i, j, k);
1448
1449 mculeptonnu.setMu(Muw);
1450
1451 switch (mculeptonnu.getOrder()) {
1452 case NNLO:
1453 case NLO:
1454 case LO:
1455 mculeptonnu.setCoeff(0, -4. * GF * myCKM(j, i).conjugate() / sqrt(2.)*(mP * mP * zetal.conjugate()*(md * zetad + mu * zetau) / mHp2 / (md + mu)), LO);
1456 break;
1457 default:
1458 std::stringstream out;
1459 out << mculeptonnu.getOrder();
1460 throw std::runtime_error("GeneralTHDMMatching::CMdiujleptonknu(): order " + out.str() + "not implemented");
1461 }
1462 vmculeptonnu.push_back(mculeptonnu);
1463 return (vmculeptonnu);
1464
1465
1466}
QCD::meson
meson
An enum type for mesons.
Definition: QCD.h:336
QCD::K_P
@ K_P
Definition: QCD.h:340
QCD::D_P
@ D_P
Definition: QCD.h:342
QCD::B_C
@ B_C
Definition: QCD.h:347
QCD::B_P
@ B_P
Definition: QCD.h:345
QCD::D_S
@ D_S
Definition: QCD.h:343
QCD::P_P
@ P_P
Definition: QCD.h:338
QCD::UP
@ UP
Definition: QCD.h:324
QCD::DOWN
@ DOWN
Definition: QCD.h:325
QCD::STRANGE
@ STRANGE
Definition: QCD.h:327
QCD::CHARM
@ CHARM
Definition: QCD.h:326

◆ CMgminus2mu()

std::vector< WilsonCoefficient > & GeneralTHDMMatching::CMgminus2mu ( )
virtual

Wilson coefficient for \( (g-2)_{\mu} \).

Returns

Definition at line 917 of file GeneralTHDMMatching.cpp.

917 {
918
919 /* gslpp::complex su = myGTHDM.getNu_11();
920 gslpp::complex sd = myGTHDM.getNd_11();
921 gslpp::complex sl = myGTHDM.getNl_11();*/
922
923 /* double Imlambda5=myGTHDM.getImlambda5();
924 double Imlambda6=myGTHDM.getMyGTHDMCache()->Imlambda6;
925 double Imlambda7=myGTHDM.getMyGTHDMCache()->Imlambda7;
926 double sinalpha2=myGTHDM.getsinalpha2();
927 double sinalpha3=myGTHDM.getsinalpha3();*/
928
929 vmcgminus2mu = StandardModelMatching::CMgminus2mu();
930
931 double gminus2muLOvalue = gminus2muLO();
932 double gminus2muNLOvalue = gminus2muNLO();
933 //gminus2muNLOvalue=0.; //For the moment we set the NLO contribution to zero
934
935// std::cout<<"\033[1;33m gminus2muLOvalue = \033[0m "<< gminus2muLOvalue << std::endl;
936// std::cout<<"\033[1;33m gminus2muNLOvalue = \033[0m "<< gminus2muNLOvalue << std::endl;
937
938
939
940
941 switch (mcgminus2mu.getOrder()) {
942 case LO:
943 mcgminus2mu.setCoeff(0, gminus2muLOvalue, LO); //g-2_muR
944 mcgminus2mu.setCoeff(1, 0., LO); //g-2_muL
945 break;
946 case NLO:
947 mcgminus2mu.setCoeff(0, gminus2muLOvalue + gminus2muNLOvalue, NLO); //g-2_muR
948 mcgminus2mu.setCoeff(1, 0., NLO); //g-2_muL
949 break;
950 case NNLO:
951 default:
952 std::stringstream out;
953 out << mcgminus2mu.getOrder();
954 throw std::runtime_error("GeneralTHDMMatching::CMgminus2mu(): order " + out.str() + " not implemented.\nOnly leading order (LO) or next-to-leading order (NLO) are allowed.");
955 }
956
957 vmcgminus2mu.push_back(mcgminus2mu);
958 return (vmcgminus2mu);
959
960}
GeneralTHDMMatching::gminus2muLO
virtual double gminus2muLO()
Calculates amplitudes for at one loop from 1502.04199, before was used.
Definition: GeneralTHDMMatching.cpp:129
GeneralTHDMMatching::gminus2muNLO
virtual double gminus2muNLO()
Calculates amplitudes for at approximate two-loop from .
Definition: GeneralTHDMMatching.cpp:525

◆ CPboxBll()

gslpp::complex GeneralTHDMMatching::CPboxBll ( double  xt,
double  xHp,
gslpp::complex  su,
gslpp::complex  sd,
gslpp::complex  sl 
)
virtual
Returns
Box CP Wilson coefficient of \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1083 of file GeneralTHDMMatching.cpp.

1083 {
1084
1085
1086 gslpp::complex CPboxU = -xt / (8 * (xHp - xt))*(-sl * su.conjugate()*(xt / (xt - 1) * log(xt) - xHp / (xHp - 1) * log(xHp))
1087 + su * sl.conjugate()*(1 - (xHp - xt * xt) / ((xHp - xt)*(xt - 1)) * log(xt)
1088 - (xHp * (xt - 1)) / ((xHp - xt)*(xHp - 1)) * log(xHp)) + 2 * sd * sl.conjugate()*(log(xt) - log(xHp)));
1089
1090 return CPboxU;
1091}

◆ CphiU()

gslpp::complex GeneralTHDMMatching::CphiU ( double  xHp,
double  xt,
double  vev,
double  xphi,
double  mu,
double  Ri1,
double  Ri2,
double  Ri3,
double  mHi_2,
double  lambda3,
double  Relambda7,
double  Imlambda7,
gslpp::complex  su,
gslpp::complex  sd 
)
virtual

Definition at line 1260 of file GeneralTHDMMatching.cpp.

1260 {
1261 gslpp::complex i = gslpp::complex::i();
1262 gslpp::complex CphiU = xt * ((1 / (2 * xphi))*(su - sd)*(1 + su.conjugate() * sd)*
1263 (Ri2 + i * Ri3) * mu + (vev * vev / mHi_2) * lambdaHHphi(lambda3, Relambda7, Imlambda7, Ri1, Ri2, Ri3) * g0(xHp, xt, su, sd) + Ri1 * ((1 / (2.0 * xphi)) * g1a(xHp, xt, su, sd)) + Ri2 * (1 / (2.0 * xphi) * g2a(xHp, xt, su, sd)) + i * Ri3 * (1 / (2.0 * xphi)) * g3a(xHp, xt, su, sd));
1264
1265 return CphiU;
1266}
GeneralTHDMMatching::g3a
virtual gslpp::complex g3a(double xHp, double xt, gslpp::complex su, gslpp::complex sd)
Definition: GeneralTHDMMatching.cpp:1246
GeneralTHDMMatching::lambdaHHphi
virtual gslpp::complex lambdaHHphi(double lambda3, double Relambda7, double Imlambda7, double Ri1, double Ri2, double Ri3)
Definition: GeneralTHDMMatching.cpp:1255
GeneralTHDMMatching::g1a
virtual gslpp::complex g1a(double xHp, double xt, gslpp::complex su, gslpp::complex sd)
Definition: GeneralTHDMMatching.cpp:1228
GeneralTHDMMatching::g0
virtual gslpp::complex g0(double xHp, double xt, gslpp::complex su, gslpp::complex sd)
Definition: GeneralTHDMMatching.cpp:1220
GeneralTHDMMatching::g2a
virtual gslpp::complex g2a(double xHp, double xt, gslpp::complex su, gslpp::complex sd)
Definition: GeneralTHDMMatching.cpp:1237

◆ CPZUBll()

gslpp::complex GeneralTHDMMatching::CPZUBll ( double  xt,
double  xHp,
double  sW2,
gslpp::complex  su,
gslpp::complex  sd 
)
virtual
Returns
Z-penguin in the unitary gauge CP Wilson coefficient of \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1093 of file GeneralTHDMMatching.cpp.

1093 {
1094
1095
1096 // CPZF. Z-penguins diagrams. Eq. (52)
1097
1098 gslpp::complex CPZF = (xt / (4 * (xHp - xt)*(xHp - xt)))*(sd * su.conjugate()*(-((xt + xHp) / 2)
1099 + ((xt * xHp) / (xHp - xt))*(log(xHp) - log(xt))) +
1100 su.abs2()*(1 / (6 * (xHp - xt)))*((xHp * xHp - 8 * xHp * xt - 17 * xt * xt) / 6
1101 + ((xt * xt * (3 * xHp + xt)) / (xHp - xt))*(log(xHp) - log(xt))))
1102 + sW2 * (xt / (6 * (xHp - xt)*(xHp - xt)))*(sd * su.conjugate()*((5 * xt - 3 * xHp) / 2
1103 + ((xHp * (2 * xHp - 3 * xt)) / (xHp - xt))*(log(xHp) - log(xt))) -
1104 su.abs2()*(1 / (6 * (xHp - xt)))*(((4 * xHp * xHp * xHp - 12 * xHp * xHp * xHp * xt + 9 * xHp * xt * xt + 3 * xt * xt * xt) / (xHp - xt))* (log(xHp) - log(xt))
1105 - (17 * xHp * xHp - 64 * xHp * xt + 71.0 * xt * xt) / 6));
1106
1107 //CPGBF. Goldstone penguin diagrams. Eq. (53)
1108
1109 gslpp::complex CPGBF = su.abs2()*(1 - sW2)*(xt * xt / (4 * (xHp - xt)*(xHp - xt)))*
1110 (xHp * (log(xHp) - log(xt)) + xt - xHp);
1111
1112 //CPZU. Z-penguin diagrams in the unitary gauge. Eq. (54)
1113
1114 gslpp::complex CPZU = CPZF + CPGBF;
1115
1116 return CPZU;
1117
1118}

◆ CSboxBll()

gslpp::complex GeneralTHDMMatching::CSboxBll ( double  xt,
double  xHp,
gslpp::complex  su,
gslpp::complex  sd,
gslpp::complex  sl 
)
virtual
Returns
Box CS Wilson coefficient of \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1073 of file GeneralTHDMMatching.cpp.

1073 {
1074
1075 gslpp::complex CSboxU = xt / (8 * (xHp - xt))*(sl * su.conjugate()*(xt / (xt - 1) * log(xt)
1076 - xHp / (xHp - 1) * log(xHp))
1077 + su * sl.conjugate()*(1 - (xHp - xt * xt) / ((xHp - xt)*(xt - 1)) * log(xt)
1078 - (xHp * (xt - 1)) / ((xHp - xt)*(xHp - 1)) * log(xHp)) + 2 * sd * sl.conjugate()*(log(xt) - log(xHp)));
1079
1080 return (CSboxU);
1081}

◆ f1()

double GeneralTHDMMatching::f1 ( double  xHp,
double  xt 
)
virtual
Returns
f1 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1120 of file GeneralTHDMMatching.cpp.

1120 {
1121
1122
1123 double f1 = (1.0 / (2.0 * (xHp - xt)))*(-xHp + xt + xHp * log(xHp) -
1124 xt * log(xt));
1125 return f1;
1126
1127}
GeneralTHDMMatching::f1
virtual double f1(double xHp, double xt)
Definition: GeneralTHDMMatching.cpp:1120

◆ f10()

double GeneralTHDMMatching::f10 ( double  xHp,
double  xt 
)
virtual
Returns
f10 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1211 of file GeneralTHDMMatching.cpp.

1211 {
1212 double f10 = (1.0 / (8 * (xHp - xt)))*
1213 ((xHp - xt) / ((xHp - 1.0)*(xt - 1.0)) + ((xt * (xt - 2.0)) / (xt - 1.0)*(xt - 1.0)) *
1214 log(xt) - ((xHp * (xHp - 2.0)) / (xHp - 1.0)*(xHp - 1.0)) * log(xHp));
1215
1216 return f10;
1217
1218}
GeneralTHDMMatching::f10
virtual double f10(double xHp, double xt)
Definition: GeneralTHDMMatching.cpp:1211

◆ F1oneloopgm2()

double GeneralTHDMMatching::F1oneloopgm2 ( const double  ratio_sq)

Loop \( (g-2)_{\mu} \) at one loop from 1502.04199. The complete loop function is \int_0^1 dx \frac{x^2(2-x)}{ratio_sq*x^2-x+1} but for small values of ratio_sq we can expand and take only the first term for which the functions becomes (-7.0/6.0-log(ratio_sq))

One-loop function for the contribution to the g-2

Returns

Definition at line 37 of file GeneralTHDMMatching.cpp.

37 {
38
39
40 double final_result;
41
42
43 if(ratio_sq>0.25){
44 final_result = (ratio_sq*(-2 + 3*ratio_sq) - 2*(-1 + ratio_sq)*sqrt(-1 + 4*ratio_sq)*
45 (M_PI/2 - atan(sqrt(-1 + 4*ratio_sq))) - 2*(-1 + ratio_sq)*sqrt(-1 + 4*ratio_sq)*
46 atan((-1 + 2*ratio_sq)/sqrt(-1 + 4*ratio_sq)) + (-1 + 3*ratio_sq)*log(ratio_sq))
47 /(2.*ratio_sq*ratio_sq*ratio_sq);
48 }
49 else if(ratio_sq==0.25){
50 final_result=-10 + 4*log(16);
51 }
52 else if(0.001<ratio_sq and ratio_sq<0.25){
53
54 double arccoth_ratio = 0.5 * log(((-1 + 2*ratio_sq)/sqrt(1 - 4*ratio_sq) + 1.0) / ((-1 + 2*ratio_sq)/sqrt(1 - 4*ratio_sq) - 1.0));
55 //ArcCoth((-1 + 2*ratio_sq)/sqrt(1 - 4*ratio_sq))
56 final_result = (2*arccoth_ratio + 2*sqrt(1 -
57 4*ratio_sq)*(ratio_sq*(-2 + 3*ratio_sq) - log(ratio_sq)) +
58 ratio_sq*(10*atanh((1 + sqrt(1 - 4*ratio_sq) - pow(ratio_sq,2))/(2 +
59 (-2 + ratio_sq)*ratio_sq)) + log((1 + sqrt(1 - 4*ratio_sq) -
60 2*ratio_sq)/8.) + 6*sqrt(1 - 4*ratio_sq)*log(ratio_sq) +
61 ratio_sq*log((16*pow(ratio_sq,8))/pow(1 + sqrt(1 - 4*ratio_sq) - 2*(2
62 + sqrt(1 - 4*ratio_sq) - ratio_sq)*ratio_sq,4)) + 2*log(1 + sqrt(1 -
63 4*ratio_sq) + 2*ratio_sq*(-2 - sqrt(1 - 4*ratio_sq) +
64 ratio_sq))))/(4.*sqrt(1 - 4*ratio_sq)*pow(ratio_sq,3));
65 }
66 else{
67 final_result = (-7.0 / 6.0 - log(ratio_sq)) + 0.25*ratio_sq*(-13 - 12*log(ratio_sq));
68 }
69
70 return final_result;
71}

◆ F1tildetwoloopgm2()

double GeneralTHDMMatching::F1tildetwoloopgm2 ( const double  ratio_sq)

Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{ratio_sq}{2}*\int_0^1 dx \frac{1}{ratio_sq-x(1-x)}\log(\frac{ratio_sq}{x(1-x)})

Two-loop (Barr-Zee) function for the contribution to the g-2

Returns

Definition at line 402 of file GeneralTHDMMatching.cpp.

402 {
403
404 gslpp::complex rsqc = ratio_sq; //We need the functions to be evaluated in the complex space so we define the ratio to be complex
405 double final_result;
406
407 if (ratio_sq == 0.25) {
408 final_result = log(2);
409 } else {
410 //Let's expand when the function is small to avoid numerical instabilities with the complete function
411 if (ratio_sq < 0.001) {
412 final_result = (ratio_sq * (pow(M_PI, 2) + 3 * pow(log(ratio_sq), 2))) / 6. + (pow(ratio_sq, 2)*(-6 + pow(M_PI, 2)
413 + 6 * log(ratio_sq) + 3 * pow(log(ratio_sq), 2))) / 3. + (pow(ratio_sq, 3)*(-11 + 2 * pow(M_PI, 2)
414 + 14 * log(ratio_sq) + 6 * pow(log(ratio_sq), 2))) / 2.;
415 } else {
416 final_result = ((rsqc * (-PolyLog.Li2(-0.5 * (1 + sqrt(1 - 4 * rsqc) - 2 * rsqc) / rsqc) +
417 PolyLog.Li2((-1 + sqrt(1 - 4 * rsqc) + 2 * rsqc) / (2. * rsqc)))) / sqrt(1 - 4 * rsqc)).real();
418 }
419
420 }
421
422 //std::cout<<"\033[1;31m final_result = \033[0m "<< final_result <<std::endl;
423 return (final_result);
424}
GeneralTHDMMatching::PolyLog
const Polylogarithms PolyLog
Definition: GeneralTHDMMatching.h:377
Polylogarithms::Li2
gslpp::complex Li2(const double x) const
The dilogarithm with a real argument, .
Definition: Polylogarithms.cpp:22

◆ F1twoloopgm2()

double GeneralTHDMMatching::F1twoloopgm2 ( const double  ratio_sq)

Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{ratio_sq}{2}*\int_0^1 dx \frac{2x(1-x)-1}{ratio_sq-x(1-x)}\log(\frac{ratio_sq}{x(1-x)})

Two-loop (Barr-Zee) function for the contribution to the g-2

Returns

Definition at line 373 of file GeneralTHDMMatching.cpp.

373 {
374
375 gslpp::complex rsqc = ratio_sq; //We need the functions to be evaluated in the complex space so we define the ratio to be complex
376 double final_result;
377
378 if (ratio_sq == 0.25) {
379 final_result = -2 * ratio_sq;
380 } else {
381 //Let's expand when the function is small to avoid numerical instabilities with the complete function
382 if (ratio_sq < 0.001) {
383 final_result = -(1 / 6.)*(ratio_sq * (12 + M_PI * M_PI - 12 * ratio_sq - 9 * ratio_sq * ratio_sq +
384 2 * M_PI * M_PI * ratio_sq * ratio_sq + 6 * (1 + 2 * ratio_sq + 3 * ratio_sq * ratio_sq) * log(ratio_sq) +
385 (3 + 6 * ratio_sq * ratio_sq) * log(ratio_sq) * log(ratio_sq)));
386 } else {
387 final_result = ((rsqc * ((pow(M_PI, 2) - 12 * sqrt(1 - 4 * rsqc) - 6 * sqrt(1 - 4 * rsqc) * log(rsqc) +
388 3 * pow(log((-2 * rsqc) / (-1 + sqrt(1 - 4 * rsqc) + 2 * rsqc)), 2) -
389 2 * rsqc * (pow(M_PI, 2) + 6 * pow(log(1 - sqrt(1 - 4 * rsqc)), 2) -
390 6 * pow(log(1 + sqrt(1 - 4 * rsqc)), 2) +
391 2 * arctanh(sqrt(1 - 4 * rsqc)) * log(4096 * pow(rsqc, 6)) +
392 3 * pow(log((-2 * rsqc) / (-1 + sqrt(1 - 4 * rsqc) + 2 * rsqc)), 2))) / 6. +
393 (2 - 4 * rsqc) * PolyLog.Li2(-0.5 * (1 + sqrt(1 - 4 * rsqc) - 2 * rsqc) / rsqc))) / sqrt(1 - 4 * rsqc)).real();
394 }
395
396 }
397
398 //std::cout<<"\033[1;31m final_result = \033[0m "<< final_result <<std::endl;
399 return (final_result);
400}

◆ f2()

double GeneralTHDMMatching::f2 ( double  xHp,
double  xt 
)
virtual
Returns
f2 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1129 of file GeneralTHDMMatching.cpp.

1129 {
1130
1131 double f2 = (1.0 / (2.0 * (xHp - xt)))*(xt - ((xHp * xt) / (xHp - xt))*(log(xHp) - log(xt)));
1132
1133 return f2;
1134
1135}
GeneralTHDMMatching::f2
virtual double f2(double xHp, double xt)
Definition: GeneralTHDMMatching.cpp:1129

◆ F2oneloopgm2()

double GeneralTHDMMatching::F2oneloopgm2 ( const double  ratio_sq)

Loop \( (g-2)_{\mu} \) at one loop from 1502.04199. The complete loop function is \int_0^1 dx \frac{-x^3}{ratio_sq*x^2-x+1} but for small values of ratio_sq we can expand and take only the first term for which the functions becomes (11.0/6.0+log(ratio_sq))

One-loop function for the contribution to the g-2

Returns

Definition at line 73 of file GeneralTHDMMatching.cpp.

73 {
74
75 double final_result;
76
77
78 if(ratio_sq>0.25){
79 final_result = ((-2 + 6*ratio_sq)*(M_PI/2 - atan(sqrt(-1 + 4*ratio_sq))) + (-2 +
80 6*ratio_sq)*atan((-1 + 2*ratio_sq)/sqrt(-1 + 4*ratio_sq)) + sqrt(-1
81 + 4*ratio_sq)*(-(ratio_sq*(2 + ratio_sq)) + (-1 +
82 ratio_sq)*log(ratio_sq)))/(2.*pow(ratio_sq,3)*sqrt(-1 + 4*ratio_sq));
83 }
84 else if(ratio_sq==0.25){
85 final_result= -34 + 24*log(4);
86 }
87 else if(0.001<ratio_sq and ratio_sq<0.25){
88
89 final_result = -0.5*(ratio_sq*(2 + ratio_sq) - (-1 + ratio_sq)*log(ratio_sq) +
90 (3*ratio_sq*log((2*ratio_sq)/(1 + sqrt(1 - 4*ratio_sq) - 2*ratio_sq))
91 + log((-2*ratio_sq)/(-1 + sqrt(1 - 4*ratio_sq) + 2*ratio_sq)))/sqrt(1
92 - 4*ratio_sq))/pow(ratio_sq,3);
93 }
94 else{
95 final_result = ((11 + 6*log(ratio_sq))/6. + (ratio_sq*(89 + 60*log(ratio_sq)))/12.);
96 }
97
98 return final_result;
99}

◆ F2twoloopgm2()

double GeneralTHDMMatching::F2twoloopgm2 ( const double  ratio_sq)

Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{1}{2}*\int_0^1 dx \frac{x(x-1)}{ratio_sq-x(1-x)}\log(\frac{ratio_sq}{x(1-x)})

Two-loop (Barr-Zee) function for the contribution to the g-2

Returns

Definition at line 426 of file GeneralTHDMMatching.cpp.

426 {
427
428 gslpp::complex rsqc = ratio_sq; //We need the functions to be evaluated in the complex space so we define the ratio to be complex
429 double final_result;
430
431 if (ratio_sq == 0.25) {
432 final_result = (1 - log(4)) / 4;
433 } else {
434 //Let's expand when the function is small to avoid numerical instabilities with the complete function
435 if (ratio_sq < 0.001) {
436 final_result = (ratio_sq * (2 + log(ratio_sq))) / 2. + (pow(ratio_sq, 2)*(-pow(M_PI, 2) - 3 * pow(log(ratio_sq), 2))) / 6.
437 + (pow(ratio_sq, 3)*(6 - pow(M_PI, 2) - 6 * log(ratio_sq) - 3 * pow(log(ratio_sq), 2))) / 3.;
438 } else {
439 final_result = ((rsqc * (6 + 3 * log(rsqc) + (rsqc * (pow(M_PI, 2) + 6 * arctanh(sqrt(1 - 4 * rsqc)) * log((-2 * rsqc) / (-1 + sqrt(1 - 4 * rsqc)
440 + 2 * rsqc)) + 12 * PolyLog.Li2(-0.5 * (1 + sqrt(1 - 4 * rsqc) - 2 * rsqc) / rsqc))) / sqrt(1 - 4 * rsqc))) / 6.).real();
441 }
442
443 }
444
445 //std::cout<<"\033[1;31m final_result = \033[0m "<< final_result <<std::endl;
446 return (final_result);
447}

◆ f3()

double GeneralTHDMMatching::f3 ( double  xHp,
double  xt 
)
virtual
Returns
f3 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1137 of file GeneralTHDMMatching.cpp.

1137 {
1138
1139 double f3 = (1.0 / (2.0 * (xHp - xt)))*
1140 (xHp - (xHp * xHp * log(xHp)) / (xHp - xt) + (xt * (2 * xHp - xt) * log(xt)) /
1141 (xHp - xt));
1142
1143 return f3;
1144
1145}
GeneralTHDMMatching::f3
virtual double f3(double xHp, double xt)
Definition: GeneralTHDMMatching.cpp:1137

◆ F3oneloopgm2()

double GeneralTHDMMatching::F3oneloopgm2 ( const double  ratio_sq)

Loop \( (g-2)_{\mu} \) at one loop from 1502.04199. The complete loop function is \int_0^1 dx \frac{x^2(1-x)}{ratio_sq*x*(1-x)-x} but for small values of ratio_sq we can expand and take only the first term for which the functions becomes (-1/12. - ratio_sq/60.)

One-loop function for the contribution to the g-2

Returns

Definition at line 101 of file GeneralTHDMMatching.cpp.

101 {
102
103
104
105 double final_result;
106
107
108 if(ratio_sq>1.){
109 throw std::runtime_error("The implemented F3oneloopgm2 only works for values of the ratio"
110 " below one. Add a table interpolation if you need higher values"
111 );
112 }
113 else if(ratio_sq==1.){
114 final_result= -0.5;
115 }
116 else if(0.001<ratio_sq and ratio_sq<1.){
117
118 final_result = ((-2 + ratio_sq)*ratio_sq + 2*(-1 + ratio_sq)*log(1 - ratio_sq))/(2.*ratio_sq*ratio_sq*ratio_sq);
119 }
120 else{
121 final_result = (-(1/6) - ratio_sq/12);
122 }
123
124
125
126 return final_result;
127}

◆ F3twoloopgm2()

double GeneralTHDMMatching::F3twoloopgm2 ( const double  ratio_sq)

Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{1}{2}*\int_0^1 dx \frac{x(3x(4x-1)+10)ratio_sq-x(1-x)}{ratio_sq-x(1-x)}\log(\frac{ratio_sq}{x(1-x)})

Two-loop (Barr-Zee) function for the contribution to the g-2

Returns

Definition at line 449 of file GeneralTHDMMatching.cpp.

449 {
450
451 gslpp::complex rsqc = ratio_sq; //We need the functions to be evaluated in the complex space so we define the ratio to be complex
452 double final_result;
453
454 if (ratio_sq == 0.25) {
455 final_result = 19. / 16;
456 } else {
457 //Let's expand when the function is small to avoid numerical instabilities with the complete function
458 if (ratio_sq < 0.001) {
459 final_result = (ratio_sq * (2 + log(ratio_sq))) / 2. + (pow(ratio_sq, 3)*(-51 + pow(M_PI, 2) + 51 * log(ratio_sq) + 3 * pow(log(ratio_sq), 2))) / 3.
460 + (pow(ratio_sq, 2)*(180 + 17 * pow(M_PI, 2) + 90 * log(ratio_sq) + 51 * pow(log(ratio_sq), 2))) / 12.;
461 } else {
462 final_result = ((rsqc * (12 * sqrt(1 - 4 * rsqc)*(1 + 15 * rsqc) + pow(M_PI, 2) * rsqc * (-17 + 30 * rsqc) + 6 * sqrt(1 - 4 * rsqc) * log(rsqc)
463 + 3 * rsqc * ((17 - 30 * rsqc) * log(1 + sqrt(1 - 4 * rsqc)) * log((1 + sqrt(1 - 4 * rsqc)) / pow(rsqc, 2)) + 30 * sqrt(1 - 4 * rsqc) * log(rsqc)
464 + (-17 + 30 * rsqc)*(8 * arctanh(sqrt(1 - 4 * rsqc)) * log(2) + log(1 - sqrt(1 - 4 * rsqc))*(3 * log(1 - sqrt(1 - 4 * rsqc))
465 - 2 * log((1 + sqrt(1 - 4 * rsqc)) * rsqc)))) + 12 * rsqc * (-17 + 30 * rsqc) * PolyLog.Li2(-0.5 * (1 + sqrt(1 - 4 * rsqc) - 2 * rsqc) / rsqc))) / (12. * sqrt(1 - 4 * rsqc))).real();
466 }
467
468 }
469
470 //std::cout<<"\033[1;31m final_result = \033[0m "<< final_result <<std::endl;
471 return (final_result);
472}

◆ f4()

double GeneralTHDMMatching::f4 ( double  xHp,
double  xt 
)
virtual
Returns
f4 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1147 of file GeneralTHDMMatching.cpp.

1147 {
1148
1149
1150 double f4 = (1.0 / (4 * (xHp - xt)*(xHp - xt)))*((xt * (3.0 * xHp - xt)) / 2.0 -
1151 ((xHp * xHp * xt) / (xHp - xt))*(log(xHp) - log(xt)));
1152 return f4;
1153
1154}
GeneralTHDMMatching::f4
virtual double f4(double xHp, double xt)
Definition: GeneralTHDMMatching.cpp:1147

◆ F4twoloopgm2()

double GeneralTHDMMatching::F4twoloopgm2 ( const double  ratio_sq)

Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. In the notation of [1502.04199] this corresponds to (Q_t*x+Q_b*(1-x))x(1+x)G(ratio_sq,0), i.e. we neglect mb/mHp and mb/mW which is a great approximation The complete loop function is \int_0^1 dx (2/3 x-1/3 *(1-x))x(1+x)\frac{\log{\frac{ratio_sq*x}{x(1-x)}}}{x(1-x-ratio_sq*x)})

Two-loop (Barr-Zee) function for the contribution to the g-2

Returns

◆ f5()

double GeneralTHDMMatching::f5 ( double  xHp,
double  xt 
)
virtual
Returns
f5 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1156 of file GeneralTHDMMatching.cpp.

1156 {
1157
1158
1159 double f5 = (1.0 / (4 * (xHp - xt)*(xHp - xt)))*((xt * (xHp - 3.0 * xt)) / 2.0 -
1160 ((xHp * xt * (xHp - 2.0 * xt)) / (xHp - xt))*(log(xHp) - log(xt)));
1161 return f5;
1162
1163}
GeneralTHDMMatching::f5
virtual double f5(double xHp, double xt)
Definition: GeneralTHDMMatching.cpp:1156

◆ F5twoloopgm2()

double GeneralTHDMMatching::F5twoloopgm2 ( const double  ratio_sq)

Loop \( (g-2)_{\mu} \) at two loops (Barr-Zee) from 1502.04199. In the notation of [1502.04199] this corresponds to (Q_t*x+Q_b*(1-x))x(1-x)G(ratio_sq,0), i.e. we neglect mb/mHp and mb/mW which is a great approximation The complete loop function is \int_0^1 dx (2/3 x-1/3 *(1-x))x(1-x)\frac{\log{\frac{ratio_sq*x}{x(1-x)}}}{x(1-x-ratio_sq*x)})

Two-loop (Barr-Zee) function for the contribution to the g-2

Returns

◆ f6()

double GeneralTHDMMatching::f6 ( double  xHp,
double  xt 
)
virtual
Returns
f6 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1165 of file GeneralTHDMMatching.cpp.

1165 {
1166
1167
1168 double f6 = (1.0 / (2.0 * (xHp - xt)))*
1169 ((xt * (xt * xt - 3.0 * xHp * xt + 9 * xHp - 5 * xt - 2.0)) / (4 * (xt - 1.0)*(xt - 1.0)) +
1170 ((xHp * (xHp * xt - 3.0 * xHp + 2.0 * xt)) / (2.0 * (xHp - 1.0)*(xHp - xt))) *
1171 log(xHp) + ((xHp * xHp * (-2.0 * xt * xt * xt + 6 * xt * xt - 9 * xt + 2.0) +
1172 3.0 * xHp * xt * xt * (xt * xt - 2.0 * xt + 3.0) - xt * xt * (2.0 * xt * xt * xt - 3.0 * xt * xt +
1173 3.0 * xt + 1.0)) / (2.0 * (xt - 1.0)*(xt - 1.0)*(xt - 1.0)*(xHp - xt))) * log(xt));
1174
1175 return f6;
1176
1177}
GeneralTHDMMatching::f6
virtual double f6(double xHp, double xt)
Definition: GeneralTHDMMatching.cpp:1165

◆ f7()

double GeneralTHDMMatching::f7 ( double  xHp,
double  xt 
)
virtual
Returns
f7 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1179 of file GeneralTHDMMatching.cpp.

1179 {
1180
1181
1182 double f7 = (1.0 / (2.0 * (xHp - xt)))*
1183 (((xt * xt + xt - 8)*(xHp - xt)) / (4 * (xt - 1.0)*(xt - 1.0)) -
1184 ((xHp * (xHp + 2.0)) / (2.0 * (xHp - 1.0))) * log(xHp) +
1185 ((xHp * (xt * xt * xt - 3.0 * xt * xt + 3.0 * xt + 2.0) + 3.0 * (xt - 2.0) * xt * xt) /
1186 (2.0 * (xt - 1.0)*(xt - 1.0)*(xt - 1.0))) * log(xt));
1187 return f7;
1188
1189}
GeneralTHDMMatching::f7
virtual double f7(double xHp, double xt)
Definition: GeneralTHDMMatching.cpp:1179

◆ f8()

double GeneralTHDMMatching::f8 ( double  xHp,
double  xt 
)
virtual
Returns
f8 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1191 of file GeneralTHDMMatching.cpp.

1191 {
1192
1193 double f8 = (1.0 / (4 * (xHp - xt)))*((xt * log(xt)) / (xt - 1.0) -
1194 (xHp * log(xHp)) / (xHp - 1.0));
1195
1196 return f8;
1197
1198}
GeneralTHDMMatching::f8
virtual double f8(double xHp, double xt)
Definition: GeneralTHDMMatching.cpp:1191

◆ f9()

double GeneralTHDMMatching::f9 ( double  xHp,
double  xt 
)
virtual
Returns
f9 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1200 of file GeneralTHDMMatching.cpp.

1200 {
1201
1202
1203 double f9 = (1.0 / (8 * (xHp - xt)))*(xHp / (xHp - 1.0) +
1204 (xt * xt * log(xt)) / ((xt - 1.0)*(xHp - xt)) -
1205 ((xHp * (xHp * xt + xHp - 2.0 * xt)) / ((xHp - 1.0)*(xHp - 1.0)*(xHp - xt))) *
1206 log(xHp));
1207 return f9;
1208
1209}
GeneralTHDMMatching::f9
virtual double f9(double xHp, double xt)
Definition: GeneralTHDMMatching.cpp:1200

◆ g0()

gslpp::complex GeneralTHDMMatching::g0 ( double  xHp,
double  xt,
gslpp::complex  su,
gslpp::complex  sd 
)
virtual
Returns
g0 needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1220 of file GeneralTHDMMatching.cpp.

1220 {
1221
1222 gslpp::complex g0 = (1.0 / (4.0 * (xHp - xt)))*((sd * su.conjugate()*(xt / (xHp - xt)*(log(xHp) - log(xt)) - 1.0))
1223 + su.abs2()*(xt * xt / (2.0 * (xHp - xt)*(xHp - xt))*(log(xHp) - log(xt)) + (xHp - 3.0 * xt) / (4.0 * (xHp - xt))));
1224
1225 return g0;
1226}

◆ g1a()

gslpp::complex GeneralTHDMMatching::g1a ( double  xHp,
double  xt,
gslpp::complex  su,
gslpp::complex  sd 
)
virtual
Returns
g1a needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1228 of file GeneralTHDMMatching.cpp.

1228 {
1229
1230 gslpp::complex g1a = -(3.0 / 4.0) + sd * su.conjugate()*(xt / (xHp - xt))*(1.0 - (xHp / (xHp - xt))*(log(xHp) - log(xt)))
1231 + su.abs2()*(xt / (2.0 * (xHp - xt)*(xHp - xt)))*((xHp + xt) / 2.0 - ((xHp * xt) / (xHp - xt))*(log(xHp) - log(xt)));
1232
1233 return g1a;
1234
1235}

◆ g2a()

gslpp::complex GeneralTHDMMatching::g2a ( double  xHp,
double  xt,
gslpp::complex  su,
gslpp::complex  sd 
)
virtual
Returns
g2a needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1237 of file GeneralTHDMMatching.cpp.

1237 {
1238
1239 gslpp::complex g2a = sd * sd * su.conjugate() * f1(xHp, xt) + sd * su.conjugate() * su.conjugate() * f2(xHp, xt) + sd * su.abs2() * f3(xHp, xt)
1240 + su * su.abs2() * f4(xHp, xt) - su.conjugate() * su.abs2() * f5(xHp, xt) + su * f6(xHp, xt) - su.conjugate() * f7(xHp, xt) + sd * f1(xHp, xt);
1241
1242 return g2a;
1243
1244}

◆ g3a()

gslpp::complex GeneralTHDMMatching::g3a ( double  xHp,
double  xt,
gslpp::complex  su,
gslpp::complex  sd 
)
virtual
Returns
g3a needed to calculate \( B_q \to l \bar{l}\) according to [Li:2014fea]

Definition at line 1246 of file GeneralTHDMMatching.cpp.

1246 {
1247
1248 gslpp::complex g3a = sd * sd * su.conjugate() * f1(xHp, xt) - sd * su.conjugate() * su.conjugate() * f2(xHp, xt)
1249 + sd * su.abs2() * f3(xHp, xt) + su * su.abs2() * f4(xHp, xt) + su.conjugate() * su.abs2() * f5(xHp, xt)
1250 + su * f6(xHp, xt) + su.conjugate() * f7(xHp, xt) + sd * f1(xHp, xt);
1251
1252 return g3a;
1253}

◆ getPolyLog()

const Polylogarithms GeneralTHDMMatching::getPolyLog ( ) const
inline

Definition at line 359 of file GeneralTHDMMatching.h.

360 {
361 return PolyLog;
362 }

◆ gminus2muLO()

double GeneralTHDMMatching::gminus2muLO ( )
virtual

Calculates amplitudes for \( (g-2)_{\mu} \) at one loop from 1502.04199, before [Broggio:2014mna] was used.

Calculates the muon g-2 at LO

Returns

Definition at line 129 of file GeneralTHDMMatching.cpp.

129 {
130
131 //We'll include the contributions from the paper [1502.04199] of V. Ilisie
132
133 updateGTHDMParameters();
134 double gminus2muLO;
135
136 //add something to note that this is only valid in the aligned case and in the CP-conserving limit
137
138
139
140 double GF = myGTHDM.getGF();
141 double mMU = myGTHDM.getLeptons(StandardModel::MU).getMass();
142
143
144 //std::cout<<"\033[1;31m F3oneloopgm2(2.) = \033[0m "<< F3oneloopgm2(2.) << std::endl;
145 //std::cout<<"\033[1;31m F3oneloopgm2(0.250000001) = \033[0m "<< F3oneloopgm2(0.250000001) << std::endl;
146 //std::cout<<"\033[1;31m F3oneloopgm2(0.25) = \033[0m "<< F3oneloopgm2(0.25) << std::endl;
147 //std::cout<<"\033[1;31m F3oneloopgm2(0.249999999) = \033[0m "<< F3oneloopgm2(0.249999999) << std::endl;
148 //std::cout<<"\033[1;31m F3oneloopgm2(0.00100000000000001) = \033[0m "<< F3oneloopgm2(0.00100000000000001) << std::endl;
149 //std::cout<<"\033[1;31m F3oneloopgm2(0.001) = \033[0m "<< F3oneloopgm2(0.001) << std::endl;
150
151
152
153 /*Mass of the physical scalars*/
154 double mH1_2 = myGTHDM.getMyGTHDMCache()->mH1sq;
155 double mH2_2 = myGTHDM.getMyGTHDMCache()->mH2sq;
156 double mH3_2 = myGTHDM.getMyGTHDMCache()->mH3sq;
157 double mHp2 = myGTHDM.getMyGTHDMCache()->mHp2;
158
159 //gslpp::complex yl1 = myGTHDM.getyl1();
160 gslpp::complex yl1 = myGTHDM.getMyGTHDMCache()->yl1;
161 gslpp::complex yl2 = myGTHDM.getMyGTHDMCache()->yl2;
162 gslpp::complex yl3 = myGTHDM.getMyGTHDMCache()->yl3;
163 gslpp::complex sl = myGTHDM.getMyGTHDMCache()->sl;
164
165 //double alpha1 = myGTHDM.getalpha1();
167 //OLD part of the code, this doesn't make sense anymore
168 //double eta = 0.0;
169 //if(myGTHDM.getSMHiggs()){
170 // eta = alpha1;
171 //}
172 //else{
173 // eta = alpha1 - pi/2.0;
174 //}
175
176 //This was obviously WRONG, here you're expanding for small values of alpha but in the
177 //non "heavy" case this doesn't make any sense! Let's forget about the expansion and consider
178 //the most general case. We will take the Yukawas from the cache
179 //double yl_h = 1.0 + eta*sl.real();
180 //double yl_H = -sl.real() + eta;
181 //double yl_A = - sl.real();
183
184
185
186 double rmu_hSM, rmu_h, rmu_H, rmu_A, rmu_Hp;
187 double part_hSM, part_h, part_H, part_A, part_Hp;
188
189 if (mH1_2 < 100.0 || mH2_2 < 100.0 || mH3_2 < 100.0 || mHp2 < 100.0) {
190 throw std::runtime_error("The implemented approximation for g-2_\\mu only works for Higgs masses above 10 GeV.");
191 }
192
193 if (!myGTHDM.getATHDMflag()) {
194 throw std::runtime_error("(g-2) is only available in the A2HDM.");
195 }
196
197 if (!myGTHDM.getCPconservationflag()) {
198 //Now the CP-violating case is being included, check again if we can remove this
199 throw std::runtime_error("(g-2) is only available in the CP-conserving limit.");
200 }
201
202
203 rmu_hSM = mMU * mMU / mH1_2;
204 rmu_h = mMU * mMU / mH1_2;
205 rmu_H = mMU * mMU / mH2_2;
206 rmu_A = mMU * mMU / mH3_2;
207 rmu_Hp = mMU * mMU / mHp2;
208
209 part_hSM = rmu_hSM * (F1oneloopgm2(rmu_hSM));
210 part_h = rmu_h * (yl1.real() * yl1.real() * F1oneloopgm2(rmu_h) + yl1.imag() * yl1.imag() * F2oneloopgm2(rmu_h));
211 part_H = rmu_H * (yl2.real() * yl2.real() * F1oneloopgm2(rmu_H) + yl2.imag() * yl2.imag() * F2oneloopgm2(rmu_H));
212 part_A = rmu_A * (yl3.real() * yl3.real() * F1oneloopgm2(rmu_A) + yl3.imag() * yl3.imag() * F2oneloopgm2(rmu_A));
213 part_Hp = rmu_Hp * sl.abs2() * F3oneloopgm2(rmu_Hp);
214
215
216 // We're interested only in the NP contribution so we need to remove the SM part
217 gminus2muLO = GF * mMU * mMU / (4.0 * M_PI * M_PI * sqrt(2.0)) * (-part_hSM + part_h + part_H + part_A + part_Hp);
218
219 return (gminus2muLO);
220}
GeneralTHDMMatching::F3oneloopgm2
double F3oneloopgm2(const double ratio_sq)
Loop at one loop from 1502.04199. The complete loop function is \int_0^1 dx \frac{x^2(1-x)}{ratio_sq...
Definition: GeneralTHDMMatching.cpp:101
GeneralTHDMMatching::F1oneloopgm2
double F1oneloopgm2(const double ratio_sq)
Loop at one loop from 1502.04199. The complete loop function is \int_0^1 dx \frac{x^2(2-x)}{ratio_sq...
Definition: GeneralTHDMMatching.cpp:37
GeneralTHDMMatching::F2oneloopgm2
double F2oneloopgm2(const double ratio_sq)
Loop at one loop from 1502.04199. The complete loop function is \int_0^1 dx \frac{-x^3}{ratio_sq*x^2...
Definition: GeneralTHDMMatching.cpp:73
QCD::MU
@ MU
Definition: QCD.h:314

◆ gminus2muNLO()

double GeneralTHDMMatching::gminus2muNLO ( )
virtual

Calculates amplitudes for \( (g-2)_{\mu} \) at approximate two-loop from [Broggio:2014mna].

Calculates the NLO contribution to the muon g-2

Returns

Definition at line 525 of file GeneralTHDMMatching.cpp.

525 {
526
527 updateGTHDMParameters();
528
529// if (!myGTHDM.getSMHiggs()) {
530// //throw std::runtime_error("The NLO computation of g-2 is only implemented for SM Higgs since the coupling of charged scalars to neutral scalars is included only for that case.");
531// throw std::runtime_error("The NLO computation of g-2 for the heavy higgs scenario must be reviewed although in principle is implemented");
532//
533// }
534
535 //myGTHDM.getMyGTHDMCache()->mH1sq;
536
537
538
539// std::cout<<"\033[1;31m ip_integral_x_1mx2_G_variable_set_0(0.000004) = \033[0m "<< myGTHDM.getMyGTHDMCache()->ip_integral_x_1mx2_G_variable_set_0(0.0000041) << std::endl;
540// std::cout<<"\033[1;31m ip_integral_x_1mx2_G_variable_set_0(1.) = \033[0m "<< myGTHDM.getMyGTHDMCache()->ip_integral_x_1mx2_G_variable_set_0(1.) << std::endl;
541// std::cout<<"\033[1;31m ip_integral_x_1mx2_G_variable_set_0(500.) = \033[0m "<< myGTHDM.getMyGTHDMCache()->ip_integral_x_1mx2_G_variable_set_0(500.) << std::endl;
542// std::cout<<"\033[1;31m ip_integral_x_1mx2_G_variable_set_0(40000) = \033[0m "<< myGTHDM.getMyGTHDMCache()->ip_integral_x_1mx2_G_variable_set_0(40000.) << std::endl;
543
544
545
546
547 //SM constants
548 double Ale = myGTHDM.getAle();
549 double vev = myGTHDM.v();
550 double sW2 = myGTHDM.sW2();
551 double Vtb2 = (myGTHDM.getVCKM()(2, 2)).abs2();
552 double Nc = 3.;
553 double Qu = 2./3.;
554 double Qd = -1./3.;
555// double Ql = -1;
556
557 //SM masses
558 double mMU = myGTHDM.getLeptons(StandardModel::MU).getMass();
559 double mTAU = myGTHDM.getLeptons(StandardModel::TAU).getMass();
560 double mt = myGTHDM.getQuarks(QCD::TOP).getMass();
561 double mb = myGTHDM.getQuarks(QCD::BOTTOM).getMass();
562 double Mw = myGTHDM.Mw();
563
564 //std::cout<<"\033[1;31m Vtb = \033[0m "<< Vtb <<std::endl;
565
566 //GTHDM Yukawa couplings
567 gslpp::complex yu1 = myGTHDM.getMyGTHDMCache()->yu1;
568 double yu1R = yu1.real();
569 double yu1I = yu1.imag();
570 gslpp::complex yu2 = myGTHDM.getMyGTHDMCache()->yu2;
571 double yu2R = yu2.real();
572 double yu2I = yu2.imag();
573 gslpp::complex yu3 = myGTHDM.getMyGTHDMCache()->yu3;
574 double yu3R = yu3.real();
575 double yu3I = yu3.imag();
576
577 gslpp::complex yd1 = myGTHDM.getMyGTHDMCache()->yd1;
578 double yd1R = yd1.real();
579 double yd1I = yd1.imag();
580 gslpp::complex yd2 = myGTHDM.getMyGTHDMCache()->yd2;
581 double yd2R = yd2.real();
582 double yd2I = yd2.imag();
583 gslpp::complex yd3 = myGTHDM.getMyGTHDMCache()->yd3;
584 double yd3R = yd3.real();
585 double yd3I = yd3.imag();
586
587 gslpp::complex yl1 = myGTHDM.getMyGTHDMCache()->yl1;
588 double yl1R = yl1.real();
589 double yl1I = yl1.imag();
590 gslpp::complex yl2 = myGTHDM.getMyGTHDMCache()->yl2;
591 double yl2R = yl2.real();
592 double yl2I = yl2.imag();
593 gslpp::complex yl3 = myGTHDM.getMyGTHDMCache()->yl3;
594 double yl3R = yl3.real();
595 double yl3I = yl3.imag();
596
597
598
599// std::cout<<"\033[1;31m yu1R = \033[0m "<< yu1R <<std::endl;
600// std::cout<<"\033[1;31m yu2R = \033[0m "<< yu2R <<std::endl;
601// std::cout<<"\033[1;31m yu3R = \033[0m "<< yu3R <<std::endl;
602
603// std::cout<<"\033[1;31m yd1R = \033[0m "<< yd1R <<std::endl;
604// std::cout<<"\033[1;31m yd2R = \033[0m "<< yd2R <<std::endl;
605// std::cout<<"\033[1;31m yd3R = \033[0m "<< yd3R <<std::endl;
606
607// std::cout<<"\033[1;31m yl1R = \033[0m "<< yl1R <<std::endl;
608// std::cout<<"\033[1;31m yl2R = \033[0m "<< yl2R <<std::endl;
609// std::cout<<"\033[1;31m yl3R = \033[0m "<< yl3R <<std::endl;
610
611
612// std::cout<<"\033[1;31m yu1I = \033[0m "<< yu1I <<std::endl;
613// std::cout<<"\033[1;31m yu2I = \033[0m "<< yu2I <<std::endl;
614// std::cout<<"\033[1;31m yu3I = \033[0m "<< yu3I <<std::endl;
615
616// std::cout<<"\033[1;31m yd1I = \033[0m "<< yd1I <<std::endl;
617// std::cout<<"\033[1;31m yd2I = \033[0m "<< yd2I <<std::endl;
618// std::cout<<"\033[1;31m yd3I = \033[0m "<< yd3I <<std::endl;
619
620// std::cout<<"\033[1;31m yl1I = \033[0m "<< yl1I <<std::endl;
621// std::cout<<"\033[1;31m yl2I = \033[0m "<< yl2I <<std::endl;
622// std::cout<<"\033[1;31m yl3I = \033[0m "<< yl3I <<std::endl;
623
624 gslpp::complex su = myGTHDM.getMyGTHDMCache()->su;
625 gslpp::complex sd = myGTHDM.getMyGTHDMCache()->sd;
626 gslpp::complex sl = myGTHDM.getMyGTHDMCache()->sl;
627
628
629// std::cout<<"\033[1;31m su = \033[0m "<< su <<std::endl;
630// std::cout<<"\033[1;31m sd = \033[0m "<< sd <<std::endl;
631// std::cout<<"\033[1;31m sl = \033[0m "<< sl <<std::endl;
632
633
634
635 //GTHDM masses
636 double mH1_2 = myGTHDM.getMyGTHDMCache()->mH1sq; //NOTE THAT mH1 MUST BE THE SM HIGGS BOSON MASS
637 double mH2_2 = myGTHDM.getMyGTHDMCache()->mH2sq;
638 double mH3_2 = myGTHDM.getMyGTHDMCache()->mH3sq;
639 double mHp2 = myGTHDM.getMyGTHDMCache()->mHp2;
640
641 //GTHDM Rotation Couplings
642 double R11 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(0,0);
643 double R12 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(0,1);
644 double R13 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(0,2);
645 double R21 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(1,0);
646 double R22 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(1,1);
647 double R23 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(1,2);
648 double R31 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(2,0);
649 double R32 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(2,1);
650 double R33 = myGTHDM.getMyGTHDMCache()->Rij_GTHDM(2,2);
651
652 //GTHDM lambda couplings
653 double lambda3 = myGTHDM.getlambda3();
654 double Relambda7 = myGTHDM.getRelambda7();
655 double Imlambda7 = myGTHDM.getImlambda7();
656
657 //GTHDM coupling of charged scalar with neutral scalars
658 double lambdaHph;
659 double lambdaHpH;
660 double lambdaHpA;
661
662 lambdaHph = R11 * lambda3 + R12 * Relambda7 - R13*Imlambda7;
663 lambdaHpH = R21 * lambda3 + R22 * Relambda7 - R23*Imlambda7;
664 lambdaHpA = R31 * lambda3 + R32 * Relambda7 - R33*Imlambda7;
665
666
667 //std::cout<<"\033[1;31m lambdaHph = \033[0m "<< lambdaHph <<std::endl;
668 //std::cout<<"\033[1;31m lambdaHpH = \033[0m "<< lambdaHpH <<std::endl;
669 //std::cout<<"\033[1;31m lambdaHpA = \033[0m "<< lambdaHpA <<std::endl;
670
671
672
673// std::cout<<"\033[1;30m Mw = \033[0m "<< Mw <<std::endl;
674
675 //top-quark ratios square
676 double rsqt_h = mt * mt / mH1_2;
677 double rsqt_H = mt * mt / mH2_2;
678 double rsqt_A = mt * mt / mH3_2;
679 double rsqt_Hp = mt * mt / mHp2;
680 double rsqt_W = mt * mt / (Mw * Mw);
681 //bottom-quark ratios square
682 double rsqb_h = mb * mb / mH1_2;
683 double rsqb_H = mb * mb / mH2_2;
684 double rsqb_A = mb * mb / mH3_2;
685 double rsqb_Hp = mb * mb / mHp2;
686 double rsqb_W = mb * mb / (Mw * Mw);
687// double rsqb_Hp = mb * mb / mHp2;
688 //tau-lepton ratios square
689 double rsqtau_h = mTAU * mTAU / mH1_2;
690 double rsqtau_H = mTAU * mTAU / mH2_2;
691 double rsqtau_A = mTAU * mTAU / mH3_2;
692 double rsqtau_Hp = mTAU * mTAU / mHp2;
693 double rsqtau_W = mTAU * mTAU / (Mw * Mw);
694
695 //charged-higgs ratios square
696 double rsqHp_h = mHp2 / mH1_2;
697 double rsqHp_H = mHp2 / mH2_2;
698 double rsqHp_A = mHp2 / mH3_2;
699 double rsqHp_W = mHp2 / (Mw*Mw);
700
701 //neutral scalars ratios square
702 double rsqh_Hp = mH1_2 / mHp2;
703 double rsqH_Hp = mH2_2 / mHp2;
704 double rsqA_Hp = mH3_2 / mHp2;
705 double rsqh_W = mH1_2 / (Mw * Mw);
706 double rsqH_W = mH2_2 / (Mw * Mw);
707 double rsqA_W = mH3_2 / (Mw * Mw);
708
709 //W-boson ratios square
710 double rsqW_h = Mw * Mw / mH1_2;
711 double rsqW_H = Mw * Mw / mH2_2;
712 double rsqW_A = Mw * Mw / mH3_2;
713 double rsqW_Hp = Mw * Mw / mHp2;
714
715
716 //Let's define here the integrals
717 double intgl_x2_1mx_G_rsq_tHp_bHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G(rsqt_Hp,rsqb_Hp);
718 double intgl_x2_1mx_G_rsq_tW_bW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G(rsqt_W,rsqb_W);
719
720 double intgl_x2_1px_G_rsq_tHp_bHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1px_G(rsqt_Hp,rsqb_Hp);
721 double intgl_x2_1px_G_rsq_tW_bW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1px_G(rsqt_W,rsqb_W);
722
723 double intgl_x_1mx2_G_rsq_tHp_bHp = myGTHDM.getMyGTHDMCache()->ip_integral_x_1mx2_G(rsqt_Hp,rsqb_Hp);
724 double intgl_x_1mx2_G_rsq_tW_bW = myGTHDM.getMyGTHDMCache()->ip_integral_x_1mx2_G(rsqt_W,rsqb_W);
725
726 double intgl_x_1mx_1px_G_rsq_tHp_bHp = myGTHDM.getMyGTHDMCache()->ip_integral_x_1mx_1px_G(rsqt_Hp,rsqb_Hp);
727 double intgl_x_1mx_1px_G_rsq_tW_bW = myGTHDM.getMyGTHDMCache()->ip_integral_x_1mx_1px_G(rsqt_W,rsqb_W);
728
729 double intgl_x_1mx2_G_variable_set_0_rsq_tauHp = myGTHDM.getMyGTHDMCache()->ip_integral_x_1mx2_G_variable_set_0(rsqtau_Hp);
730 double intgl_x_1mx2_G_variable_set_0_rsq_tauW = myGTHDM.getMyGTHDMCache()->ip_integral_x_1mx2_G_variable_set_0(rsqtau_W);
731
732
733 double intgl_x2_1mx_G_rsq_WHp_hHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G(rsqW_Hp,rsqh_Hp);
734 double intgl_x2_1mx_G_variable_set_1_rsq_hW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G_variable_set_1(rsqh_W);
735 double intgl_x2_G_rsq_WHp_hHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_G(rsqW_Hp,rsqh_Hp);
736 double intgl_x2_G_variable_set_1_rsq_hW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_G_variable_set_1(rsqh_W);
737
738 double intgl_x2_1mx_G_rsq_HpW_hW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G(rsqHp_W,rsqh_W);
739 double intgl_x2_1mx_G_variable_set_1_rsq_hHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G_variable_set_1(rsqh_Hp);
740
741
742
743 double intgl_x2_1mx_G_rsq_WHp_HHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G(rsqW_Hp,rsqH_Hp);
744 double intgl_x2_1mx_G_variable_set_1_rsq_HW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G_variable_set_1(rsqH_W);
745 double intgl_x2_G_rsq_WHp_HHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_G(rsqW_Hp,rsqH_Hp);
746 double intgl_x2_G_variable_set_1_rsq_HW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_G_variable_set_1(rsqH_W);
747
748 double intgl_x2_1mx_G_rsq_HpW_HW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G(rsqHp_W,rsqH_W);
749 double intgl_x2_1mx_G_variable_set_1_rsq_HHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G_variable_set_1(rsqH_Hp);
750
751
752
753 double intgl_x2_1mx_G_rsq_WHp_AHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G(rsqW_Hp,rsqA_Hp);
754 double intgl_x2_1mx_G_variable_set_1_rsq_AW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G_variable_set_1(rsqA_W);
755 double intgl_x2_G_rsq_WHp_AHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_G(rsqW_Hp,rsqA_Hp);
756 double intgl_x2_G_variable_set_1_rsq_AW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_G_variable_set_1(rsqA_W);
757
758 double intgl_x2_1mx_G_rsq_HpW_AW = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G(rsqHp_W,rsqA_W);
759 double intgl_x2_1mx_G_variable_set_1_rsq_AHp = myGTHDM.getMyGTHDMCache()->ip_integral_x2_1mx_G_variable_set_1(rsqA_Hp);
760
761
762
763
764 //Computation a1
765 double a1_const = (Ale * mMU * mMU / (4 * M_PI * M_PI * M_PI * vev * vev));
766 double a1_SM;
767 double a1_h;
768 double a1_H;
769 double a1_A;
770 double a1_total;
771
772
773 a1_SM = a1_const * (Nc * Qu * Qu * F1twoloopgm2(rsqt_h) + Nc * Qd * Qd * F1twoloopgm2(rsqb_h) + F1twoloopgm2(rsqtau_h));
774
775 a1_h = a1_const * (Nc * Qu * Qu * (yu1R * yl1R * F1twoloopgm2(rsqt_h) + yu1I * yl1I * F1tildetwoloopgm2(rsqt_h))
776 + Nc * Qd * Qd * (yd1R * yl1R * F1twoloopgm2(rsqb_h) + yd1I * yl1I * F1tildetwoloopgm2(rsqb_h))
777 +(yl1R * yl1R * F1twoloopgm2(rsqtau_h) + yl1I * yl1I * F1tildetwoloopgm2(rsqtau_h)));
778
779 a1_H = a1_const * (Nc * Qu * Qu * (yu2R * yl2R * F1twoloopgm2(rsqt_H) + yu2I * yl2I * F1tildetwoloopgm2(rsqt_H))
780 + Nc * Qd * Qd * (yd2R * yl2R * F1twoloopgm2(rsqb_H) + yd2I * yl2I * F1tildetwoloopgm2(rsqb_H))
781 +(yl2R * yl2R * F1twoloopgm2(rsqtau_H) + yl2I * yl2I * F1tildetwoloopgm2(rsqtau_H)));
782
783 a1_A = a1_const * (Nc * Qu * Qu * (yu3R * yl3R * F1twoloopgm2(rsqt_A) + yu3I * yl3I * F1tildetwoloopgm2(rsqt_A))
784 + Nc * Qd * Qd * (yd3R * yl3R * F1twoloopgm2(rsqb_A) + yd3I * yl3I * F1tildetwoloopgm2(rsqb_A))
785 +(yl3R * yl3R * F1twoloopgm2(rsqtau_A) + yl3I * yl3I * F1tildetwoloopgm2(rsqtau_A)));
786
787 a1_total = -a1_SM + a1_h + a1_H + a1_A;
788
789
790
791// std::cout<<"\033[1;34m Nc * Qu * Qu = \033[0m "<< Nc * Qu * Qu <<std::endl;
792
793
794// std::cout<<"\033[1;34m Nc * Qu * Qu * (yu2R * yu2R * F1twoloopgm2(rsqt_H)) = \033[0m "<< Nc * Qu * Qu * (yu2R * yu2R * F1twoloopgm2(rsqt_H)) <<std::endl;
795
796
797// std::cout<<"\033[1;34m F1twoloopgm2(rsqt_H) = \033[0m "<< F1twoloopgm2(rsqt_H) <<std::endl;
798// std::cout<<"\033[1;34m F1twoloopgm2(rsqb_H) = \033[0m "<< F1twoloopgm2(rsqb_H) <<std::endl;
799// std::cout<<"\033[1;34m F1twoloopgm2(rsqtau_H) = \033[0m "<< F1twoloopgm2(rsqtau_H) <<std::endl;
800
801
802// std::cout<<"\033[1;33m a1_SM = \033[0m "<< a1_SM <<std::endl;
803// std::cout<<"\033[1;33m a1_h = \033[0m "<< a1_h <<std::endl;
804// std::cout<<"\033[1;33m a1_H = \033[0m "<< a1_H <<std::endl;
805// std::cout<<"\033[1;33m a1_A = \033[0m "<< a1_A <<std::endl;
806
807
808 //Computation a2
809 double a2_const = (Ale * mMU * mMU / (8 * M_PI * M_PI * M_PI));
810 double a2_SM = 0.; //The SM contribution for this part is zero
811 double a2_h;
812 double a2_H;
813 double a2_A;
814 double a2_total;
815
816
817 a2_h = a2_const * (yl1R / mH1_2) * lambdaHph * F2twoloopgm2(rsqHp_h);
818 a2_H = a2_const * (yl1R / mH2_2) * lambdaHpH * F2twoloopgm2(rsqHp_H);
819 a2_A = a2_const * (yl1R / mH3_2) * lambdaHpA * F2twoloopgm2(rsqHp_A);
820
821 a2_total = -a2_SM + a2_h + a2_H + a2_A;
822
823
824
825 //Computation a3
826 double a3_const = (Ale * mMU * mMU / (8 * M_PI * M_PI * M_PI * vev * vev));
827 double a3_SM;
828 double a3_h;
829 double a3_H;
830 double a3_A;
831 double a3_total;
832
833 a3_SM = a3_const * F3twoloopgm2(rsqW_h);
834 a3_h = a3_const * yl1R * R11 * F3twoloopgm2(rsqW_h);
835 a3_H = a3_const * yl2R * R21 * F3twoloopgm2(rsqW_H);
836 a3_A = a3_const * yl3R * R31 * F3twoloopgm2(rsqW_A);
837
838 a3_total = -a3_SM + a3_h + a3_H + a3_A;
839
840 //Computation a4
841 double a4_const = (Ale * mMU * mMU / (32 * M_PI * M_PI * M_PI * vev * vev * sW2 * (mHp2 - Mw * Mw)));
842 double a4_Hp_tb;
843 double a4_Hp_tt;
844 double a4_Hp_bb;
845 double a4_Hp_bt;
846 double a4_Hp_tau;
847 double a4_total;
848
849 a4_Hp_tb = a4_const * Nc * Vtb2 * (Qu * (sd * (sl.conjugate())).real()*mb*mb)*(intgl_x2_1mx_G_rsq_tHp_bHp - intgl_x2_1mx_G_rsq_tW_bW) ;
850 a4_Hp_tt = a4_const * Nc * Vtb2 * (Qu * (su * (sl.conjugate())).real()*mt*mt)*(intgl_x2_1px_G_rsq_tHp_bHp - intgl_x2_1px_G_rsq_tW_bW) ;
851
852 a4_Hp_bb = a4_const * Nc * Vtb2 * (Qd * (sd * (sl.conjugate())).real()*mb*mb)*(intgl_x_1mx2_G_rsq_tHp_bHp - intgl_x_1mx2_G_rsq_tW_bW) ;
853 a4_Hp_bt = a4_const * Nc * Vtb2 * (Qd * (su * (sl.conjugate())).real()*mt*mt)*(intgl_x_1mx_1px_G_rsq_tHp_bHp - intgl_x_1mx_1px_G_rsq_tW_bW) ;
854
855 a4_Hp_tau = a4_const * (-1 * (sl * (sl.conjugate())).real()*mTAU*mTAU)*(intgl_x_1mx2_G_variable_set_0_rsq_tauHp - intgl_x_1mx2_G_variable_set_0_rsq_tauW) ;
856
857 a4_total = a4_Hp_tb + a4_Hp_tt + a4_Hp_bb + a4_Hp_bt + a4_Hp_tau;
858
859
860
861 //Computation a5
862 gslpp::complex img_i = gslpp::complex::i();
863 double a5_const = (Ale * mMU * mMU / (64 * M_PI * M_PI * M_PI * vev * vev * sW2 * (mHp2 - Mw * Mw)));
864 double a5_h;
865 double a5_H;
866 double a5_A;
867 double a5_total;
868
869 a5_h = a5_const*(sl.conjugate()*R11*(R12 - img_i*R13)).real()*((mHp2+Mw*Mw-mH1_2)*(intgl_x2_1mx_G_rsq_WHp_hHp -
870 intgl_x2_1mx_G_variable_set_1_rsq_hW) - 4*Mw*Mw*(intgl_x2_G_rsq_WHp_hHp - intgl_x2_G_variable_set_1_rsq_hW) );
871
872 a5_H = a5_const*(sl.conjugate()*R21*(R22 - img_i*R23)).real()*((mHp2+Mw*Mw-mH2_2)*(intgl_x2_1mx_G_rsq_WHp_HHp -
873 intgl_x2_1mx_G_variable_set_1_rsq_HW) - 4*Mw*Mw*(intgl_x2_G_rsq_WHp_HHp - intgl_x2_G_variable_set_1_rsq_HW) );
874
875 a5_A = a5_const*(sl.conjugate()*R31*(R32 - img_i*R33)).real()*((mHp2+Mw*Mw-mH3_2)*(intgl_x2_1mx_G_rsq_WHp_AHp -
876 intgl_x2_1mx_G_variable_set_1_rsq_AW) - 4*Mw*Mw*(intgl_x2_G_rsq_WHp_AHp - intgl_x2_G_variable_set_1_rsq_AW) );
877
878 a5_total = a5_h + a5_H + a5_A;
879
880
881
882
883 //Computation a6
884 double a6_const = (Ale * mMU * mMU / (64 * M_PI * M_PI * M_PI * sW2 * (mHp2 - Mw * Mw)));
885 double a6_h;
886 double a6_H;
887 double a6_A;
888 double a6_total;
889 //lambdaHph
890
891 a6_h = a6_const * (sl.conjugate()*(R12 - img_i*R13)).real() * lambdaHph *
892 (- intgl_x2_1mx_G_variable_set_1_rsq_hHp + intgl_x2_1mx_G_rsq_HpW_hW);
893
894 a6_H = a6_const * (sl.conjugate()*(R22 - img_i*R23)).real() * lambdaHpH *
895 (- intgl_x2_1mx_G_variable_set_1_rsq_HHp + intgl_x2_1mx_G_rsq_HpW_HW);
896
897 a6_A = a6_const * (sl.conjugate()*(R32 - img_i*R33)).real() * lambdaHpA *
898 (- intgl_x2_1mx_G_variable_set_1_rsq_AHp + intgl_x2_1mx_G_rsq_HpW_AW);
899
900 a6_total = a6_h + a6_H + a6_A;
901
902
903 double gminus2muNLO = a1_total + a2_total + a3_total + a4_total + a5_total + a6_total;
904
905
906// std::cout<<"\033[1;32m a1_total = \033[0m "<< a1_total <<std::endl;
907// std::cout<<"\033[1;32m a2_total = \033[0m "<< a2_total <<std::endl;
908// std::cout<<"\033[1;32m a3_total = \033[0m "<< a3_total <<std::endl;
909// std::cout<<"\033[1;32m a4_total = \033[0m "<< a4_total <<std::endl;
910// std::cout<<"\033[1;32m a5_total = \033[0m "<< a5_total <<std::endl;
911// std::cout<<"\033[1;32m a6_total = \033[0m "<< a6_total <<std::endl;
912
913
914 return (gminus2muNLO);
915}
GeneralTHDMMatching::F2twoloopgm2
double F2twoloopgm2(const double ratio_sq)
Loop at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{1}{2}*\int_0^1 dx ...
Definition: GeneralTHDMMatching.cpp:426
GeneralTHDMMatching::F1twoloopgm2
double F1twoloopgm2(const double ratio_sq)
Loop at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{ratio_sq}{2}*\int_...
Definition: GeneralTHDMMatching.cpp:373
GeneralTHDMMatching::F3twoloopgm2
double F3twoloopgm2(const double ratio_sq)
Loop at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{1}{2}*\int_0^1 dx ...
Definition: GeneralTHDMMatching.cpp:449
GeneralTHDMMatching::F1tildetwoloopgm2
double F1tildetwoloopgm2(const double ratio_sq)
Loop at two loops (Barr-Zee) from 1502.04199. The complete loop function is \frac{ratio_sq}{2}*\int_...
Definition: GeneralTHDMMatching.cpp:402
Mw
An observable class for the -boson mass.
Definition: Mw.h:22
QCD::TAU
@ TAU
Definition: QCD.h:316

◆ lambdaHHphi()

gslpp::complex GeneralTHDMMatching::lambdaHHphi ( double  lambda3,
double  Relambda7,
double  Imlambda7,
double  Ri1,
double  Ri2,
double  Ri3 
)
virtual

Definition at line 1255 of file GeneralTHDMMatching.cpp.

1255 {
1256 gslpp::complex lambdaHHphi = lambda3 * Ri1 + Relambda7 * Ri2 - Imlambda7*Ri3;
1257 return lambdaHHphi;
1258}

◆ neglog()

gslpp::complex GeneralTHDMMatching::neglog ( gslpp::complex  argument)
virtual

Calculates the log of a negative number.

Calculates the log of a negative number

Returns
log of a negative number

Definition at line 296 of file GeneralTHDMMatching.cpp.

296 {
297
298 gslpp::complex result;
299
300 result = log(argument.abs()) + gslpp::complex::i() * argument.arg();
301
302
303 return result;
304}

◆ negpow()

gslpp::complex GeneralTHDMMatching::negpow ( double  basis,
double  exp 
)
virtual

Calculates the power root of a negative number.

Calculates the power root of a negative number

Returns
power root of a negative number

Definition at line 286 of file GeneralTHDMMatching.cpp.

286 {
287 gslpp::complex result;
288 if (basis > 0) {
289 result = pow(basis, exp);
290 } else {
291 result = pow(-basis, exp)*(cos(M_PI * exp) + gslpp::complex::i() * sin(M_PI * exp));
292 }
293 return result;
294}

◆ negsquareroot()

gslpp::complex GeneralTHDMMatching::negsquareroot ( double  x)
virtual

Calculates amplitudes for \( (g-2)_{\mu} \) at approximate two-loop from [Cherchiglia:2016eui].

Calculates the bosonic NLO contribution to the muon g-2

Returns
Calculates the bosonic NLO contribution to the muon g-2

Calculates amplitudes for \( (g-2)_{\mu} \) at approximate two-loop from [Cherchiglia:2016eui].

Returns
Calculates the square root of a negative number

Calculates the square root of a negative number

Returns
square root of a negative number

Definition at line 276 of file GeneralTHDMMatching.cpp.

276 {
277 gslpp::complex result;
278 if (x > 0) {
279 result = sqrt(x);
280 } else {
281 result = gslpp::complex::i() * sqrt(-x);
282 }
283 return result;
284}

◆ setWCbsg()

gslpp::complex GeneralTHDMMatching::setWCbsg ( int  i,
gslpp::complex  sigmau,
gslpp::complex  sigmad,
double  mt,
double  mhp,
double  mu,
orders  order 
)
Parameters
iint, flag for the caching
sigmau
sigmad
mttop mass
mhpcharged Higgs mass
mumatching scale
order
Returns
return the value of the Wilson coefficients for \( B \rightarrow X_{s} \gamma, l^{+} l^{-} \)

Definition at line 1517 of file GeneralTHDMMatching.cpp.

1517 {
1518 if (su.abs() == sigmau.abs() && su.arg() == sigmau.arg() &&
1519 sd.abs() == sigmad.abs() && sd.arg() == sigmad.arg() &&
1520 mtbsg == mt && mhpbsg == mhp && mubsg == mu) {
1521 switch (order) {
1522 case NNLO:
1523 return ( CWbsgArrayNNLO[i]);
1524 case NLO:
1525 return ( CWbsgArrayNLO[i]);
1526 break;
1527 case LO:
1528 return ( CWbsgArrayLO[i]);
1529 break;
1530 default:
1531 std::stringstream out;
1532 out << order;
1533 throw std::runtime_error("order" + out.str() + "not implemeted");
1534 }
1535 }
1536
1537 su = sigmau;
1538 sd = sigmad;
1539 mtbsg = mt;
1540 mhpbsg = mhp;
1541 mubsg = mu;
1542
1543 double x = mt * mt / mhp / mhp;
1544
1545 double x2 = x*x;
1546 double x3 = x2*x;
1547 double x4 = x3*x;
1548 double x5 = x4*x;
1549 double xm = 1. - x;
1550 double xm2 = xm*xm;
1551 double xm3 = xm2*xm;
1552 double xm4 = xm3*xm;
1553 double xm6 = xm4*xm2;
1554 double xm8 = xm4*xm4;
1555 double xo = 1. - 1. / x;
1556 double xo2 = xo*xo;
1557 double xo4 = xo * xo2*xo;
1558 double xo6 = xo4*xo2;
1559 double xo8 = xo4*xo4;
1560 double Lx = log(x);
1561 double Lx2 = Lx*Lx;
1562 double Lx3 = Lx2*Lx;
1563 double Li2 = gsl_sf_dilog(1. - 1. / x);
1564
1565 double abssu = su.abs();
1566 double abssu2 = abssu*abssu;
1567 gslpp::complex susd = su.conjugate() * sd;
1568
1569 double lstmu = 2. * log(mu / mt);
1570
1571 double n70ct = -((7. - 5. * x - 8. * x2) / (36. * xm3) + (2. * x - 3. * x2) * Lx / (6. * xm4)) * x / 2.;
1572 double n70fr = ((3. * x - 5. * x2) / (6. * xm2) + (2. * x - 3. * x2) * Lx / (3. * xm3)) / 2.;
1573
1574 double n80ct = -((2. + 5. * x - x2) / (12. * xm3) + (x * Lx) / (2. * xm4)) * x / 2.;
1575 double n80fr = ((3. * x - x2) / (2. * xm2) + (x * Lx) / xm3) / 2.;
1576
1577 double n41ct = (-16. * x + 29. * x2 - 7. * x3) / (36. * xm3) + (-2. * x + 3. * x2) * Lx / (6. * xm4);
1578
1579 double n71ct = (797. * x - 5436. * x2 + 7569. * x3 - 1202. * x4) / (486. * xm4) +
1580 (36. * x2 - 74. * x3 + 16. * x4) * Li2 / (9. * xm4) +
1581 ((7. * x - 463. * x2 + 807. * x3 - 63. * x4) * Lx) / (81. * xm3 * xm2);
1582 double cd71ct = (-31. * x - 18. * x2 + 135. * x3 - 14. * x4) / (27. * xm4) + (-28. * x2 + 46. * x3 + 6. * x4) * Lx / (9. * xm3 * xm2);
1583 double n71fr = (28. * x - 52. * x2 + 8. * x3) / (3. * xm3) + (-48. * x + 112. * x2 - 32. * x3) * Li2 / (9. * xm3) +
1584 (66. * x - 128. * x2 + 14. * x3) * Lx / (9. * xm4);
1585 double cd71fr = (42. * x - 94. * x2 + 16. * x3) / (9. * xm3) + (32. * x - 56. * x2 - 12. * x3) * Lx / (9. * xm4);
1586
1587 double n81ct = (1130. * x - 18153. * x2 + 7650. * x3 - 4451. * x4) / (1296. * xm4) +
1588 (30. * x2 - 17. * x3 + 13. * x4) * Li2 / (6. * xm4) +
1589 (-2. * x - 2155. * x2 + 321. * x3 - 468. * x4) * Lx / (216. * xm3 * xm2);
1590 double cd81ct = (-38. * x - 261. * x2 + 18. * x3 - 7. * x4) / (36. * xm4) + (-31. * x2 - 17. * x3) * Lx / (6. * xm3 * xm2);
1591 double n81fr = (143. * x - 44. * x2 + 29. * x3) / (8. * xm3) + (-36. * x + 25. * x2 - 17. * x3) * Li2 / (6. * xm3) +
1592 (165. * x - 7. * x2 + 34. * x3) * Lx / (12. * xm4);
1593 double cd81fr = (81. * x - 16. * x2 + 7. * x3) / (6. * xm3) + (19. * x + 17. * x2) * Lx / (3. * xm4);
1594
1595 double n32ct = (10. * x4 + 30. * x2 - 20. * x) / (27. * xm4) * Li2 +
1596 (30. * x3 - 66. * x2 - 56. * x) / (81. * xm4) * Lx + (6. * x3 - 187. * x2 + 213. * x) / (-81. * xm3);
1597 double cd32ct = (-30. * x2 + 20. * x) / (27. * xm4) * Lx + (-35. * x3 + 145. * x2 - 80. * x) / (-81. * xm3);
1598
1599 double n42ct = (515. * x4 - 906. * x3 + 99. * x2 + 182. * x) / (54. * xm4) * Li2 +
1600 (1030. * x4 - 2763. * x3 - 15. * x2 + 980. * x) / (-108. * xm3 * xm2) * Lx +
1601 (-29467. * x4 + 68142. * x3 - 6717. * x2 - 18134. * x) / (1944. * xm4);
1602 double cd42ct = (-375. * x3 - 95. * x2 + 182. * x) / (-54. * xm3 * xm2) * Lx +
1603 (133. * x4 - 108. * x3 + 4023. * x2 - 2320. * x) / (324. * xm4);
1604
1605 double cd72ct = -(x * (67930. * x4 - 470095. * x3 + 1358478. * x2 - 700243. * x + 54970.)) / (-2187. * xm3 * xm2) +
1606 (x * (10422. * x4 - 84390. * x3 + 322801. * x2 - 146588. * x + 1435.)) / (729. * xm6) * Lx +
1607 (2. * x2 * (260. * x3 - 1515. * x2 + 3757. * x - 1446.)) / (-27. * xm3 * xm2) * Li2;
1608 double ce72ct = (x * (-518. * x4 + 3665. * x3 - 17397. * x2 + 3767. * x + 1843.)) / (-162. * xm3 * xm2) +
1609 (x2 * (-63. * x3 + 532. * x2 + 2089. * x - 1118.)) / (27. * xm6) * Lx;
1610 double cd72fr = -(x * (3790. * x3 - 22511. * x2 + 53614. * x - 21069.)) / (81. * xm4) -
1611 (2. * x * (-1266. * x3 + 7642. * x2 - 21467. * x + 8179.)) / (-81. * xm3 * xm2) * Lx +
1612 (8. * x * (139. * x3 - 612. * x2 + 1103. * x - 342.)) / (27. * xm4) * Li2;
1613 double ce72fr = -(x * (284. * x3 - 1435. * x2 + 4304. * x - 1425.)) / (27. * xm4) -
1614 (2. * x * (63. * x3 - 397. * x2 - 970. * x + 440.)) / (-27. * xm3 * xm2) * Lx;
1615
1616 double cd82ct = -(x * (522347. * x4 - 2423255. * x3 + 2706021. * x2 - 5930609. * x + 148856)) / (-11664. * xm3 * xm2) +
1617 (x * (51948. * x4 - 233781. * x3 + 48634. * x2 - 698693. * x + 2452.)) / (1944. * xm6) * Lx +
1618 (x2 * (481. * x3 - 1950. * x2 + 1523. * x - 2550.)) / (-18. * xm3 * xm2) * Li2;
1619 double ce82ct = (x * (-259. * x4 + 1117. * x3 + 2925. * x2 + 28411. * x + 2366.)) / (-216. * xm3 * xm2) -
1620 (x2 * (139. * x2 + 2938. * x + 2683.)) / (36. * xm6) * Lx;
1621 double cd82fr = -(x * (1463. * x3 - 5794. * x2 + 5543. * x - 15036.)) / (27. * xm4) -
1622 (x * (-1887. * x3 + 7115. * x2 + 2519. * x + 19901.)) / (-54. * xm3 * xm2) * Lx -
1623 (x * (-629. * x3 + 2178. * x2 - 1729. * x + 2196.)) / (18. * xm4) * Li2;
1624 double ce82fr = -(x * (259. * x3 - 947. * x2 - 251. * x - 5973.)) / (36. * xm4)-
1625 (x * (139. * x2 + 2134. * x + 1183.)) / (-18. * xm3 * xm2) * Lx;
1626
1627 double n72ct = 0.;
1628 if (mhp < 50.)
1629 n72ct = -274.2 / x5 - 72.13 * Lx / x5 + 24.41 / x4 - 168.3 * Lx / x4 + 79.15 / x3 -
1630 103.8 * Lx / x3 + 47.09 / x2 - 38.12 * Lx / x2 + 15.35 / x - 8.753 * Lx / x + 3.970;
1631 else if (mhp < mt)
1632 n72ct = 1.283 + 0.7158 * xo + 0.4119 * xo2 + 0.2629 * xo * xo2 + 0.1825 * xo4 +
1633 0.1347 * xo * xo4 + 0.1040 * xo6 + 0.0831 * xo * xo6 + 0.06804 * xo8 +
1634 0.05688 * xo * xo8 + 0.04833 * xo2 * xo8 + 0.04163 * xo * xo2 * xo8 + 0.03625 * xo4 * xo8 +
1635 0.03188 * xo * xo4 * xo8 + 0.02827 * xo6 * xo8 + 0.02525 * xo * xo6 * xo8 + 0.02269 * xo8 * xo8;
1636 else if (mhp < 520.)
1637 n72ct = 1.283 - 0.7158 * xm - 0.3039 * xm2 - 0.1549 * xm3 - 0.08625 * xm4 -
1638 0.05020 * xm3 * xm2 - 0.02970 * xm6 - 0.01740 * xm3 * xm4 - 0.009752 * xm8 -
1639 0.004877 * xm3 * xm6 - 0.001721 * xm2 * xm8 + 0.0003378 * xm3 * xm8 + 0.001679 * xm4 * xm8 +
1640 0.002542 * xm * xm4 * xm8 + 0.003083 * xm6 * xm8 + 0.003404 * xm * xm6 * xm8 + 0.003574 * xm8 * xm8;
1641 else n72ct = -823.0 * x5 + 42.30 * x5 * Lx3 - 412.4 * x5 * Lx2 - 3362 * x5 * Lx -
1642 1492 * x4 - 23.26 * x4 * Lx3 - 541.4 * x4 * Lx2 - 2540 * x4 * Lx -
1643 1158 * x3 - 34.50 * x3 * Lx3 - 348.2 * x3 * Lx2 - 1292 * x3 * Lx -
1644 480.9 * x2 - 20.73 * x2 * Lx3 - 112.4 * x2 * Lx2 - 396.1 * x2 * Lx -
1645 8.278 * x + 0.9225 * x * Lx2 + 4.317 * x*Lx;
1646
1647 double n72fr = 0.;
1648 if (mhp < 50.)
1649 n72fr = -(194.3 / x5 + 101.1 * Lx / x5 - 24.97 / x4 + 168.4 * Lx / x4 - 78.90 / x3 +
1650 106.2 * Lx / x3 - 49.32 / x2 + 38.43 * Lx / x2 - 12.91 / x + 9.757 * Lx / x + 8.088);
1651 else if (mhp < mt)
1652 n72fr = -(12.82 - 1.663 * xo - 0.8852 * xo2 - 0.4827 * xo * xo2 - 0.2976 * xo4 -
1653 0.2021 * xo * xo4 - 0.1470 * xo6 - 0.1125 * xo * xo6 - 0.08931 * xo8 -
1654 0.07291 * xo * xo8 - 0.06083 * xo2 * xo8 - 0.05164 * xo * xo2 * xo8 - 0.04446 * xo4 * xo8 -
1655 0.03873 * xo * xo4 * xo8 - 0.03407 * xo6 * xo8 - 0.03023 * xo * xo6 * xo8 - 0.02702 * xo8 * xo8);
1656 else if (mhp < 400.)
1657 n72fr = -(12.82 + 1.663 * xm + 0.7780 * xm2 + 0.3755 * xm3 + 0.1581 * xm4 +
1658 0.03021 * xm3 * xm2 - 0.04868 * xm6 - 0.09864 * xm3 * xm4 - 0.1306 * xm8 -
1659 0.1510 * xm3 * xm6 - 0.1637 * xm2 * xm8 - 0.1712 * xm3 * xm8 - 0.1751 * xm4 * xm8 -
1660 0.1766 * xm * xm4 * xm8 - 0.1763 * xm6 * xm8 - 0.1748 * xm * xm6 * xm8 - 0.1724 * xm8 * xm8);
1661 else n72fr = -(2828 * x5 - 66.63 * x5 * Lx3 + 469.4 * x5 * Lx2 + 1986 * x5 * Lx +
1662 1480 * x4 + 36.08 * x4 * Lx3 + 323.2 * x4 * Lx2 + 169.9 * x4 * Lx +
1663 166.7 * x3 + 19.73 * x3 * Lx3 - 46.61 * x3 * Lx2 - 826.2 * x3 * Lx -
1664 524.1 * x2 - 8.889 * x2 * Lx3 - 195.7 * x2 * Lx2 - 870.3 * x2 * Lx -
1665 572.2 * x - 20.94 * x * Lx3 - 123.5 * x * Lx2 - 453.5 * x * Lx);
1666
1667 double n82ct = 0.;
1668 if (mhp < 50.)
1669 n82ct = 826.2 / x5 - 300.7 * Lx / x5 + 96.35 / x4 + 91.89 * Lx / x4 - 66.39 / x3 +
1670 78.58 * Lx / x3 - 39.76 / x2 + 20.02 * Lx / x2 - 5.214 / x + 2.278;
1671 else if (mhp < mt)
1672 n82ct = 1.188 + 0.4078 * xo + 0.2002 * xo2 + 0.1190 * xo * xo2 + 0.07861 * xo4 +
1673 0.05531 * xo * xo4 + 0.04061 * xo6 + 0.03075 * xo * xo6 + 0.02386 * xo8 +
1674 0.01888 * xo * xo8 + 0.01520 * xo2 * xo8 + 0.01241 * xo * xo2 * xo8 + 0.01026 * xo4 * xo8 +
1675 0.008575 * xo * xo4 * xo8 + 0.007238 * xo6 * xo8 + 0.006164 * xo * xo6 * xo8 + 0.005290 * xo8 * xo8;
1676 else if (mhp < 600.)
1677 n82ct = 1.188 - 0.4078 * xm - 0.2076 * xm2 - 0.1265 * xm3 - 0.08570 * xm4 -
1678 0.06204 * xm3 * xm2 - 0.04689 * xm6 - 0.03652 * xm3 * xm4 - 0.02907 * xm8 -
1679 0.02354 * xm3 * xm6 - 0.01933 * xm2 * xm8 - 0.01605 * xm3 * xm8 - 0.01345 * xm4 * xm8 -
1680 0.01137 * xm * xm4 * xm8 - 0.009678 * xm6 * xm8 - 0.008293 * xm * xm6 * xm8 - 0.007148 * xm8 * xm8;
1681 else n82ct = -19606 * x5 - 226.7 * x5 * Lx3 - 5251 * x5 * Lx2 - 26090 * x5 * Lx -
1682 9016 * x4 - 143.4 * x4 * Lx3 - 2244 * x4 * Lx2 - 10102 * x4 * Lx -
1683 3357 * x3 - 66.32 * x3 * Lx3 - 779.6 * x3 * Lx2 - 3077 * x3 * Lx -
1684 805.5 * x2 - 22.98 * x2 * Lx3 - 169.1 * x2 * Lx2 - 602.7 * x2 * Lx +
1685 0.7437 * x + 0.6908 * x * Lx2 + 3.238 * x*Lx;
1686
1687 double n82fr = 0.;
1688 if (mhp < 50.)
1689 n82fr = -(-1003 / x5 + 476.9 * Lx / x5 - 205.7 / x4 - 71.62 * Lx / x4 + 62.26 / x3 -
1690 110.7 * Lx / x3 + 63.74 / x2 - 35.42 * Lx / x2 + 10.89 / x - 3.174);
1691 else if (mhp < mt)
1692 n82fr = -(-0.6110 - 1.095 * xo - 0.4463 * xo2 - 0.2568 * xo * xo2 - 0.1698 * xo4 -
1693 0.1197 * xo * xo4 - 0.08761 * xo6 - 0.06595 * xo * xo6 - 0.05079 * xo8 -
1694 0.03987 * xo * xo8 - 0.03182 * xo2 * xo8 - 0.02577 * xo * xo2 * xo8 - 0.02114 * xo4 * xo8 -
1695 0.01754 * xo * xo4 * xo8 - 0.01471 * xo6 * xo8 - 0.01244 * xo * xo6 * xo8 - 0.01062 * xo8 * xo8);
1696 else if (mhp < 520.)
1697 n82fr = -(-0.6110 + 1.095 * xm + 0.6492 * xm2 + 0.4596 * xm3 + 0.3569 * xm4 +
1698 0.2910 * xm3 * xm2 + 0.2438 * xm6 + 0.2075 * xm3 * xm4 + 0.1785 * xm8 +
1699 0.1546 * xm3 * xm6 + 0.1347 * xm2 * xm8 + 0.1177 * xm3 * xm8 + 0.1032 * xm4 * xm8 +
1700 0.09073 * xm * xm4 * xm8 + 0.07987 * xm6 * xm8 + 0.07040 * xm * xm6 * xm8 + 0.06210 * xm8 * xm8);
1701 else n82fr = -(-15961 * x5 + 1003 * x5 * Lx3 - 2627 * x5 * Lx2 - 29962 * x5 * Lx -
1702 11683 * x4 + 54.66 * x4 * Lx3 - 2777 * x4 * Lx2 - 17770 * x4 * Lx -
1703 6481 * x3 - 40.68 * x3 * Lx3 - 1439 * x3 * Lx2 - 7906 * x3 * Lx -
1704 2943 * x2 - 31.83 * x2 * Lx3 - 612.6 * x2 * Lx2 - 2770 * x2 * Lx -
1705 929.8 * x - 19.80 * x * Lx3 - 174.7 * x * Lx2 - 658.4 * x * Lx);
1706
1707
1708 switch (order) {
1709 case NNLO:
1710 CWbsgArrayNNLO[2] = (n32ct + cd32ct * lstmu) * abssu2;
1711 CWbsgArrayNNLO[3] = (n42ct + cd42ct * lstmu) * abssu2;
1712 CWbsgArrayNNLO[4] = 2. / 15. * n41ct * abssu2 - CWbsgArrayNNLO[2] / 10.;
1713 CWbsgArrayNNLO[5] = 0.25 * n41ct * abssu2 - CWbsgArrayNNLO[2] * 3. / 16.;
1714 CWbsgArrayNNLO[6] = -(n72ct + cd72ct * lstmu + ce72ct * lstmu * lstmu) * abssu2
1715 + (n72fr + cd72fr * lstmu + ce72fr * lstmu * lstmu) * susd
1716 - 1. / 3. * CWbsgArrayNNLO[2] - 4. / 9. * CWbsgArrayNNLO[3] - 20. / 3. * CWbsgArrayNNLO[4] - 80. / 9. * CWbsgArrayNNLO[5];
1717 CWbsgArrayNNLO[7] = -(n82ct + cd82ct * lstmu + ce82ct * lstmu * lstmu) * abssu2
1718 + (n82fr + cd82fr * lstmu + ce82fr * lstmu * lstmu) * susd
1719 + CWbsgArrayNNLO[2] - 1. / 6. * CWbsgArrayNNLO[3] - 20. * CWbsgArrayNNLO[4] - 10. / 3. * CWbsgArrayNNLO[5];
1720 case NLO:
1721 CWbsgArrayNLO[3] = n41ct * abssu2;
1722 CWbsgArrayNLO[6] = (n71ct + cd71ct * lstmu) * abssu2 - (n71fr + cd71fr * lstmu) * susd - 4. / 9. * CWbsgArrayNLO[3];
1723 CWbsgArrayNLO[7] = (n81ct + cd81ct * lstmu) * abssu2 - (n81fr + cd81fr * lstmu) * susd - 1. / 6. * CWbsgArrayNLO[3];
1724 case LO:
1725 CWbsgArrayLO[6] = n70ct * abssu2 - n70fr * susd;
1726 CWbsgArrayLO[7] = n80ct * abssu2 - n80fr * susd;
1727 break;
1728 default:
1729 std::stringstream out;
1730 out << order;
1731 throw std::runtime_error("order" + out.str() + "not implemeted");
1732 }
1733
1734 /*std::cout << "CWbsgArrayLO[6] = " << CWbsgArrayLO[6] << std::endl;
1735 std::cout << "CWbsgArrayLO[7] = " << CWbsgArrayLO[7] << std::endl;
1736 std::cout << "CWbsgArrayNLO[3] = " << CWbsgArrayNLO[3] << std::endl;
1737 std::cout << "CWbsgArrayNLO[6] = " << CWbsgArrayNLO[6] << std::endl;
1738 std::cout << "CWbsgArrayNLO[7] = " << CWbsgArrayNLO[7] << std::endl;
1739 std::cout << "CWbsgArrayNNLO[2] = " << CWbsgArrayNNLO[2] << std::endl;
1740 std::cout << "CWbsgArrayNNLO[3] = " << CWbsgArrayNNLO[3] << std::endl;
1741 std::cout << "CWbsgArrayNNLO[4] = " << CWbsgArrayNNLO[4] << std::endl;
1742 std::cout << "CWbsgArrayNNLO[5] = " << CWbsgArrayNNLO[5] << std::endl;
1743 std::cout << "CWbsgArrayNNLO[6] = " << CWbsgArrayNNLO[6] << std::endl;
1744 std::cout << "CWbsgArrayNNLO[7] = " << CWbsgArrayNNLO[7] << std::endl;*/
1745
1746
1747 switch (order) {
1748 case NNLO:
1749 return ( CWbsgArrayNNLO[i]);
1750 case NLO:
1751 return ( CWbsgArrayNLO[i]);
1752 break;
1753 case LO:
1754 return ( CWbsgArrayLO[i]);
1755 break;
1756 default:
1757 std::stringstream out;
1758 out << order;
1759 throw std::runtime_error("order" + out.str() + "not implemeted");
1760 }
1761}
GeneralTHDMMatching::CWbsgArrayLO
gslpp::complex CWbsgArrayLO[8]
Definition: GeneralTHDMMatching.h:373
GeneralTHDMMatching::CWbsgArrayNNLO
gslpp::complex CWbsgArrayNNLO[8]
Definition: GeneralTHDMMatching.h:373
GeneralTHDMMatching::CWbsgArrayNLO
gslpp::complex CWbsgArrayNLO[8]
Definition: GeneralTHDMMatching.h:373
Li2
cd Li2(cd x)
Definition: hpl.h:1011

◆ updateGTHDMParameters()

void GeneralTHDMMatching::updateGTHDMParameters ( )

Definition at line 28 of file GeneralTHDMMatching.cpp.

28 {
29 GF = myGTHDM.getGF();
30 mMU = myGTHDM.getLeptons(StandardModel::MU).getMass();
31}

Member Data Documentation

◆ CWbsgArrayLO

gslpp::complex GeneralTHDMMatching::CWbsgArrayLO[8]
private

Definition at line 373 of file GeneralTHDMMatching.h.

◆ CWbsgArrayNLO

gslpp::complex GeneralTHDMMatching::CWbsgArrayNLO[8]
private

Definition at line 373 of file GeneralTHDMMatching.h.

◆ CWbsgArrayNNLO

gslpp::complex GeneralTHDMMatching::CWbsgArrayNNLO[8]
private

Definition at line 373 of file GeneralTHDMMatching.h.

◆ GF

double GeneralTHDMMatching::GF
private

Definition at line 372 of file GeneralTHDMMatching.h.

◆ mcBMll

WilsonCoefficient GeneralTHDMMatching::mcBMll
private

Definition at line 370 of file GeneralTHDMMatching.h.

◆ mcbsg

WilsonCoefficient GeneralTHDMMatching::mcbsg
private

Definition at line 370 of file GeneralTHDMMatching.h.

◆ mcbsmm

WilsonCoefficient GeneralTHDMMatching::mcbsmm
private

Definition at line 370 of file GeneralTHDMMatching.h.

◆ mcdbs2

WilsonCoefficient GeneralTHDMMatching::mcdbs2
private

Definition at line 370 of file GeneralTHDMMatching.h.

◆ mcgminus2mu

WilsonCoefficient GeneralTHDMMatching::mcgminus2mu
private

Definition at line 370 of file GeneralTHDMMatching.h.

◆ mculeptonnu

WilsonCoefficient GeneralTHDMMatching::mculeptonnu
private

Definition at line 370 of file GeneralTHDMMatching.h.

◆ mhpbsg

double GeneralTHDMMatching::mhpbsg
private

Definition at line 374 of file GeneralTHDMMatching.h.

◆ mMU

double GeneralTHDMMatching::mMU
private

Definition at line 372 of file GeneralTHDMMatching.h.

◆ mtbsg

double GeneralTHDMMatching::mtbsg
private

Definition at line 374 of file GeneralTHDMMatching.h.

◆ mubsg

double GeneralTHDMMatching::mubsg
private

Definition at line 374 of file GeneralTHDMMatching.h.

◆ myCKM

gslpp::matrix<gslpp::complex> GeneralTHDMMatching::myCKM
private

Definition at line 369 of file GeneralTHDMMatching.h.

◆ myGTHDM

const GeneralTHDM& GeneralTHDMMatching::myGTHDM
private

Definition at line 367 of file GeneralTHDMMatching.h.

◆ PolyLog

const Polylogarithms GeneralTHDMMatching::PolyLog
private

Definition at line 377 of file GeneralTHDMMatching.h.

◆ sd

gslpp::complex GeneralTHDMMatching::sd
private

Definition at line 375 of file GeneralTHDMMatching.h.

◆ sl

gslpp::complex GeneralTHDMMatching::sl
private

Definition at line 375 of file GeneralTHDMMatching.h.

◆ su

gslpp::complex GeneralTHDMMatching::su
private

Definition at line 375 of file GeneralTHDMMatching.h.

The documentation for this class was generated from the following files: