11 #include <gsl/gsl_integration.h>
25 mcdbd2Chi0Chi0(5,
NDR,
NLO),
29 mcdbd2ChiChiT(5,
NDR,
NLO),
30 mcdbd2Chi0Chi0T(5,
NDR,
NLO),
36 mcdbs2Chi0Chi0(5,
NDR,
NLO),
40 mcdbs2ChiChiT(5,
NDR,
NLO),
41 mcdbs2Chi0Chi0T(5,
NDR,
NLO),
47 mcdk2Chi0Chi0(5,
NDR,
NLO),
52 mcdk2Chi0Chi0T(5,
NDR,
NLO),
58 mcdd2Chi0Chi0(5,
NDR,
NLO),
63 mcdd2Chi0Chi0T(5,
NDR,
NLO),
67 mcmueconv(8,
NDR,
LO),
68 mcgminus2mu(2,
NDR,
LO),
125 AmpTauA1LN(4, 6, 0.),
126 AmpTauA1RN(4, 6, 0.),
127 AmpTauA1LC(2, 3, 0.),
128 AmpTauA1RC(2, 3, 0.),
151 gamULCKMMD(6, 3, 0.),
152 gamURCKMMU(6, 3, 0.),
155 gamULgamULdag(6, 6, 0.),
156 gamDRgamDRdag(6, 6, 0.),
161 Lambda0EpsYCache(3,3,0.),
162 DeltaDL_Cache(3,3,0.),
165 myCKM_cache(3, 3, 0.),
166 VUDHH_cache(6, 6, 0.),
167 DeltaMd_cache(3, 3, 0.),
224 for (
int i = 0; i < 6; i++) {
225 for (
int j = 0; j < 3; j++) {
226 gslpp::complex sum_k_gamULCKM = 0.,sum_k_gamURCKMMU = 0., sum_k_gamULCKMMD = 0., sum_k_gamDRMD = 0., sum_k_gamDLMD = 0.;
227 for (
int k = 0; k < 3; k++) {
229 for (
int l = 0; l < 3; l++) {
233 sum_k_gamULCKM +=
myRu(i, k) *
myCKM(k, j);
234 sum_k_gamURCKMMU +=
myRu(i, k + 3) * temp_QmassUCKM;
235 sum_k_gamULCKMMD +=
myRu(i, k) * temp_CKMQmassD;
239 gamULCKM.assign(i, j, sum_k_gamULCKM);
242 gamDRMD.assign(i, j, sum_k_gamDRMD);
243 gamDLMD.assign(i, j, sum_k_gamDLMD);
247 for (
int i = 0; i < 6; i++) {
248 for (
int j = 0; j < 6; j++) {
250 for (
int k = 0; k < 3; k++) {
251 sum_k_gamULgamULdag +=
myRu(i, k) *
myRu(j, k).conjugate();
252 sum_k_gamDRgamDRdag +=
myRd(i, k + 3) *
myRd(j, k + 3).conjugate();
259 for (
int i = 0; i < 6; i++) {
260 for (
int j = 0; j < 3; j++) {
261 for (
int k = 0; k < 2; k++) {
268 for (
int i = 0; i < 6; i++) {
269 for (
int j = 0; j < 3; j++) {
270 for (
int k = 0; k < 4; k++) {
285 if ((fabs(z) < SUSYLEPS3) || (fabs(t) < SUSYLEPS3)) {
289 return (
sqrt(z * t) *
Dk(x, y, z, t, 0));
295 return ( (
log(a) -
log(b)) / (4. * (-a + b)) );
303 return ( (a * (-b + c) *
log(a) + b * (a - c) *
log(b)
304 + (-a + b) * c *
log(c)) / (4. * (a - b)*(-a + c)*(-b + c)));
308 return (((-b + c) *
log(a) + (a - c) *
log(b) + (-a + b) *
log(c)) /
309 ((a - b)*(-a + c)*(-b + c)));
313 throw std::runtime_error(
"Error in DL0(a,b,c,k) in SUSYMatching.cpp ");
316 return (EXIT_FAILURE);
321 if (fabs(a) < SUSYLEPS) {
323 throw std::runtime_error(
"Error in DL(a,b,c,0) in SUSYMatching.cpp because the limit DL(0,b,c, k = 0) is singular");
325 else if (fabs(b) < SUSYLEPS) {
327 return ( (-a + c + a *
log(a) - a *
log(c)) / (a * (a - c) * (a - c)));
329 else if (fabs(c) < SUSYLEPS) {
331 return ( (-a + b + a *
log(a) - a *
log(b)) / (a * (a - b) * (a - b)));
334 return( (a*a-b*c)/(a-b)/(a-b)/(a-c)/(a-c)*
log(a) +
335 b/(a-b)/(a-b)/(c-b)*
log(b) - c/(a-c)/(a-c)/(c-b)*
log(c)
340 if (fabs(a) < SUSYLEPS) {
342 return ( (-
log(b) +
log(c)) / (4. * (b - c)));
344 else if (fabs(b) < SUSYLEPS) {
346 return ( (-a + c + c *
log(a) - c *
log(c)) / (4. * (a - c) * (a - c)));
348 else if (fabs(c) < SUSYLEPS) {
350 return ( (-a + b + b *
log(a) - b *
log(b)) / (4. * (a - b) * (a - b)));
353 return ( a*(-2.*b*c + a*(b + c))/4./(a - b)/(a - b)/(a - c)/(a - c)*
log(a)
354 + b*b/4./(a - b)/(a - b)/(c - b)*
log(b) - c*c/4./(a - c)/(a - c)/
355 (c - b)*
log(c) + a/4./(a - b)/(c - a) );
359 throw std::runtime_error(
"Error in DL(a,b,c,k) in SUSYMatching.cpp ");
362 return (EXIT_FAILURE);
368 if (fabs(a) < SUSYLEPS) {
370 throw std::runtime_error(
"Error in DLL(a,b,k = 0) in SUSYMatching.cpp because the limit DLL(0,b, k = 0) is singular");
373 if (fabs(b) < SUSYLEPS) {
374 return (1. / (2. * a * a));
377 return ( -(-a * a + b * b + 2 * a * b *
log(a) - 2 * a * b *
log(b)) /
378 (2. * a * (a - b) * (a - b) * (a - b)));
382 if (fabs(a) < SUSYLEPS) {
384 throw std::runtime_error(
"Error in DLL(a,b,k = 2) in SUSYMatching.cpp because the limit DLL(0,b, k = 2) is singular");
387 if (fabs(b) < SUSYLEPS) {
389 return ( -1. / (8. * a));
392 return (a * a - 4 * a * b + 3 * b * b + 2 * b * b *
log(a)
393 - 2 * b * b *
log(b)) / (8. * (-a + b) * (-a + b) * (-a + b));
397 throw std::runtime_error(
"Error in DLL(a,b,k) in SUSYMatching.cpp ");
400 return (EXIT_FAILURE);
406 if ( (fabs(a) < SUSYLEPS) || (fabs(b) < SUSYLEPS) ) {
408 std::cout <<
"MChargini = " <<
mySUSY.
getMch() << std::endl;
409 std::cout <<
"MNeutralini = " <<
mySUSY.
getMneu() << std::endl;
410 std::cout <<
"MD2squarks = " <<
mySUSY.
getMsd2() << std::endl;
411 std::cout <<
"MU2squarks = " <<
mySUSY.
getMsu2() << std::endl;
412 std::cout <<
"Mgluino = " <<
mySUSY.
getM3() << std::endl;
415 std::cout <<
"Mup(Q_S) = " <<
mySUSYMQ(0) << std::endl;
416 std::cout <<
"Mdown(Q_S) = " <<
mySUSYMQ(1) << std::endl;
417 std::cout <<
"Mc(Q_S) = " <<
mySUSYMQ(2) << std::endl;
418 std::cout <<
"Ms(Q_S) = " <<
mySUSYMQ(3) << std::endl;
419 std::cout <<
"Mtop(Q_S) = " <<
mySUSYMQ(4) << std::endl;
420 std::cout <<
"Mb(Q_S) = " <<
mySUSYMQ(5) << std::endl;
423 throw std::runtime_error(
"Error in DLLp function, because the limits D0(0,0,b,b) and D0(a,a,0,0) are singular ");
425 return (-2 * a + 2 * b + (a + b) *
log(a) - (a + b) *
log(b))
426 / ((a - b)*(a - b)*(a - b));
430 if ( fabs(a) < SUSYLEPS ){
432 return (-1. / (4. * b ));
434 else if ( fabs(b) < SUSYLEPS ){
436 return ( -1. / (4. * a) );
439 return (-a * a + b * b + 2 * a * b *
log(a) - 2 * a * b *
log(b))
440 / (4. * (a - b)*(a - b)*(a - b));
444 throw std::runtime_error(
"Error in DLLp(a,b,k) in SUSYMatching.cpp ");
447 return (EXIT_FAILURE);
453 return (1. / (6. * a * a));
456 return (-1. / (12. * a));
460 throw std::runtime_error(
"Error in DLLL(a,k) in SUSYMatching.cpp ");
463 return (EXIT_FAILURE);
472 if ((fabs(1. - y / x) < SUSYLEPS) && (fabs(1. - z / x) < SUSYLEPS) &&
473 (fabs(1. - t / x) < SUSYLEPS)) {
478 if ((fabs(1. - y / x) < SUSYLEPS) && (fabs(1. - z / x) < SUSYLEPS)) {
481 if ((fabs(1. - y / x) < SUSYLEPS) && (fabs(1. - t / x) < SUSYLEPS)) {
484 if ((fabs(1. - z / x) < SUSYLEPS) && (fabs(1. - t / x) < SUSYLEPS)) {
487 if ((fabs(1. - z / y) < SUSYLEPS) && (fabs(1. - t / y) < SUSYLEPS)) {
492 if (fabs(1. - y / x) < SUSYLEPS) {
493 if ((fabs(1. - t / z) < SUSYLEPS)) {
494 return DLLp(x, z, k);
497 return DL(x, z, t, k);
499 if (fabs(1. - z / x) < SUSYLEPS) {
500 if (fabs(1. - t / y) < SUSYLEPS) {
501 return DLLp(x, y, k);
504 return DL(x, y, t, k);
506 if (fabs(1. - t / x) < SUSYLEPS) {
507 if (fabs(1. - z / y) < SUSYLEPS) {
508 return DLLp(x, y, k);
511 return DL(x, z, y, k);
513 if ((fabs(1. - z / y) < SUSYLEPS)) {
514 return DL(y, x, t, k);
516 if ((fabs(1. - t / y) < SUSYLEPS)) {
517 return DL(y, z, x, k);
519 if ((fabs(1. - z / t) < SUSYLEPS)) {
521 return (
DL(x, y, z, k));
526 if (fabs(x) < SUSYLEPS) {
528 return (
DL0(y, z, t, k));
530 if (fabs(y) < SUSYLEPS) {
532 return (
DL0(x, z, t, k));
534 if (fabs(z) < SUSYLEPS) {
536 return (
DL0(x, y, t, k));
538 if (fabs(t) < SUSYLEPS) {
540 return (
DL0(x, y, z, k));
546 return x *
log(x) / ((t - x)*(z - x)*(y - x)) + y *
log(y)
547 / ((t - y)*(z - y)*(x - y)) + z *
log(z)
548 / ((t - z)*(x - z)*(y - z)) +
549 t *
log(t) / ((x - t)*(z - t)*(y - t));
552 return (((t * t *
log(t)) / ((-t + x)*(-t + y)*(-t + z))
553 + (x * x *
log(x)) / ((t - x)*(-x + y)*(-x + z))
554 + (y * y *
log(y)) / ((t - y)*(x - y)*(-y + z))
555 + (z * z *
log(z)) / ((t - z)*(x - z)*(y - z)))/4.);
559 throw std::runtime_error(
"Error in Dk(x,y,z,t,k) in SUSYMatching.cpp ");
562 return (EXIT_FAILURE);
572 return ((a - b - b *
log(a) + b *
log(b)) / ((a - b)*(a - b)));
577 return (-
log(
mu2R) + (a * (a - b) + a * (a - 2 * b) *
log(a) + b * b *
log(b)) /
582 throw std::runtime_error(
"Error in CL(a,b,k) in SUSYMatching.cpp ");
585 return (EXIT_FAILURE);
591 return (1. / (2. * a));
600 throw std::runtime_error(
"Error in CLL(a,k) in SUSYMatching.cpp ");
603 return (EXIT_FAILURE);
612 if ((fabs(1. - y / x) < SUSYLEPS)&(fabs(1. - z / x) < SUSYLEPS)) {
617 else if ((fabs(1. - y / x) < SUSYLEPS)) {
622 else if ((fabs(1. - z / x) < SUSYLEPS)) {
627 else if ((fabs(1. - z / y) < SUSYLEPS)) {
634 return ((y *
log(y / x)) / ((x - y)*(-y + z)) + (z *
log(z / x)) /
638 return (-
log(
mu2R) +
log(x) + y * y *
log(y / x) / ((x - y)*(-y + z))
639 + z * z *
log(z / x) / ((x - z)*(y - z)));
643 throw std::runtime_error(
"Error in Ck(x,y,z,k) in SUSYMatching.cpp ");
646 return (EXIT_FAILURE);
656 return (1. +
log(a));
663 throw std::runtime_error(
"Error in BL(a,k) in SUSYMatching.cpp ");
666 return (EXIT_FAILURE);
672 if((fabs(x) < SUSYLEPS2)&&(k==0))
return (
log(y));
673 if((fabs(y) < SUSYLEPS2)&&(k==0))
return (
log(x));
676 if (fabs((1. - y / x)) < SUSYLEPS2) {
681 return (
log(x) + (y *
log(x / y)) / (x - y));
685 return (1. / 4. + 0.5 *
Ck(x, y, y, 2));
689 throw std::runtime_error(
"Error in Bk(x,y,k) in SUSYMatching.cpp ");
692 return (EXIT_FAILURE);
720 for (J = 0; J < 3; J++) {
721 for (I = 0; I < 3; I++) {
725 for (k = 0; k < 6; k++) {
727 temp += -2. / 3. *
Als / M_PI *
Mg *
myRd(k, J + 3).conjugate()
728 *
myRd(k, I) *
Bk(
Mg *
Mg, myMD2Squarks(k), 0);
732 for (k = 0; k < 6; k++) {
733 for (l = 0; l < 4; l++) {
742 for (k = 0; k < 6; k++) {
743 for (l = 0; l < 2; l++) {
765 for (
int J = 0; J < 3; J++) {
782 for( I = 0;I < 3;I++){
783 for( J = 0;J < 3;J++){
812 for (I = 0; I < 3; I++) {
813 for (J = 0; J < 3; J++) {
818 (mdJ * mdJ - mdI * mdI));
851 for (I = 0; I < 3; I++) {
852 for (J = 0; J < 3; J++) {
854 for (l = 0; l < 3; l++) {
861 Delta_CKM_IJ.
assign(0.,0.,0);
881 for (b = 0; b < 3; b++) {
882 for (k = 0; k < 6; k++) {
883 for (j = 0; j < 2; j++) {
886 for (l = 0; l < 3; l++) {
887 VdUCL_bkj += -
gW *
myRu(k, l) *
myCKM(l, b) *
myV(j, 0).conjugate()
893 VdUCL_bkj.
assign(0., 0., 0);
916 for (b = 0; b < 3; b++) {
917 for (k = 0; k < 6; k++) {
918 for (j = 0; j < 2; j++) {
919 for (l = 0; l < 3; l++) {
922 VdUCR_bkj += Ydb *
myCKM(l, b) *
myRu(k, l) *
myU(j, 1);
924 for (p = 0; p < 3; p++) {
933 VdUCR_bkj += Ydb *
myRu(k, l) *
myU(j, 1)
937 else{
throw std::runtime_error(
"Wrong flag assigned to vertex VdUCR_bkj å");}
942 VdUCR_bkj.
assign(0., 0., 0);
968 double sW2 = 1.0 - cW2;
969 double CosThetaW =
sqrt(cW2);
970 double SinThetaW =
sqrt(sW2);
973 for (b = 0; b < 3; b++) {
974 for (k = 0; k < 6; k++) {
975 for (j = 0; j < 4; j++) {
979 VdDNL_bkj += (-
gW /
sqrt(2.) *
myRd(k, b) * (1. / 3.
980 * SinThetaW / CosThetaW *
myN(j, 0).conjugate() -
984 *
myRd(k, b + 3) *
myN(j, 2).conjugate());
986 for (l = 0; l < 3; l++) {
990 VdDNL_bkj += (-
gW /
sqrt(2.) *
myRd(k, l) * (1. / 3.
991 * SinThetaW / CosThetaW *
myN(j, 0).conjugate() -
995 *
myRd(k, l + 3) *
myN(j, 2).conjugate())
1002 VdDNL_bkj += -
gW /
sqrt(2.) *
myRd(k, b) * (1. / 3.
1003 * SinThetaW / CosThetaW *
myN(j, 0).conjugate() -
1006 *
myRd(k, b + 3) *
myN(j, 2).conjugate();
1011 VdDNL_bkj.
assign(0., 0., 0);
1018 if(flag != 0 && flag != 1)
throw std::runtime_error(
"Error in Comp_VdDNL(flag) in SUSYMatching.cpp ");
1040 double sW2 = 1.0 - cW2;
1042 double CosThetaW =
sqrt(cW2);
1043 double SinThetaW =
sqrt(sW2);
1046 for (b = 0; b < 3; b++) {
1047 for (k = 0; k < 6; k++) {
1048 for (j = 0; j < 4; j++) {
1052 VdDNR_bkj += -
sqrt(2.) / 3. *
gW * SinThetaW / CosThetaW *
1058 for (l = 0; l < 3; l++) {
1062 VdDNR_bkj += (-
sqrt(2.) / 3. *
gW * SinThetaW / CosThetaW *
1071 VdDNR_bkj += -
sqrt(2.) / 3. *
gW * SinThetaW / CosThetaW *
1077 VdDNR_bkj.
assign(0., 0., 0);
1083 if(flag != 0 && flag != 1)
throw std::runtime_error(
"Error in Comp_VdDNR(flag) in SUSYMatching.cpp ");
1106 for (b = 0; b < 3; b++) {
1107 for (k = 0; k < 6; k++) {
1108 for (j = 0; j < 2; j++) {
1110 for (I = 0; I < 3; I++) {
1112 VuDCL_bkj += -(
gW *
myRd(k, I) *
myU(j, 0).conjugate()
1113 - YdI *
myRd(k, I + 3) *
myU(j, 1).conjugate()) *
1114 myCKM(b, I).conjugate();
1118 VuDCL_bkj.
assign(0., 0., 0);
1143 for (b = 0; b < 3; b++) {
1144 for (k = 0; k < 6; k++) {
1145 for (j = 0; j < 2; j++) {
1149 for (I = 0; I < 3; I++) {
1151 VuDCR_bkj += mySUSYCKM(b, I).conjugate() * Yub *
myRd(k, I) *
myV(j, 1) ;
1156 VuDCR_bkj.
assign(0., 0., 0);
1171 if (Dmixingflag == 0) {
1173 return (
VdUCL(b, k, j));
1175 else if (Dmixingflag == 1) {
1177 return (
VuDCL(b, k, j));
1181 throw std::runtime_error(
"Error in VdUCL(b,k,j,flag,Dmixingflag) in SUSYMatching.cpp ");
1187 if (Dmixingflag == 0) {
1189 return (
VdUCR(b, k, j, flag));
1191 else if (Dmixingflag == 1) {
1193 return (
VuDCR(b, k, j));
1197 throw std::runtime_error(
"Error in VdUCR(b,k,j,flag,Dmixingflag) in SUSYMatching.cpp ");
1213 for (b = 0; b < 3; b++) {
1214 for (k = 0; k < 6; k++) {
1215 for (j = 0; j < 4; j++) {
1219 VuUNL_bkj += -1. /
sqrt(2.) *
gW *
myRu(k, b) * (1. / (3. * TanThetaW) *
1220 myN(j, 0).conjugate() +
myN(j, 1).conjugate())
1221 - Yub *
myRu(k, b + 3) *
myN(j, 3).conjugate();
1223 VuUNR_bkj += 2. *
sqrt(2.) / 3. *
gW * TanThetaW *
myRu(k, b + 3) *
1227 VuUNL_bkj.
assign(0., 0., 0);
1229 VuUNR_bkj.
assign(0.,0.,0);
1241 if (chirality.compare(
"L") == 0) {
1246 else if (chirality.compare(
"R") == 0) {
1253 std::cout <<
" VuUN error in SUSYMatching.cpp" << std::endl;
1257 return (EXIT_FAILURE);
1266 if (Dmixingflag == 0) {
1268 return (
VdDNL(b, k, j, flag));
1271 else if (Dmixingflag == 1) {
1273 return (
VuUN(b, k, j,
"L"));
1277 throw std::runtime_error(
"Error in VdDNL(b,k,j,flag,Dmixingflag) in SUSYMatching.cpp ");
1280 return (EXIT_FAILURE);
1285 if (Dmixingflag == 0) {
1287 return (
VdDNR(b, k, j, flag));
1290 else if (Dmixingflag == 1) {
1292 return (
VuUN(b, k, j,
"R"));
1296 throw std::runtime_error(
"Error in VdDNR(b,k,j,flag,Dmixingflag) in SUSYMatching.cpp ");
1299 return (EXIT_FAILURE);
1324 for (i = 0; i < 3; i++) {
1325 for (j = 0; j < 3; j++) {
1332 for (k = 0; k < 3; k++) {
1373 for (i = 0; i < 6; i++) {
1374 for (j = 0; j < 6; j++) {
1377 for (I = 0; I < 3; I++) {
1378 for (J = 0; J < 3; J++) {
1383 (1 +
Eps_J(I) *
tanb) * mySUSYCKM(J, I) *
myRd(j, I + 3).conjugate() *
1390 * ZH(1, 0)) * mySUSYCKM(J, I) *
myRd(j, I).conjugate() *
myRu(i, J)
1395 ZH(1,0)*(myTU(J,0)*mySUSYCKM(0,I) + myTU(J,1)*mySUSYCKM(1,I)
1396 + myTU(J,2)*mySUSYCKM(2,I)))
1397 *
myRu(i, J + 3) *
myRd(j, I).conjugate()
1399 + (ZH(0,0) * (myTD(I,0).conjugate() * mySUSYCKM(J,0)
1400 + myTD(I,1).conjugate() * mySUSYCKM(J,1)
1401 + myTD(I,2).conjugate() * mySUSYCKM(J,2))
1404 myRd(j, I + 3).conjugate();
1413 VUDHijH.
assign(0., 0., 0);
1442 for (m = 0; m < 6; m++) {
1443 for (l = 0; l < 6; l++) {
1445 DFHL_ji +=
VUDHH(m,l) *(-2 *
Als / (3. * M_PI) *
Mg
1446 *
myRu(m, j + 3).conjugate() *
myRd(l, i) *
1447 Ck(
Mg *
Mg, myMU2Squarks(m), myMD2Squarks(l), 0)
1449 - 1. / (16. * M_PI * M_PI) * Yuj * Ydi *
1450 myRu(m, j).conjugate() *
myRd(l, i + 3) *
1464 for (i = 0; i < 3; i++) {
1465 for (j = 0; j < 3; j++) {
1488 throw std::runtime_error(
"Error in PLRk(j,i,k) in SUSYMatching.cpp ");
1491 return (EXIT_FAILURE);
1503 throw std::runtime_error(
"Error in PRLk(j,i,k) in SUSYMatching.cpp ");
1512 if (Dmixingflag == 0) {
1513 return (
PRLk(j, i, k));
1516 else if (Dmixingflag == 1) {
1517 return (
PLRk(i, j, k).conjugate());
1520 throw std::runtime_error(
"Error in PRLk(j,i,k,Dmixingflag) in SUSYMatching.cpp ");
1528 if (Dmixingflag == 0) {
1529 return (
PLRk(j, i, k));
1532 else if (Dmixingflag == 1) {
1533 return (
PRLk(i, j, k).conjugate());
1536 throw std::runtime_error(
"Error in PLRk(j,i,k,Dmixingflag) in SUSYMatching.cpp ");
1539 return (EXIT_FAILURE);
1548 double tan2alpha = 2. *
tanb / (1 -
tanb *
tanb) * (M2A + M2Z) / (M2A - M2Z);
1549 double tana = (-1. +
sqrt(1. + tan2alpha * tan2alpha)) / tan2alpha;
1550 double sina = tana /
sqrt(1. + tana * tana);
1551 double cosa = 1. /
sqrt(1 + tana * tana);
1569 throw std::runtime_error(
"Error in xdS(S) in SUSYMatching.cpp ");
1573 return (EXIT_FAILURE);
1581 double tan2alpha = 2. *
tanb / (1 -
tanb *
tanb) * (M2A + M2Z) / (M2A - M2Z);
1582 double tana = (-1. +
sqrt(1. + tan2alpha * tan2alpha)) / tan2alpha;
1583 double sina = tana /
sqrt(1. + tana * tana);
1584 double cosa = 1. /
sqrt(1 + tana * tana);
1601 throw std::runtime_error(
"Error in xuS(S) in SUSYMatching.cpp ");
1605 return (EXIT_FAILURE);
1617 return (
XLRS(I, J,
S).conjugate());
1620 throw std::runtime_error(
"Error in XRLS(J,I,S) in SUSYMatching.cpp ");
1623 return (EXIT_FAILURE);
1649 throw std::runtime_error(
"Error in XLRS(J,I,S) in SUSYMatching.cpp ");
1652 return (EXIT_FAILURE);
1671 if (Dmixingflag == 0) {
1672 for (i = 0; i < 6; i++) {
1676 else if (Dmixingflag == 1) {
1677 for (i = 0; i < 3; i++) {
1679 MQuarks(2 * i) +=
mySUSYMQ(2 * i + 1);
1680 MQuarks(2 * i + 1) +=
mySUSYMQ(2 * i);
1685 throw std::runtime_error(
"Error in Dmixingflag in SUSYMatching.cpp. Flag can be either 0 or 1 ");
1696 int I, J, k, l,
S, O;
1708 int D = Dmixingflag;
1710 for (O = 1; O < 9; O++) {
1715 for (I = 0; I < 3; I++) {
1716 M2I = MQuarks(2 * I);
1719 for (J = 0; J < 3; J++) {
1721 M2J = MQuarks(2 * J);
1723 CLO +=
gW *
gW / (32. * M_PI * M_PI) *
1725 PRLk(I, b, 0, D) *
D0N(M2W, M2W, M2I, M2J)
1729 Dk(M2W, M2W, M2I, M2J, 2)
1733 Dk(M2W, M2W, M2I, M2J, 2);
1738 for (I = 1; I < 3; I++) {
1739 M2I = MQuarks(2 * I);
1742 for (J = 1; J < 3; J++) {
1744 M2J = MQuarks(2 * J);
1747 for (k = 0; k < 2; k++) {
1748 for (l = 0; l < 2; l++) {
1750 CLO += -1. / (32. * M_PI * M_PI) *
PLRk(I, q, l, D).
conjugate() *
1752 D0N(M2Hk(k), M2Hk(l), M2I, M2J);
1761 for (I = 0; I < 3; I++) {
1762 M2I = MQuarks(2 * I);
1765 for (J = 0; J < 3; J++) {
1767 M2J = MQuarks(2 * J);
1770 for (k = 0; k < 2; k++) {
1773 CLO +=
gW *
gW / (8. * M_PI * M_PI) *
1776 Dk(M2W, M2Hk(k), M2I, M2J, 2);
1780 for (l = 0; l < 2; l++) {
1782 CLO += -1. / (16. * M_PI * M_PI) *
PLRk(I, q, l, D).
conjugate() *
1784 PLRk(J, b, l, D) *
D0N(M2Hk(k), M2Hk(l), M2I, M2J);
1793 for (
S = 0;
S < 3;
S++) {
1804 for (I = 0; I < 3; I++) {
1805 M2I = MQuarks(2 * I);
1808 for (J = 0; J < 3; J++) {
1810 M2J = MQuarks(2 * J);
1813 for (k = 0; k < 2; k++) {
1814 for (l = 0; l < 2; l++) {
1818 PLRk(J, b, l, D) *
Dk(M2Hk(k), M2Hk(l), M2I, M2J, 2);
1825 for (I = 0; I < 3; I++) {
1826 M2I = MQuarks(2 * I);
1829 for (J = 0; J < 3; J++) {
1831 M2J = MQuarks(2 * J);
1834 for (k = 0; k < 2; k++) {
1835 for (l = 0; l < 2; l++) {
1837 CLO += -1. / (32. * M_PI * M_PI) *
PLRk(I, q, k, D).
conjugate() *
1839 PLRk(J, b, l, D) *
Dk(M2Hk(k), M2Hk(l), M2I, M2J, 2);
1847 for (I = 0; I < 3; I++) {
1848 M2I = MQuarks(2 * I);
1851 for (J = 0; J < 3; J++) {
1853 M2J = MQuarks(2 * J);
1856 for (k = 0; k < 2; k++) {
1857 for (l = 0; l < 2; l++) {
1859 CLO += -1. / (32. * M_PI * M_PI) *
PRLk(I, q, l, D).
conjugate() *
1861 *
D0N(M2Hk(k), M2Hk(l), M2I, M2J);
1869 VCLO.assign(O - 1,CLO);
1888 if (Dmixingflag == 0) {
1892 else if (Dmixingflag == 1) {
1909 for (O = 1; O < 9; O++) {
1914 for (h = 0; h < 6; h++) {
1915 for (k = 0; k < 6; k++) {
1918 myR(h, b) * myR(k, b)
1919 * myR(h, q).conjugate() * myR(k, q).conjugate() *
1921 Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 0) +
1922 11. / 9. *
Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 2));
1927 for (h = 0; h < 6; h++) {
1928 for (k = 0; k < 6; k++) {
1930 CLO += -
Als *
Als * 17. / 18. * M2g * myR(h, b)
1931 * myR(k, b) * myR(h, q + 3).conjugate()
1932 * myR(k, q + 3).conjugate() *
1933 Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 0);
1938 for (h = 0; h < 6; h++) {
1939 for (k = 0; k < 6; k++) {
1941 CLO +=
Als *
Als * 1. / 6. * M2g * myR(h, b)
1942 * myR(k, b) * myR(h, q + 3).conjugate()
1943 * myR(k, q + 3).conjugate() *
1944 Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 0);
1949 for (h = 0; h < 6; h++) {
1950 for (k = 0; k < 6; k++) {
1952 CLO += -
Als *
Als * 7. / 3. * M2g * myR(h, b)
1953 * myR(k, b + 3) * myR(h, q).conjugate()
1954 * myR(k, q + 3).conjugate() *
1955 Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 0) +
1957 Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 2)
1958 * myR(h, b) * myR(k, b + 3) *
1959 (6. * myR(h, q).conjugate() * myR(k, q + 3).conjugate()
1960 + 11. * myR(k, q).conjugate() * myR(h, q + 3).conjugate())
1966 for (h = 0; h < 6; h++) {
1967 for (k = 0; k < 6; k++) {
1969 CLO += -
Als *
Als * 1. / 9. * M2g * myR(h, b)
1970 * myR(k, b + 3) * myR(h, q).conjugate()
1971 * myR(k, q + 3).conjugate() *
1972 Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 0) +
1974 Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 2)
1975 * myR(h, b) * myR(k, b + 3) *
1976 (3. * myR(k, q).conjugate() * myR(h, q + 3).conjugate()
1977 - 2. * myR(h, q).conjugate() * myR(k, q + 3).conjugate())
1983 for (h = 0; h < 6; h++) {
1984 for (k = 0; k < 6; k++) {
1986 CLO += -
Als *
Als * myR(h, b + 3) * myR(k, b + 3)
1987 * myR(h, q + 3).conjugate() * myR(k, q + 3).conjugate() *
1988 (1. / 9. * M2g *
Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 0) -
1989 11. / 9. *
Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 2));
1994 for (h = 0; h < 6; h++) {
1995 for (k = 0; k < 6; k++) {
1997 CLO += -
Als *
Als * 17. / 18. * M2g * myR(h, b + 3)
1998 * myR(k, b + 3) * myR(h, q).conjugate()
1999 * myR(k, q).conjugate() *
2000 Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 0);
2005 for (h = 0; h < 6; h++) {
2006 for (k = 0; k < 6; k++) {
2008 CLO +=
Als *
Als * 1. / 6. * M2g * myR(h, b + 3)
2009 * myR(k, b + 3) * myR(h, q).conjugate()
2010 * myR(k, q).conjugate() *
2011 Dk(myM2Squarks(h), myM2Squarks(k), M2g, M2g, 0);
2017 VCLO.assign(O - 1, CLO);
2039 if (Dmixingflag == 0) {
2042 else if (Dmixingflag == 1) {
2046 int D = Dmixingflag;
2048 for (O = 1; O < 9; O++) {
2053 for (i = 0; i < 2; i++) {
2054 for (j = 0; j < 2; j++) {
2055 for (h = 0; h < 6; h++) {
2056 for (k = 0; k < 6; k++) {
2058 CLO += -1. / (32. * M_PI * M_PI) *
Dk(myM2Squarks(k),
2069 for (i = 0; i < 2; i++) {
2070 for (j = 0; j < 2; j++) {
2071 for (h = 0; h < 6; h++) {
2072 for (k = 0; k < 6; k++) {
2074 CLO += -1. / (32. * M_PI * M_PI) *
2087 for (i = 0; i < 2; i++) {
2088 for (j = 0; j < 2; j++) {
2089 for (h = 0; h < 6; h++) {
2090 for (k = 0; k < 6; k++) {
2094 VdUCR(b, h, j, 1, D) *
Dk(myM2Squarks(k),
2103 for (i = 0; i < 2; i++) {
2104 for (j = 0; j < 2; j++) {
2105 for (h = 0; h < 6; h++) {
2106 for (k = 0; k < 6; k++) {
2108 CLO += -1. / (16. * M_PI * M_PI) *
VdUCL(q, k, j, D).
conjugate() *
2110 VdUCR(b, h, j, 1, D) *
D0N(myM2Squarks(k),
2119 for (i = 0; i < 2; i++) {
2120 for (j = 0; j < 2; j++) {
2121 for (h = 0; h < 6; h++) {
2122 for (k = 0; k < 6; k++) {
2124 CLO += -1. / (32. * M_PI * M_PI) *
Dk(myM2Squarks(k),
2128 VdUCR(b, k, j, 1, D);
2135 for (i = 0; i < 2; i++) {
2136 for (j = 0; j < 2; j++) {
2137 for (h = 0; h < 6; h++) {
2138 for (k = 0; k < 6; k++) {
2140 CLO += -1. / (32. * M_PI * M_PI) *
VdUCL(q, k, j, D).
conjugate() *
2142 VdUCR(b, h, j, 1, D) *
2143 D0N(myM2Squarks(k), myM2Squarks(h),
2151 VCLO.assign(O - 1, CLO);
2171 if (Dmixingflag == 0) {
2174 else if (Dmixingflag == 1) {
2179 int D = Dmixingflag;
2181 for (O = 1; O < 9; O++) {
2186 for (i = 0; i < 4; i++) {
2187 for (j = 0; j < 4; j++) {
2188 for (h = 0; h < 6; h++) {
2189 for (k = 0; k < 6; k++) {
2192 (1./ 32. / M_PI / M_PI *
VdDNL(b, h, j, 1, D)
2194 *
Dk(myM2Squarks(k),myM2Squarks(h),
MChi0(i) *
MChi0(i),
2196 / 64. / M_PI / M_PI *
Dk(myM2Squarks(k),
2207 for (i = 0; i < 4; i++) {
2208 for (j = 0; j < 4; j++) {
2209 for (h = 0; h < 6; h++) {
2210 for (k = 0; k < 6; k++) {
2212 CLO +=
MChi0(i) *
MChi0(j) / 32. / M_PI / M_PI *
2213 Dk(myM2Squarks(k), myM2Squarks(h),
2224 for (i = 0; i < 4; i++) {
2225 for (j = 0; j < 4; j++) {
2226 for (h = 0; h < 6; h++) {
2227 for (k = 0; k < 6; k++) {
2229 CLO += -
MChi0(i) *
MChi0(j) / 32. / M_PI / M_PI *
2230 Dk(myM2Squarks(k), myM2Squarks(h),
MChi0(i) *
MChi0(i),
2242 for (i = 0; i < 4; i++) {
2243 for (j = 0; j < 4; j++) {
2244 for (h = 0; h < 6; h++) {
2245 for (k = 0; k < 6; k++) {
2247 CLO += 1. / 8. / M_PI / M_PI *
Dk(myM2Squarks(k),
2260 for (i = 0; i < 4; i++) {
2261 for (j = 0; j < 4; j++) {
2262 for (h = 0; h < 6; h++) {
2263 for (k = 0; k < 6; k++) {
2265 CLO += -1. / 8. / M_PI / M_PI *
Dk(myM2Squarks(k),
2270 M_PI *
Dk(myM2Squarks(k), myM2Squarks(h),
2281 for (i = 0; i < 4; i++) {
2282 for (j = 0; j < 4; j++) {
2283 for (h = 0; h < 6; h++) {
2284 for (k = 0; k < 6; k++) {
2287 (1./ 32. / M_PI / M_PI *
VdDNR(b, h, j, 1, D)
2289 *
Dk(myM2Squarks(k),myM2Squarks(h),
MChi0(i) *
MChi0(i),
2291 / 64. / M_PI / M_PI *
Dk(myM2Squarks(k),
2302 for (i = 0; i < 4; i++) {
2303 for (j = 0; j < 4; j++) {
2304 for (h = 0; h < 6; h++) {
2305 for (k = 0; k < 6; k++) {
2307 CLO +=
MChi0(i) *
MChi0(j) / 32. / M_PI / M_PI *
2308 Dk(myM2Squarks(k), myM2Squarks(h),
2320 for (i = 0; i < 4; i++) {
2321 for (j = 0; j < 4; j++) {
2322 for (h = 0; h < 6; h++) {
2323 for (k = 0; k < 6; k++) {
2325 CLO += -
MChi0(i) *
MChi0(j) / 32. / M_PI / M_PI *
2326 Dk(myM2Squarks(k), myM2Squarks(h),
MChi0(i) *
MChi0(i),
2337 VCLO.assign(O - 1, CLO);
2362 if (Dmixingflag == 0) {
2366 else if (Dmixingflag == 1) {
2374 int D = Dmixingflag;
2377 for (O = 1; O < 9; O++) {
2382 for (i = 0; i < 4; i++) {
2383 for (h = 0; h < 6; h++) {
2384 for (k = 0; k < 6; k++) {
2386 CLO += -
Als * 2. / 3. / 4. / M_PI *
Dk(myM2Squarks(k), myM2Squarks(h),
2387 MChi0(i) *
MChi0(i), M2g, 2) * myR(k, q).conjugate()
2391 , M2g, 0) * (myR(h, q).conjugate() * myR(k, q).conjugate()
2392 *
VdDNL(b, h, i, 1, D) *
VdDNL(b, k, i, 1, D)
2401 for (i = 0; i < 4; i++) {
2402 for (h = 0; h < 6; h++) {
2403 for (k = 0; k < 6; k++) {
2406 Dk(myM2Squarks(k), myM2Squarks(h),
MChi0(i)
2407 *
MChi0(i), M2g, 0) * (3. * myR(h, b) *
2408 myR(k, q + 3).conjugate() *
VdDNL(b, k, i, 1, D) *
2410 + myR(k, b) * myR(h, b) *
2412 myR(k, q + 3).conjugate() * myR(h, q + 3).conjugate() *
2413 VdDNL(b, k, i, 1, D) *
VdDNL(b, h, i, 1, D));
2420 for (i = 0; i < 4; i++) {
2421 for (h = 0; h < 6; h++) {
2422 for (k = 0; k < 6; k++) {
2424 CLO += -
Als / 4. / M_PI *
MChi0(i) *
Mg / 3. *
2425 Dk(myM2Squarks(k), myM2Squarks(h),
MChi0(i)
2426 *
MChi0(i), M2g, 0) * (
2427 myR(h, b) * myR(k, q + 3).conjugate() *
VdDNL(b, k, i, 1, D) *
2431 myR(h, q + 3).conjugate() *
VdDNL(b, k, i, 1, D) *
2432 VdDNL(b, h, i, 1, D));
2439 for (i = 0; i < 4; i++) {
2440 for (h = 0; h < 6; h++) {
2441 for (k = 0; k < 6; k++) {
2443 CLO += -
Als / 4. / M_PI * 2. / 3. *
2444 Dk(myM2Squarks(k), myM2Squarks(h),
MChi0(i)
2445 *
MChi0(i), M2g, 2) * (
2446 myR(h, b) * myR(k, q).conjugate() *
VdDNR(b, k, i, 1, D) *
2449 myR(k, q + 3).conjugate() *
VdDNL(b, k, i, 1, D) *
2451 - myR(h, b) * myR(k, b + 3) *
2453 - myR(h, q).conjugate() * myR(k, q + 3).conjugate() *
2455 - 3. * myR(k, b + 3) * myR(h, b) *
2457 - 3. * myR(k, q + 3).conjugate() *
2458 myR(h, q).conjugate() *
VdDNL(b, h, i, 1, D) *
VdDNR(b, k, i, 1, D)) +
2460 Dk(myM2Squarks(k), myM2Squarks(h),
MChi0(i)*
MChi0(i), M2g, 0) *
2461 (myR(h, b) * myR(k, q + 3).conjugate() *
2463 myR(h, b + 3) * myR(k, q).conjugate() *
2471 for (i = 0; i < 4; i++) {
2472 for (h = 0; h < 6; h++) {
2473 for (k = 0; k < 6; k++) {
2475 CLO +=
Als / 4. / M_PI * 2. / 3. *
2476 Dk(myM2Squarks(k), myM2Squarks(h),
MChi0(i)
2477 *
MChi0(i), M2g, 2)*(
2478 3. * myR(h, b) * myR(k, q).conjugate() *
VdDNR(b, k, i, 1, D) *
2480 myR(k, q + 3).conjugate() *
VdDNL(b, k, i, 1, D) *
2482 - 3. * myR(h, b) * myR(k, b + 3) *
2484 - 3. * myR(h, q).conjugate() * myR(k, q + 3).conjugate() *
2485 VdDNL(b, k, i, 1, D) *
VdDNR(b, h, i, 1, D) - myR(k, b + 3) *
2488 myR(h, q).conjugate() *
VdDNL(b, h, i, 1, D) *
VdDNR(b, k, i, 1, D)) -
2490 Dk(myM2Squarks(k), myM2Squarks(h),
MChi0(i)
2491 *
MChi0(i), M2g, 0)* (myR(h, b) * myR(k, q + 3).conjugate() *
2493 myR(h, b + 3) * myR(k, q).conjugate() *
2501 for (i = 0; i < 4; i++) {
2502 for (h = 0; h < 6; h++) {
2503 for (k = 0; k < 6; k++) {
2505 CLO += -
Als * 2. / 3. / 4. / M_PI *
Dk(myM2Squarks(k), myM2Squarks(h),
2506 MChi0(i) *
MChi0(i), M2g, 2) * myR(k, q + 3).conjugate()
2510 , M2g, 0) * (myR(h, q + 3).conjugate() * myR(k, q + 3).conjugate()
2511 *
VdDNR(b, h, i, 1, D) *
VdDNR(b, k, i, 1, D)
2512 + myR(h, b + 3) * myR(k, b + 3) *
VdDNR(q, h, i, 1, D).
conjugate() *
2520 for (i = 0; i < 4; i++) {
2521 for (h = 0; h < 6; h++) {
2522 for (k = 0; k < 6; k++) {
2525 Dk(myM2Squarks(k), myM2Squarks(h),
MChi0(i)
2526 *
MChi0(i), M2g, 0) * (3. * myR(h, b + 3) *
2527 myR(k, q).conjugate() *
VdDNR(b, k, i, 1, D) *
2529 + myR(k, b + 3) * myR(h, b + 3) *
2531 myR(k, q).conjugate() * myR(h, q).conjugate() *
2532 VdDNR(b, k, i, 1, D) *
VdDNR(b, h, i, 1, D));
2539 for (i = 0; i < 4; i++) {
2540 for (h = 0; h < 6; h++) {
2541 for (k = 0; k < 6; k++) {
2543 CLO += -
Als / 4. / M_PI *
MChi0(i) *
Mg / 3. *
2544 Dk(myM2Squarks(k), myM2Squarks(h),
MChi0(i)
2545 *
MChi0(i), M2g, 0) * (
2546 myR(h, b + 3) * myR(k, q).conjugate() *
VdDNR(b, k, i, 1, D) *
2550 myR(h, q).conjugate() *
VdDNR(b, k, i, 1, D) *
2551 VdDNR(b, h, i, 1, D));
2558 VCLO.assign(O - 1, CLO);
2594 for (i = 0; i < 5; i++) {
2604 throw std::runtime_error(
"Error in SUSYMatching::CMdbd2()");
2612 for (i = 0; i < 5; i++) {
2622 throw std::runtime_error(
"Error in SUSYMatching::CMdbd2()");
2629 for (i = 0; i < 5; i++) {
2639 throw std::runtime_error(
"Error in SUSYMatching::CMdbd2()");
2646 for (i = 0; i < 5; i++) {
2656 throw std::runtime_error(
"Error in SUSYMatching::CMdbd2()");
2663 for (i = 0; i < 5; i++) {
2673 throw std::runtime_error(
"Error in SUSYMatching::CMdbd2()");
2683 for (i = 0; i < 3; i++) {
2693 throw std::runtime_error(
"Error in SUSYMatching::CMdbd2()");
2701 for(i = 0; i < 3 ; i++){
2711 throw std::runtime_error(
"Error in SUSYMatching::CMdbd2()");
2718 for (i = 0; i < 3; i++) {
2728 throw std::runtime_error(
"Error in SUSYMatching::CMdbd2()");
2736 for (i = 0; i < 3; i++) {
2746 throw std::runtime_error(
"Error in SUSYMatching::CMdbd2()");
2753 for (i = 0; i < 3; i++) {
2763 throw std::runtime_error(
"Error in SUSYMatching::CMdbd2()");
2795 for (i = 0; i < 5; i++) {
2805 throw std::runtime_error(
"Error in SUSYMatching::CMdbs2()");
2812 for(i = 0; i < 5 ; i++){
2822 throw std::runtime_error(
"Error in SUSYMatching::CMdbs2()");
2829 for (i = 0; i < 5; i++) {
2839 throw std::runtime_error(
"Error in SUSYMatching::CMdbs2()");
2846 for (i = 0; i < 5; i++) {
2856 throw std::runtime_error(
"Error in SUSYMatching::CMdbs2()");
2863 for (i = 0; i < 5; i++) {
2873 throw std::runtime_error(
"Error in SUSYMatching::CMdbs2()");
2883 for (i = 0; i < 3; i++) {
2893 throw std::runtime_error(
"Error in SUSYMatching::CMdbs2()");
2900 for(i = 0; i < 3 ; i++){
2910 throw std::runtime_error(
"Error in SUSYMatching::CMdbs2()");
2917 for (i = 0; i < 3; i++) {
2927 throw std::runtime_error(
"Error in SUSYMatching::CMdbs2()");
2934 for (i = 0; i < 3; i++) {
2944 throw std::runtime_error(
"Error in SUSYMatching::CMdbs2()");
2951 for (i = 0; i < 3; i++) {
2961 throw std::runtime_error(
"Error in SUSYMatching::CMdbs2()");
2985 vmdk2 = StandardModelMatching::CMdk2();
2994 for (i = 0; i < 5; i++) {
3004 throw std::runtime_error(
"Error in SUSYMatching::CMdk2()");
3011 for(i = 0; i < 5 ; i++){
3021 throw std::runtime_error(
"Error in SUSYMatching::CMdk2()");
3028 for (i = 0; i < 5; i++) {
3038 throw std::runtime_error(
"Error in SUSYMatching::CMdk2()");
3045 for (i = 0; i < 5; i++) {
3055 throw std::runtime_error(
"Error in SUSYMatching::CMdk2()");
3062 for (i = 0; i < 5; i++) {
3072 throw std::runtime_error(
"Error in SUSYMatching::CMdk2()");
3082 for (i = 0; i < 3; i++) {
3092 throw std::runtime_error(
"Error in SUSYMatching::CMdk2()");
3099 for(i = 0; i < 3 ; i++){
3109 throw std::runtime_error(
"Error in SUSYMatching::CMdk2()");
3116 for (i = 0; i < 3; i++) {
3126 throw std::runtime_error(
"Error in SUSYMatching::CMdk2()");
3133 for (i = 0; i < 3; i++) {
3143 throw std::runtime_error(
"Error in SUSYMatching::CMdk2()");
3150 for (i = 0; i < 3; i++) {
3160 throw std::runtime_error(
"Error in SUSYMatching::CMdk2()");
3194 for (i = 0; i < 5; i++) {
3204 throw std::runtime_error(
"Error in SUSYMatching::CMdd2()");
3211 for(i = 0; i < 5 ; i++){
3221 throw std::runtime_error(
"Error in SUSYMatching::CMdd2()");
3228 for (i = 0; i < 5; i++) {
3238 throw std::runtime_error(
"Error in SUSYMatching::CMdd2()");
3245 for (i = 0; i < 5; i++) {
3255 throw std::runtime_error(
"Error in SUSYMatching::CMdd2()");
3262 for (i = 0; i < 5; i++) {
3272 throw std::runtime_error(
"Error in SUSYMatching::CMdd2()");
3282 for (i = 0; i < 3; i++) {
3292 throw std::runtime_error(
"Error in SUSYMatching::CMdd2()");
3299 for(i = 0; i < 3 ; i++){
3309 throw std::runtime_error(
"Error in SUSYMatching::CMdd2()");
3316 for (i = 0; i < 3; i++) {
3326 throw std::runtime_error(
"Error in SUSYMatching::CMdd2()");
3333 for (i = 0; i < 3; i++) {
3343 throw std::runtime_error(
"Error in SUSYMatching::CMdd2()");
3350 for (i = 0; i < 3; i++) {
3360 throw std::runtime_error(
"Error in SUSYMatching::CMdd2()");
3390 if (fabs(x - 1.) < SUSYLEPS)
return (-5. / 48.);
3392 return ((x * (7 - 5 * x - 8 * x * x)) / (24. * (-1 + x)*(-1 + x)*(-1 + x)) +
3393 x * x * (-2 + 3 * x) *
log(x) / (4. * (-1 + x)*(-1 + x)*(-1 + x)*(-1 + x)));
3398 if (fabs(x - 1.) < SUSYLEPS)
return (-7. / 36.);
3400 return (((3 - 5 * x) * x) / (12. * (-1 + x)*(-1 + x)) +
3401 (x * (-2 + 3 * x) *
log(x)) / (6. * (-1 + x)*(-1 + x)*(-1 + x)));
3404 else throw std::runtime_error(
"Error in F7k ");
3422 *
tanb) *
F7k(m2top / M2W, 2));
3424 VCLO.assign(1, 1. / (3. *
tanb *
tanb) *
F7k(m2top / M2Hp, 1) +
3428 for (j = 0; j < 2; j++) {
3429 for (k = 0; k < 6; k++) {
3430 CLO += 1. / (3. *
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2)) *
3440 VCLO.assign(2, CLO);
3456 std::stringstream out;
3458 throw std::runtime_error(
"SUSYMatching::CMbsg(): scheme " + out.str() +
"not implemented");
3472 std::stringstream out;
3474 throw std::runtime_error(
"SUSYMatching::CMbsg(): order " + out.str() +
"not implemented");
3500 std::runtime_error(
"SUSYMatching::CMBMll(): LEPTON not implemented");
3511 std::stringstream out;
3513 throw std::runtime_error(
"StandardModel::CMBKstrall(): scheme " + out.str() +
"not implemented");
3526 std::stringstream out;
3528 throw std::runtime_error(
"StandardModelMatching::CMBKstrall(): order " + out.str() +
"not implemented");
3531 vmcBMll.push_back(
mcBMll);
3541 if (std::fabs(x - 1.) > SUSYLEPS2) {
3542 return (-7. + 5. * x + 8. * x * x) / (6. * (1. - x) * (1. - x) * (1. - x))
3543 - (2. * x - 3. * x * x) / ((1. - x) * (1. - x) * (1. - x) * (1. - x)) *
log(x);
3552 if (std::fabs(x - 1.) > SUSYLEPS2) {
3553 return (3. * x - 5. * x * x) / (2. * (1. - x) * (1. - x))
3554 + (2. * x - 3. * x * x) / ((1. - x) * (1. - x) * (1. - x)) *
log(x);
3563 if (std::fabs(x - 1.) > SUSYLEPS2) {
3564 return (2. + 5. * x - x * x) / (6. * (1. - x) * (1. - x) * (1. - x))
3565 + x / ((1. - x) * (1. - x) * (1. - x) * (1. - x)) *
log(x);
3574 if (std::fabs(x - 1.) > SUSYLEPS2) {
3575 return (1. + x) / (2. * (1. - x) * (1. - x))
3576 + x / ((1. - x) * (1. - x) * (1. - x)) *
log(x);
3585 if (std::fabs(x - 1.) > SUSYLEPS2) {
3586 return x / (1. - x) + x / ((1. - x) * (1. - x)) *
log(x);
3587 }
else return -1. / 2.;
3593 if (std::fabs(x - 1.) > SUSYLEPS2) {
3594 return (38. * x - 79. * x * x + 47. * x * x * x) / (6. * (1. - x) * (1. - x) * (1. - x))
3595 + (4. * x - 6. * x * x + 3. * x * x * x * x) / ((1. - x) * (1. - x) * (1. - x) * (1. - x)) *
log(x);
3604 if (std::fabs(x - 1.) > SUSYLEPS2) {
3605 return (52. - 101. * x + 43. * x * x) / (6. * (1. - x) * (1. - x) * (1. - x))
3606 + (6. - 9. * x + 2. * x * x * x) / ((1. - x) * (1. - x) * (1. - x) * (1. - x)) *
log(x);
3615 if (std::fabs(x - 1.) > SUSYLEPS2) {
3616 return (2. - 7. * x + 11. * x * x) / ((1. - x) * (1. - x) * (1. - x))
3617 + (6. * x * x * x) / ((1. - x) * (1. - x) * (1. - x) * (1. - x)) *
log(x);
3626 if (m1 <= 0.0 || m2 <= 0.0 || m3 <= 0.0)
throw std::runtime_error(
"SUSYMatching::bsll_C0(): Invalid argument!");
3628 if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m3) / m2) <= SUSYLEPS2)) {
3629 return -1. / (2 * m1 * m1);
3630 }
else if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m3) / m2) > SUSYLEPS2)) {
3631 return (m3 * m3 * (1. -
log(m3 * m3 / (m2 * m2))) - m2 * m2) / ((m2 * m2 - m3 * m3)*(m2 * m2 - m3 * m3));
3632 }
else if ((std::fabs((m1 - m3) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m3) / m2) > SUSYLEPS2)) {
3633 return (m2 * m2 * (1. -
log(m2 * m2 / (m3 * m3))) - m3 * m3) / ((m2 * m2 - m3 * m3)*(m2 * m2 - m3 * m3));
3634 }
else if ((std::fabs((m1 - m2) / m1) > SUSYLEPS2) && (std::fabs((m2 - m3) / m2) <= SUSYLEPS2)) {
3635 return (m1 * m1 * (1. +
log(m3 * m3 / (m1 * m1))) - m3 * m3) / ((m1 * m1 - m3 * m3)*(m1 * m1 - m3 * m3));
3636 }
else if ((std::fabs((m1 - m2) / m1) > SUSYLEPS2) && (std::fabs((m2 - m3) / m2) > SUSYLEPS2)) {
3637 return (m2 * m2 / (m1 * m1 - m2 * m2) *
log(m2 * m2 / (m1 * m1)) - m3 * m3 / (m1 * m1 - m3 * m3) *
log(m3 * m3 / (m1 * m1))) / (m2 * m2 - m3 * m3);
3639 std::cout <<
"m1, m2, m3 = " << m1 <<
", " << m2 <<
", " << m3 << std::endl;
3640 throw std::runtime_error(
"SUSYMatching::bsll_C0(): limit/case not implemented!");
3647 if (mu2 <= 0.0 || m1 <= 0.0 || m2 <= 0.0 || m3 <= 0.0)
throw std::runtime_error(
"SUSYMatching::bsll_C2(): Invalid argument!");
3649 if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m3) / m2) <= SUSYLEPS2)) {
3650 return -(1 / 4.) *
log(m3 * m3 / (mu2));
3651 }
else if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m3) / m2) > SUSYLEPS2)) {
3652 return (
pow(m2, 4.) - 4. * m2 * m2 * m3 * m3 + 3. *
pow(m3, 4.) - (2. *
pow(m2, 4.) - 4. * m2 * m2 * m3 * m3) *
log(m2 * m2 / (mu2)) - 2. *
pow(m3, 4.) *
log(m3 * m3 / (mu2))) / (8. * (m2 * m2 - m3 * m3)*(m2 * m2 - m3 * m3));
3653 }
else if ((std::fabs((m1 - m3) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m3) / m2) > SUSYLEPS2)) {
3654 return (3. *
pow(m2, 4.) - 4. * m2 * m2 * m3 * m3 +
pow(m3, 4.) + (-2. *
pow(m3, 4.) + 4. * m2 * m2 * m3 * m3) *
log(m3 * m3 / (mu2)) - 2. *
pow(m2, 4.) *
log(m2 * m2 / (mu2))) / (8. * (m2 * m2 - m3 * m3)*(m2 * m2 - m3 * m3));
3655 }
else if ((std::fabs((m1 - m2) / m1) > SUSYLEPS2) && (std::fabs((m2 - m3) / m2) <= SUSYLEPS2)) {
3656 return (3. *
pow(m1, 4.) - 4. * m1 * m1 * m3 * m3 +
pow(m3, 4.) + (-2. *
pow(m3, 4.) + 4. * m1 * m1 * m3 * m3) *
log(m3 * m3 / (mu2)) - 2. *
pow(m1, 4.) *
log(m1 * m1 / (mu2))) / (8. * (m1 * m1 - m3 * m3)*(m1 * m1 - m3 * m3));
3657 }
else if ((std::fabs((m1 - m2) / m1) > SUSYLEPS2) && (std::fabs((m2 - m3) / m2) > SUSYLEPS2)) {
3658 return 3. / 8. - (1. / 4.)*(
pow(m1, 4.) / ((m1 * m1 - m2 * m2)*(m1 * m1 - m3 * m3)) *
log(m1 * m1 / (mu2)) +
pow(m2, 4.) / ((m2 * m2 - m1 * m1)*(m2 * m2 - m3 * m3)) *
log(m2 * m2 / (mu2)) +
pow(m3, 4.) / ((m3 * m3 - m2 * m2)*(m3 * m3 - m1 * m1)) *
log(m3 * m3 / (mu2)));
3660 std::cout <<
"m1, m2, m3, mu2 = " << m1 <<
", " << m2 <<
", " << m3 <<
", " << mu2 << std::endl;
3661 throw std::runtime_error(
"SUSYMatching::bsll_C2(): limit/case not implemented!");
3668 if (m1 <= 0.0 || m2 <= 0.0 || m3 <= 0.0 || m4 <= 0.0)
throw std::runtime_error(
"SUSYMatching::bsll_D0(): Invalid argument!");
3670 if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2) && (std::fabs((m1 - m3) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m4) / m2) <= SUSYLEPS2)) {
3671 return 1. / (6. * m1 * m1 * m1 * m1);
3672 }
else if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m3) / m3) <= SUSYLEPS2)) {
3673 return (m3 * m3 * m3 * m3 - m4 * m4 * m4 * m4 + 4. * m3 * m3 * m4 * m4 *
log(m4 / m3)) / (2. * (m3 * m3 - m4 * m4)*(m3 * m3 - m4 * m4)*(m3 * m3 - m4 * m4) * m3 * m3);
3674 }
else if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m4) / m4) <= SUSYLEPS2)) {
3675 return (-m3 * m3 * m3 * m3 + m4 * m4 * m4 * m4 + 4. * m3 * m3 * m4 * m4 *
log(m3 / m4)) / (2. * (-m3 * m3 + m4 * m4)*(-m3 * m3 + m4 * m4)*(-m3 * m3 + m4 * m4) * m4 * m4);
3676 }
else if ((std::fabs((m1 - m3) / m1) <= SUSYLEPS2) && (std::fabs((m3 - m4) / m4) <= SUSYLEPS2)) {
3677 return (-m2 * m2 * m2 * m2 + m4 * m4 * m4 * m4 + 4. * m2 * m2 * m4 * m4 *
log(m2 / m4)) / (2. * (-m2 * m2 + m4 * m4)*(-m2 * m2 + m4 * m4)*(-m2 * m2 + m4 * m4) * m4 * m4);
3678 }
else if ((std::fabs((m2 - m3) / m2) <= SUSYLEPS2) && (std::fabs((m3 - m4) / m4) <= SUSYLEPS2)) {
3679 return (-m1 * m1 * m1 * m1 + m4 * m4 * m4 * m4 + 4. * m1 * m1 * m4 * m4 *
log(m1 / m4)) / (2. * (-m1 * m1 + m4 * m4)*(-m1 * m1 + m4 * m4)*(-m1 * m1 + m4 * m4) * m4 * m4);
3680 }
else if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2)) {
3681 double term1 = (m3 * m3 * (-m2 * m2 + m3 * m3 + 2. * m2 * m2 *
log(m2 / m3))) / (m2 * m2 - m3 * m3) / (m2 * m2 - m3 * m3);
3682 double term2 = (m4 * m4 * (-m2 * m2 + m4 * m4 + 2. * m2 * m2 *
log(m2 / m4))) / (m2 * m2 - m4 * m4) / (m2 * m2 - m4 * m4);
3683 return (term1 - term2) / (m2 * m2 * (m3 * m3 - m4 * m4));
3684 }
else if ((std::fabs((m1 - m3) / m1) <= SUSYLEPS2)) {
3685 double term1 = (2. * m3 * m3 *
log(m3 / m2)) / (-m2 * m2 + m3 * m3);
3686 double term2 = (2. * m4 * m4 *
log(m4 / m2)) / (-m2 * m2 + m4 * m4);
3687 double term3 = (2. * m4 * m4 *
log(m4 / m3)) / (-m3 * m3 + m4 * m4);
3688 return (-1. + term1 - term2 + term3) / (m3 * m3 - m2 * m2) / (m3 * m3 - m4 * m4);
3689 }
else if ((std::fabs((m1 - m4) / m1) <= SUSYLEPS2)) {
3690 double term1 = (2. * m3 * m3 *
log(m3 / m2)) / (-m2 * m2 + m3 * m3);
3691 double term2 = (2. * m3 * m3 *
log(m3 / m4)) / (-m4 * m4 + m3 * m3);
3692 double term3 = (2. * m4 * m4 *
log(m4 / m2)) / (-m2 * m2 + m4 * m4);
3693 return (1. + term1 - term2 - term3) / (m4 * m4 - m2 * m2) / (m3 * m3 - m4 * m4);
3694 }
else if ((std::fabs((m2 - m3) / m2) <= SUSYLEPS2)) {
3695 double term1 = (2. * m3 * m3 *
log(m3 / m1)) / (-m1 * m1 + m3 * m3);
3696 double term2 = (2. * m4 * m4 *
log(m4 / m1)) / (-m1 * m1 + m4 * m4);
3697 double term3 = (2. * m4 * m4 *
log(m4 / m3)) / (-m3 * m3 + m4 * m4);
3698 return (1. - term1 + term2 - term3) / (-m4 * m4 + m3 * m3) / (-m3 * m3 + m1 * m1);
3699 }
else if ((std::fabs((m2 - m4) / m2) <= SUSYLEPS2)) {
3700 double term1 = (2. * m3 * m3 *
log(m3 / m1)) / (-m1 * m1 + m3 * m3);
3701 double term2 = (2. * m3 * m3 *
log(m3 / m4)) / (-m4 * m4 + m3 * m3);
3702 double term3 = (2. * m4 * m4 *
log(m4 / m1)) / (-m1 * m1 + m4 * m4);
3703 return (-1. - term1 + term2 + term3) / (-m4 * m4 + m3 * m3) / (-m4 * m4 + m1 * m1);
3704 }
else if ((std::fabs((m3 - m4) / m3) <= SUSYLEPS2)) {
3705 double term1 = (-m4 * m4 + m1 * m1 + 2. * m1 * m1 *
log(m4 / m1)) / (m1 * m1 - m4 * m4) / (m1 * m1 - m4 * m4);
3706 double term2 = (-m4 * m4 + m2 * m2 + 2. * m2 * m2 *
log(m4 / m2)) / (m2 * m2 - m4 * m4) / (m2 * m2 - m4 * m4);
3707 return (term1 - term2) / (m1 * m1 - m2 * m2);
3708 }
else if ((std::fabs((m1 - m2) / m1) > SUSYLEPS2) && (std::fabs((m3 - m4) / m3) > SUSYLEPS2) && (std::fabs((m1 - m4) / m1) > SUSYLEPS2) && (std::fabs((m1 - m3) / m1) > SUSYLEPS2) && (std::fabs((m2 - m3) / m2) > SUSYLEPS2) && (std::fabs((m2 - m4) / m2) > SUSYLEPS2)) {
3709 return (
bsll_C0(m1, m3, m4) -
bsll_C0(m2, m3, m4)) / (m1 * m1 - m2 * m2);
3711 std::cout <<
"m1, m2, m3, m4 = " << m1 <<
", " << m2 <<
", " << m3 <<
", " << m4 << std::endl;
3712 throw std::runtime_error(
"SUSYMatching::bsll_D0reg(): limit/case not implemented!");
3719 if (m1 <= 0.0 || m2 <= 0.0 || m3 <= 0.0 || m4 <= 0.0)
throw std::runtime_error(
"SUSYMatching::bsll_D2(): Invalid argument!");
3721 if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2) && (std::fabs((m1 - m3) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m4) / m2) <= SUSYLEPS2)) {
3722 return -1. / (12. * m1 * m1 * m1 * m1);
3723 }
else if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m3) / m3) <= SUSYLEPS2)) {
3724 return -(m3 * m3 * m3 * m3 - 4. * m3 * m3 * m4 * m4 + 3. * m4 * m4 * m4 * m4 - 4. * m4 * m4 * m4 * m4 *
log(m4 / m3)) / (8. * (m3 * m3 - m4 * m4)*(m3 * m3 - m4 * m4)*(m3 * m3 - m4 * m4));
3725 }
else if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2) && (std::fabs((m2 - m4) / m4) <= SUSYLEPS2)) {
3726 return (m4 * m4 * m4 * m4 - 4. * m3 * m3 * m4 * m4 + 3. * m3 * m3 * m3 * m3 - 4. * m3 * m3 * m3 * m3 *
log(m3 / m4)) / (8. * (m3 * m3 - m4 * m4)*(m3 * m3 - m4 * m4)*(m3 * m3 - m4 * m4));
3727 }
else if ((std::fabs((m1 - m3) / m1) <= SUSYLEPS2) && (std::fabs((m3 - m4) / m3) <= SUSYLEPS2)) {
3728 return (m4 * m4 * m4 * m4 - 4. * m2 * m2 * m4 * m4 + 3. * m2 * m2 * m2 * m2 - 4. * m2 * m2 * m2 * m2 *
log(m2 / m4)) / (8. * (m2 * m2 - m4 * m4)*(m2 * m2 - m4 * m4)*(m2 * m2 - m4 * m4));
3729 }
else if ((std::fabs((m2 - m3) / m2) <= SUSYLEPS2) && (std::fabs((m3 - m4) / m3) <= SUSYLEPS2)) {
3730 return (m4 * m4 * m4 * m4 - 4. * m1 * m1 * m4 * m4 + 3. * m1 * m1 * m1 * m1 - 4. * m1 * m1 * m1 * m1 *
log(m1 / m4)) / (8. * (m1 * m1 - m4 * m4)*(m1 * m1 - m4 * m4)*(m1 * m1 - m4 * m4));
3731 }
else if ((std::fabs((m1 - m2) / m1) <= SUSYLEPS2)) {
3732 double num = 2. * m3 * m3 * m3 * m3 * (m2 * m2 - m4 * m4)*(m2 * m2 - m4 * m4) *
log(m3 / m2)+(m2 * m2 - m3 * m3)*(m2 * m2 * (m2 * m2 - m4 * m4)*(m3 * m3 - m4 * m4) + 2. * (-m2 * m2 + m3 * m3) * m4 * m4 * m4 * m4 *
log(m4 / m2));
3733 double den = 4. * (m2 * m2 - m3 * m3)*(m2 * m2 - m3 * m3)*(m2 * m2 - m4 * m4)*(m2 * m2 - m4 * m4)*(m3 * m3 - m4 * m4);
3734 return -(num / den);
3735 }
else if ((std::fabs((m1 - m3) / m1) <= SUSYLEPS2)) {
3736 double num = -m2 * m2 * (2. * m3 * m3 *
log(m3 / m2) / (m2 * m2 - m3 * m3) + 2. * m4 * m4 *
log(m4 / m2) / (-m2 * m2 + m4 * m4)) + m3 * m3 * (-1. + 2. * m4 * m4 *
log(m4 / m3) / (m4 * m4 - m3 * m3));
3737 double den = 4. * (-m2 * m2 + m3 * m3)*(m3 * m3 - m4 * m4);
3739 }
else if ((std::fabs((m1 - m4) / m1) <= SUSYLEPS2)) {
3740 double num = 2. * m2 * m2 * m3 * m3 *
log(m3 / m2) / (-m2 * m2 + m3 * m3) + m4 * m4 * (1. - 2. * m3 * m3 *
log(m3 / m4) / (m3 * m3 - m4 * m4) - 2. * m2 * m2 *
log(m2 / m4) / (m2 * m2 - m4 * m4));
3741 double den = 4. * (-m2 * m2 + m4 * m4)*(m3 * m3 - m4 * m4);
3743 }
else if ((std::fabs((m2 - m3) / m2) <= SUSYLEPS2)) {
3744 double num = m3 * m3 + m1 * m1 * (2. * m3 * m3 *
log(m3 / m1) / (m1 * m1 - m3 * m3) + 2. * m4 * m4 *
log(m4 / m1) / (m4 * m4 - m1 * m1)) + 2. * m3 * m3 * m4 * m4 *
log(m4 / m3) / (m3 * m3 - m4 * m4);
3745 double den = 4. * (m1 * m1 - m3 * m3)*(m3 * m3 - m4 * m4);
3747 }
else if ((std::fabs((m2 - m4) / m2) <= SUSYLEPS2)) {
3748 double num = 2. * m1 * m1 * m3 * m3 *
log(m3 / m1) / (m1 * m1 - m3 * m3) + m4 * m4 * (-1. + 2. * m3 * m3 *
log(m3 / m4) / (m3 * m3 - m4 * m4) + 2. * m1 * m1 *
log(m1 / m4) / (m1 * m1 - m4 * m4));
3749 double den = 4. * (m1 * m1 - m4 * m4)*(m3 * m3 - m4 * m4);
3751 }
else if ((std::fabs((m3 - m4) / m3) <= SUSYLEPS2)) {
3752 double num = 2. * m1 * m1 * (m1 * m1 - m4 * m4 + 2. * m1 * m1 *
log(m4 / m1)) / (m1 * m1 - m4 * m4) / (m1 * m1 - m4 * m4);
3753 num -= 2. * m2 * m2 * (m2 * m2 - m4 * m4 + 2. * m2 * m2 *
log(m4 / m2)) / (m2 * m2 - m4 * m4) / (m2 * m2 - m4 * m4);
3754 double den = 8. * (m1 * m1 - m2 * m2);
3756 }
else if ((std::fabs((m1 - m2) / m1) > SUSYLEPS2) && (std::fabs((m3 - m4) / m3) > SUSYLEPS2) && (std::fabs((m1 - m4) / m1) > SUSYLEPS2) && (std::fabs((m1 - m3) / m1) > SUSYLEPS2) && (std::fabs((m2 - m3) / m2) > SUSYLEPS2) && (std::fabs((m2 - m4) / m2) > SUSYLEPS2)) {
3757 return (m1 * m1 *
bsll_C0(m1, m3, m4) - m2 * m2 *
bsll_C0(m2, m3, m4)) / (4. * (m1 * m1 - m2 * m2));
3759 std::cout <<
"m1,m2,m3,m4 = " << m1 <<
", " << m2 <<
", " << m3 <<
", " << m4 << std::endl;
3760 throw std::runtime_error(
"SUSYMatching::bsll_D2reg(): limit/case not implemented!");
3771 double g3 =
sqrt(4. * pi *
Als);
3778 for (
int a = 0; a < 6; a++) {
3779 for (
int i = 0; i < 2; i++) {
3784 DC7_chargino *= MW * MW / (3. *
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
3789 for (
int a = 0; a < 6; a++) {
3790 for (
int i = 0; i < 4; i++) {
3795 DC7_neutralino *= -MW * MW / (3. *
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
3799 for (
int a = 0; a < 6; a++) {
3803 DC7_gluino *= 4. * MW * MW * g3 * g3 / (9. *
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
3807 return DC7_chargedHiggs + DC7_chargino + DC7_neutralino + DC7_gluino;
3818 double g3 =
sqrt(4. * pi *
Als);
3820 double tantheta =
sqrt(sw2) /
sqrt(1. - sw2);
3828 for (
int a = 0; a < 6; a++) {
3829 for (
int b = 0; b < 6; b++) {
3830 for (
int i = 0; i < 2; i++) {
3831 for (
int j = 0; j < 2; j++) {
3832 Y_chargino += (
XUL[a][1][i]).conjugate() *
XUL[b][2][j] *
3836 *
myU(i, 0) *
myU(j, 0).conjugate()) / 2.;
3845 for (
int a = 0; a < 6; a++) {
3846 for (
int b = 0; b < 6; b++) {
3847 for (
int i = 0; i < 4; i++) {
3848 for (
int j = 0; j < 4; j++) {
3849 Y_neutralino += (
ZDL[a][1][i]).conjugate() *
ZDL[b][2][j] *
3853 *(
myN(j, 2).conjugate() *
myN(i, 2) -
myN(j, 3).conjugate() *
myN(i, 3))) / 2.;
3858 Y_neutralino /= (
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
3862 for (
int a = 0; a < 6; a++) {
3863 for (
int b = 0; b < 6; b++) {
3867 Y_gluino *= g3 * g3 / (
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
3871 for (
int a = 0; a < 6; a++) {
3872 for (
int i = 0; i < 2; i++) {
3873 for (
int j = 0; j < 2; j++) {
3878 Y_chargino_box *= MW * MW / (
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
3882 for (
int a = 0; a < 6; a++) {
3883 for (
int i = 0; i < 4; i++) {
3884 for (
int j = 0; j < 4; j++) {
3885 Z_neutralino_box += (
ZDL[a][1][i]).conjugate() *
ZDL[a][2][j] / (1. - sw2)*
3891 Z_neutralino_box *= MW * MW / (
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
3894 for (
int a = 0; a < 6; a++) {
3895 for (
int i = 0; i < 4; i++) {
3896 for (
int j = 0; j < 4; j++) {
3897 Y_neutralino_box += (
ZDL[a][1][i]).conjugate() *
ZDL[a][2][j]*
3899 (tantheta *
myN(i, 0).conjugate() +
myN(i, 1).conjugate())*(tantheta *
myN(j, 0) +
myN(j, 1))
3901 *(tantheta *
myN(j, 0).conjugate() +
myN(j, 1).conjugate())*(tantheta *
myN(i, 0) +
myN(i, 1))) / 2.;
3905 Y_neutralino_box *= MW * MW / (
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
3906 Y_neutralino_box += 2. * sw2*Z_neutralino_box;
3910 return Y_chargedHiggs + Y_chargino + Y_neutralino + Y_gluino + Y_chargino_box + Y_neutralino_box;
3922 double g3 =
sqrt(4. * pi *
Als);
3933 for (
int a = 0; a < 6; a++) {
3934 for (
int b = 0; b < 6; b++) {
3935 for (
int i = 0; i < 2; i++) {
3936 for (
int j = 0; j < 2; j++) {
3937 Z_chargino += (
XUL[a][1][i]).conjugate() * (
XUL[b][2][j]) *
3941 *
myU(i, 0) *
myU(j, 0).conjugate()) / 2.;
3947 for (
int a = 0; a < 6; a++) {
3948 for (
int i = 0; i < 2; i++) {
3956 for (
int a = 0; a < 6; a++) {
3957 for (
int b = 0; b < 6; b++) {
3958 for (
int i = 0; i < 4; i++) {
3959 for (
int j = 0; j < 4; j++) {
3960 Z_neutralino += (
ZDL[a][1][i]).conjugate() *
ZDL[b][2][j] *
3964 *(
myN(j, 2).conjugate() *
myN(i, 2) -
myN(j, 3).conjugate() *
myN(i, 3))) / 2.;
3970 for (
int a = 0; a < 6; a++) {
3971 for (
int i = 0; i < 4; i++) {
3976 Z_neutralino /= (
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
3980 for (
int a = 0; a < 6; a++) {
3981 for (
int b = 0; b < 6; b++) {
3986 for (
int a = 0; a < 6; a++) {
3989 Z_gluino *= g3 * g3 / (
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
3993 for (
int a = 0; a < 6; a++) {
3994 for (
int i = 0; i < 4; i++) {
3995 for (
int j = 0; j < 4; j++) {
3996 Z_neutralino_box += (
ZDL[a][1][i]).conjugate() *
ZDL[a][2][j] / (1. - sw2)*
4002 Z_neutralino_box *= MW * MW / (
gW *
gW *
myCKM(2, 1).conjugate() *
myCKM(2, 2));
4007 return Z_chargedHiggs + Z_chargino + Z_neutralino + Z_gluino + Z_neutralino_box;
4033 for (
int a=0;a<4;a++) {
4034 Mdiag.assign(a,a,
MChi0(a));
4038 MN_tmp =
myN.transpose() * Mdiag *
myN;
4040 MN_tmp.eigensystem(cmplxONT,
MNeig);
4041 for (
int a=0;a<4;a++) {
4042 for (
int b=0;b<4;b++) {
4043 ON.assign(a,b,cmplxONT(b,a).real());
4057 double piconst = 1.0/(32.0 * pi * pi);
4058 double sw2 =
mySUSY.StandardModel::sW2(MW);
4059 double stw =
sqrt(sw2);
4060 double ctw =
sqrt(1.0 - sw2);
4061 double ttw = stw/ctw;
4067 double cdenn = MW*
cosb;
4069 double g2t = g2/
sqrt(2.0);
4074 for (
int a=0;a<4;a++) {
4075 for (
int x=0;x<6;x++) {
4077 NRlE.assign(a, x, - (g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 0) + (mE/cdenn)*
ON(a, 2)*
myRl(x, 3)));
4078 NRlMU.assign(a, x, -(g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 1) + (mMU/cdenn)*
ON(a, 2)*
myRl(x, 4)));
4079 NRlTAU.assign(a, x, -(g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 2) + (mTAU/cdenn)*
ON(a, 2)*
myRl(x, 5)));
4081 NLlE.assign(a, x, -(g2t)*((mE/cdenn)*
ON(a, 2)*
myRl(x, 0) + 2.0*
ON(a, 0)*ttw*
myRl(x, 3)));
4082 NLlMU.assign(a, x, -(g2t)*((mMU/cdenn)*
ON(a, 2)*
myRl(x, 1) + 2.0*
ON(a, 0)*ttw*
myRl(x, 4)));
4083 NLlTAU.assign(a, x, -(g2t)*((mTAU/cdenn)*
ON(a, 2)*
myRl(x, 2) + 2.0*
ON(a, 0)*ttw*
myRl(x, 5)));
4092 for (
int a=0;a<2;a++) {
4093 for (
int x=0;x<3;x++) {
4099 CLlE.assign(a, x, g2*mE/cdenc*
myU(a, 1).conjugate()*
myRn(x, 0));
4100 CLlMU.assign(a, x, g2*mMU/cdenc*
myU(a, 1).conjugate()*
myRn(x, 1));
4101 CLlTAU.assign(a, x, g2*mTAU/cdenc*
myU(a, 1).conjugate()*
myRn(x, 2));
4106 for (
int a=0;a<4;a++) {
4107 for (
int x=0;x<6;x++) {
4112 for (
int a=0;a<2;a++) {
4113 for (
int x=0;x<3;x++) {
4118 for (
int a=0;a<4;a++) {
4119 for (
int x=0;x<6;x++) {
4120 if (fabs(1.0 -
Lepty(a, x)) > 0.01) {
4126 (6.0 *
pow((1.0 -
Lepty(a,x)),4.0)) );
4138 for (
int a=0;a<2;a++) {
4139 for (
int x=0;x<3;x++) {
4140 if(fabs(1.0-
Leptz(a, x)) > 0.01) {
4146 (6.0*
pow((1.0 -
Leptz(a, x)),4.0))) );
4162 for (
int a=0;a<4;a++) {
4163 for (
int x=0;x<6;x++) {
4182 for (
int a=0;a<4;a++) {
4183 for (
int x=0;x<6;x++) {
4191 for (
int a=0;a<2;a++) {
4192 for (
int x=0;x<3;x++) {
4211 for (
int a=0;a<2;a++) {
4212 for (
int x=0;x<3;x++) {
4230 for (
int a=0;a<4;a++) {
4231 for (
int x=0;x<6;x++) {
4250 for (
int a=0;a<4;a++) {
4251 for (
int x=0;x<6;x++) {
4259 for (
int a=0;a<2;a++) {
4260 for (
int x=0;x<3;x++) {
4279 for (
int a=0;a<2;a++) {
4280 for (
int x=0;x<3;x++) {
4297 for (
int a=0;a<4;a++) {
4298 for (
int x=0;x<6;x++) {
4317 for (
int a=0;a<4;a++) {
4318 for (
int x=0;x<6;x++) {
4326 for (
int a=0;a<2;a++) {
4327 for (
int x=0;x<3;x++) {
4346 for (
int a=0;a<2;a++) {
4347 for (
int x=0;x<3;x++) {
4372 double sw2 =
mySUSY.StandardModel::sW2(MW);
4373 double stw =
sqrt(sw2);
4374 double ctw =
sqrt(1.0 - sw2);
4375 double ttw = stw/ctw;
4381 double cdenn = MW*
cosb;
4383 double g2t = g2/
sqrt(2.0);
4389 for (
int a=0;a<4;a++) {
4390 for (
int x=0;x<6;x++) {
4392 NRlE.assign(a, x, - (g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 0) + (mE/cdenn)*
ON(a, 2)*
myRl(x, 3)));
4393 NRlMU.assign(a, x, -(g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 1) + (mMU/cdenn)*
ON(a, 2)*
myRl(x, 4)));
4394 NRlTAU.assign(a, x, -(g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 2) + (mTAU/cdenn)*
ON(a, 2)*
myRl(x, 5)));
4396 NLlE.assign(a, x, -(g2t)*((mE/cdenn)*
ON(a, 2)*
myRl(x, 0) + 2.0*
ON(a, 0)*ttw*
myRl(x, 3)));
4397 NLlMU.assign(a, x, -(g2t)*((mMU/cdenn)*
ON(a, 2)*
myRl(x, 1) + 2.0*
ON(a, 0)*ttw*
myRl(x, 4)));
4398 NLlTAU.assign(a, x, -(g2t)*((mTAU/cdenn)*
ON(a, 2)*
myRl(x, 2) + 2.0*
ON(a, 0)*ttw*
myRl(x, 5)));
4407 for (
int a=0;a<2;a++) {
4408 for (
int x=0;x<3;x++) {
4414 CLlE.assign(a, x, g2*mE/cdenc*
myU(a, 1).conjugate()*
myRn(x, 0));
4415 CLlMU.assign(a, x, g2*mMU/cdenc*
myU(a, 1).conjugate()*
myRn(x, 1));
4416 CLlTAU.assign(a, x, g2*mTAU/cdenc*
myU(a, 1).conjugate()*
myRn(x, 2));
4435 for (
int a=0;a<4;a++) {
4436 for (
int b=0;b<4;b++) {
4437 for (
int x=0;x<6;x++) {
4438 for (
int t=0;t<6;t++) {
4439 B1nRMu3E = B1nRMu3E + (1.0/(4.0*pi*alph))*(0.5*
NLlMU(a,x)*
NLlE(a,t)*
NLlE(b,t)*
NLlE(b,x)
4443 B2nRMu3E = B2nRMu3E + (1.0/(4.0*pi*alph))*(0.25*(
NLlMU(a,x)*
NLlE(a,t)*
NRlE(b,t)*
NRlE(b,x)
4454 B1nLMu3E = B1nLMu3E + (1.0/(4.0*pi*alph))*(0.5*
NRlMU(a,x)*
NRlE(a,t)*
NRlE(b,t)*
NRlE(b,x)
4458 B2nLMu3E = B2nLMu3E + (1.0/(4.0*pi*alph))*(0.25*(
NRlMU(a,x)*
NRlE(a,t)*
NLlE(b,t)*
NLlE(b,x)
4484 for (
int a=0;a<2;a++) {
4485 for (
int b=0;b<2;b++) {
4486 for (
int x=0;x<3;x++) {
4487 for (
int t=0;t<3;t++) {
4488 B1cRMu3E = B1cRMu3E + (1.0/(8.0*pi*alph))*
CLlMU(a,x)*
CLlE(a,t)*
CLlE(b,t)*
CLlE(b,x)
4490 B2cRMu3E = B2cRMu3E + (1.0/(4.0*pi*alph))*(0.25*
CLlMU(a,x)*
CLlE(a,t)*
CRlE(b,t)*
CRlE(b,x)
4496 B1cLMu3E = B1cLMu3E + (1.0/(8.0*pi*alph))*(
CRlMU(a,x)*
CRlE(a,t)*
CRlE(b,t)*
CRlE(b,x)
4498 B2cLMu3E = B2cLMu3E + (1.0/(4.0*pi*alph))*(0.25*
CRlMU(a,x)*
CRlE(a,t)*
CLlE(b,t)*
CLlE(b,x)
4530 for (
int a=0;a<4;a++) {
4531 for (
int b=0;b<4;b++) {
4532 for (
int x=0;x<6;x++) {
4533 for (
int t=0;t<6;t++) {
4575 for (
int a=0;a<2;a++) {
4576 for (
int b=0;b<2;b++) {
4577 for (
int x=0;x<3;x++) {
4578 for (
int t=0;t<3;t++) {
4600 BFunctions.assign(0, B1nRTau3Mu + B1cRTau3Mu );
4601 BFunctions.assign(1, B1nLTau3Mu + B1cLTau3Mu );
4602 BFunctions.assign(2, B2nRTau3Mu + B2cRTau3Mu );
4603 BFunctions.assign(3, B2nLTau3Mu + B2cLTau3Mu );
4604 BFunctions.assign(4, B3nRTau3Mu + B3cRTau3Mu );
4605 BFunctions.assign(5, B3nLTau3Mu + B3cLTau3Mu );
4620 for (
int a=0;a<4;a++) {
4621 for (
int b=0;b<4;b++) {
4622 for (
int x=0;x<6;x++) {
4623 for (
int t=0;t<6;t++) {
4624 B1nRTau3E = B1nRTau3E + (1.0/(4.0*pi*alph))*(0.5*
NLlTAU(a,x)*
NLlE(a,t)*
NLlE(b,t)*
NLlE(b,x)
4628 B2nRTau3E = B2nRTau3E + (1.0/(4.0*pi*alph))*(0.25*(
NLlTAU(a,x)*
NLlE(a,t)*
NRlE(b,t)*
NRlE(b,x)
4639 B1nLTau3E = B1nLTau3E + (1.0/(4.0*pi*alph))*(0.5*
NRlTAU(a,x)*
NRlE(a,t)*
NRlE(b,t)*
NRlE(b,x)
4643 B2nLTau3E = B2nLTau3E + (1.0/(4.0*pi*alph))*(0.25*(
NRlTAU(a,x)*
NRlE(a,t)*
NLlE(b,t)*
NLlE(b,x)
4665 for (
int a=0;a<2;a++) {
4666 for (
int b=0;b<2;b++) {
4667 for (
int x=0;x<3;x++) {
4668 for (
int t=0;t<3;t++) {
4669 B1cRTau3E = B1cRTau3E + (1.0/(8.0*pi*alph))*
CLlTAU(a,x)*
CLlE(a,t)*
CLlE(b,t)*
CLlE(b,x)
4671 B2cRTau3E = B2cRTau3E + (1.0/(4.0*pi*alph))*(0.25*
CLlTAU(a,x)*
CLlE(a,t)*
CRlE(b,t)*
CRlE(b,x)
4677 B1cLTau3E = B1cLTau3E + (1.0/(8.0*pi*alph))*
CRlTAU(a,x)*
CRlE(a,t)*
CRlE(b,t)*
CRlE(b,x)
4679 B2cLTau3E = B2cLTau3E + (1.0/(4.0*pi*alph))*(0.25*
CRlTAU(a,x)*
CRlE(a,t)*
CLlE(b,t)*
CLlE(b,x)
4690 BFunctions.assign(0, B1nRTau3E + B1cRTau3E );
4691 BFunctions.assign(1, B1nLTau3E + B1cLTau3E );
4692 BFunctions.assign(2, B2nRTau3E + B2cRTau3E );
4693 BFunctions.assign(3, B2nLTau3E + B2cLTau3E );
4694 BFunctions.assign(4, B3nRTau3E + B3cRTau3E );
4695 BFunctions.assign(5, B3nLTau3E + B3cLTau3E );
4717 if(a == b)
return 1;
4734 double piconst = 1.0/(32.0 * pi * pi);
4735 double sw2 =
mySUSY.StandardModel::sW2(MW);
4736 double stw =
sqrt(sw2);
4737 double ctw =
sqrt(1.0 - sw2);
4738 double ttw = stw/ctw;
4758 double cdenn = MW*
cosb;
4760 double g2t = g2/
sqrt(2.0);
4796 for (
int a=0;a<4;a++) {
4797 for (
int x=0;x<6;x++) {
4799 NRlE.assign(a, x, - (g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 0) + (mE/cdenn)*
ON(a, 2)*
myRl(x, 3)));
4800 NRlMU.assign(a, x, -(g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 1) + (mMU/cdenn)*
ON(a, 2)*
myRl(x, 4)));
4801 NRlTAU.assign(a, x, -(g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 2) + (mTAU/cdenn)*
ON(a, 2)*
myRl(x, 5)));
4803 NLlE.assign(a, x, -(g2t)*((mE/cdenn)*
ON(a, 2)*
myRl(x, 0) + 2.0*
ON(a, 0)*ttw*
myRl(x, 3)));
4804 NLlMU.assign(a, x, -(g2t)*((mMU/cdenn)*
ON(a, 2)*
myRl(x, 1) + 2.0*
ON(a, 0)*ttw*
myRl(x, 4)));
4805 NLlTAU.assign(a, x, -(g2t)*((mTAU/cdenn)*
ON(a, 2)*
myRl(x, 2) + 2.0*
ON(a, 0)*ttw*
myRl(x, 5)));
4814 for (
int a=0;a<2;a++) {
4815 for (
int x=0;x<3;x++) {
4821 CLlE.assign(a, x, g2*mE/cdenc*
myU(a, 1).conjugate()*
myRn(x, 0));
4822 CLlMU.assign(a, x, g2*mMU/cdenc*
myU(a, 1).conjugate()*
myRn(x, 1));
4823 CLlTAU.assign(a, x, g2*mTAU/cdenc*
myU(a, 1).conjugate()*
myRn(x, 2));
4832 sigma1.assign(0, sina);
4833 sigma1.assign(1, -cosa);
4835 sigma2.assign(0, cosa);
4836 sigma2.assign(1, sina);
4838 sigma3.assign(0, sinapb);
4839 sigma3.assign(1, -cosapb);
4840 sigma3.assign(2, 0.);
4841 sigma4.assign(0, -sina);
4842 sigma4.assign(1, cosa);
4843 sigma4.assign(2, 0.);
4844 sigma5.assign(0, -cosbma);
4845 sigma5.assign(1, sinbma);
4849 gslpp::matrix<gslpp::complex> DL0(4, 4, 0.), DR0(4, 4, 0.), DL1(4, 4, 0.), DR1(4, 4, 0.), DL2(4, 4, 0.), DR2(4, 4, 0.);
4850 for (
int a=0;a<4;a++) {
4851 for (
int b=0;b<4;b++) {
4852 Qpp.assign(a, b, 0.5*(
ON(a,2)*(
ON(b,1)-ttw*
ON(b,0))+
ON(b,2)*(
ON(a,1)-ttw*
ON(a,0))) );
4853 Rpp.assign(a, b, (M2.
conjugate()*
ON(a,1)*
ON(b,1) +M1.
conjugate()*
ON(a,0)*
ON(b,0) -muH.
conjugate()*(
ON(a,2)*
ON(b,3)+
ON(a,3)*
ON(b,2)))/(2.0*MW) );
4854 DL0.assign(b, a, -g2/
sinb * (Qpp(a,b).conjugate()*sigma5(0) -Rpp(a,b).conjugate()*sigma2(0) +
MNeig(a)/(2.0*MW)*sigma2(0)*
delta_ab(a,b)) );
4855 DR0.assign(b, a,
DL0(b,a).conjugate() );
4856 DL1.assign(b, a, -g2/
sinb * (Qpp(a,b).conjugate()*sigma5(1) -Rpp(a,b).conjugate()*sigma2(1) +
MNeig(a)/(2.0*MW)*sigma2(1)*
delta_ab(a,b)) );
4857 DR1.assign(b, a, DL1(b,a).conjugate() );
4858 DL2.assign(b, a, -g2/
sinb * (Qpp(a,b).conjugate()*sigma5(2) -Rpp(a,b).conjugate()*sigma2(2) +
MNeig(a)/(2.0*MW)*sigma2(2)*
delta_ab(a,b)) );
4859 DR2.assign(b, a, DL2(b,a).conjugate() );
4864 gslpp::matrix<gslpp::complex> WL0(2, 2, 0.), WR0(2, 2, 0.), WL1(2, 2, 0.), WR1(2, 2, 0.), WL2(2, 2, 0.), WR2(2, 2, 0.);
4865 for (
int a=0;a<2;a++) {
4866 for (
int b=0;b<2;b++) {
4867 Qch.assign(a, b,
myU(a,1)*
myV(b,0)/
sqrt(2.0) );
4869 WR0.assign(a, b, -g2/
sinb * (Qch(a,b)*sigma5(0).conjugate() -Rch(a,b)*sigma2(0).conjugate() +
MChi(a)/(2.0*MW)*sigma2(0).conjugate()*
delta_ab(a,b)) );
4870 WL0.assign(b, a, WR0(a,b).conjugate() );
4871 WR1.assign(a, b, -g2/
sinb * (Qch(a,b)*sigma5(1).conjugate() -Rch(a,b)*sigma2(1).conjugate() +
MChi(a)/(2.0*MW)*sigma2(1).conjugate()*
delta_ab(a,b)) );
4872 WL1.assign(b, a, WR1(a,b).conjugate() );
4873 WR2.assign(a, b, -g2/
sinb * (Qch(a,b)*sigma5(2).conjugate() -Rch(a,b)*sigma2(2).conjugate() +
MChi(a)/(2.0*MW)*sigma2(2).conjugate()*
delta_ab(a,b)) );
4874 WL2.assign(b, a, WR2(a,b).conjugate() );
4882 for (
int p=0;p<3;p++) {
4883 gLLE.assign(p, MZ/ctw*sigma3(p)*(0.5-sw2) + mE*mE/(MW*
cosb)*sigma4(p));
4884 gLLMU.assign(p, MZ/ctw*sigma3(p)*(0.5-sw2) + mMU*mMU/(MW*
cosb)*sigma4(p));
4885 gLLTAU.assign(p, MZ/ctw*sigma3(p)*(0.5-sw2) + mTAU*mTAU/(MW*
cosb)*sigma4(p));
4886 gRRE.assign(p, MZ/ctw*sigma3(p)*sw2 + mE*mE/(MW*
cosb)*sigma4(p));
4887 gRRMU.assign(p, MZ/ctw*sigma3(p)*sw2 + mMU*mMU/(MW*
cosb)*sigma4(p));
4888 gRRTAU.assign(p, MZ/ctw*sigma3(p)*sw2 + mTAU*mTAU/(MW*
cosb)*sigma4(p));
4889 gLRE.assign(p, (-sigma1(p)*
TEhat(0,0)/mE*
v1/
sqrt(2.0)-sigma2(p).conjugate()*muH)*mE/(2.0*MW*
cosb));
4890 gLRMU.assign(p, (-sigma1(p)*
TEhat(1,1)/mMU*
v1/
sqrt(2.0)-sigma2(p).conjugate()*muH)*mMU/(2.0*MW*
cosb));
4891 gLRTAU.assign(p, (-sigma1(p)*
TEhat(2,2)/mTAU*
v1/
sqrt(2.0)-sigma2(p).conjugate()*muH)*mTAU/(2.0*MW*
cosb));
4892 gRLE.assign(p, gLRE(p).conjugate());
4893 gRLMU.assign(p, gLRMU(p).conjugate());
4894 gRLTAU.assign(p, gLRTAU(p).conjugate());
4895 gLLNU.assign(p, -0.5*MZ/ctw*sigma3(p));
4912 for (
int x=0;x<6;x++) {
4913 for (
int y=0;y<6;y++) {
4914 Gl0.assign(x, y, -g2*( gLLE(0)*
myRl(x,0).conjugate()*
myRl(y,0) +gRRE(0)*
myRl(x,3).conjugate()*
myRl(y,3) +gLRE(0)*
myRl(x,0).conjugate()*
myRl(y,3) +gRLE(0)*
myRl(x,3).conjugate()*
myRl(y,0)
4915 +gLLMU(0)*
myRl(x,1).conjugate()*
myRl(y,1) +gRRMU(0)*
myRl(x,4).conjugate()*
myRl(y,4) +gLRMU(0)*
myRl(x,1).conjugate()*
myRl(y,4) +gRLMU(0)*
myRl(x,4).conjugate()*
myRl(y,1)
4916 +gLLTAU(0)*
myRl(x,2).conjugate()*
myRl(y,2) +gRRTAU(0)*
myRl(x,5).conjugate()*
myRl(y,5) +gLRTAU(0)*
myRl(x,2).conjugate()*
myRl(y,5) +gRLTAU(0)*
myRl(x,5).conjugate()*
myRl(y,2)));
4917 Gl1.assign(x, y, -g2*( gLLE(1)*
myRl(x,0).conjugate()*
myRl(y,0) +gRRE(1)*
myRl(x,3).conjugate()*
myRl(y,3) +gLRE(1)*
myRl(x,0).conjugate()*
myRl(y,3) +gRLE(1)*
myRl(x,3).conjugate()*
myRl(y,0)
4918 +gLLMU(1)*
myRl(x,1).conjugate()*
myRl(y,1) +gRRMU(1)*
myRl(x,4).conjugate()*
myRl(y,4) +gLRMU(1)*
myRl(x,1).conjugate()*
myRl(y,4) +gRLMU(1)*
myRl(x,4).conjugate()*
myRl(y,1)
4919 +gLLTAU(1)*
myRl(x,2).conjugate()*
myRl(y,2) +gRRTAU(1)*
myRl(x,5).conjugate()*
myRl(y,5) +gLRTAU(1)*
myRl(x,2).conjugate()*
myRl(y,5) +gRLTAU(1)*
myRl(x,5).conjugate()*
myRl(y,2)));
4920 Gl2.assign(x, y, -g2*( gLLE(2)*
myRl(x,0).conjugate()*
myRl(y,0) +gRRE(2)*
myRl(x,3).conjugate()*
myRl(y,3) +gLRE(2)*
myRl(x,0).conjugate()*
myRl(y,3) +gRLE(2)*
myRl(x,3).conjugate()*
myRl(y,0)
4921 +gLLMU(2)*
myRl(x,1).conjugate()*
myRl(y,1) +gRRMU(2)*
myRl(x,4).conjugate()*
myRl(y,4) +gLRMU(2)*
myRl(x,1).conjugate()*
myRl(y,4) +gRLMU(2)*
myRl(x,4).conjugate()*
myRl(y,1)
4922 +gLLTAU(2)*
myRl(x,2).conjugate()*
myRl(y,2) +gRRTAU(2)*
myRl(x,5).conjugate()*
myRl(y,5) +gLRTAU(2)*
myRl(x,2).conjugate()*
myRl(y,5) +gRLTAU(2)*
myRl(x,5).conjugate()*
myRl(y,2)));
4927 for (
int x=0;x<3;x++) {
4928 Gnu0.assign(x, x, -g2*gLLNU(0) );
4929 Gnu1.assign(x, x, -g2*gLLNU(1) );
4930 Gnu2.assign(x, x, -g2*gLLNU(2) );
4933 gslpp::vector<gslpp::complex> SRE(3, 0.), SLE(3, 0.), SRMU(3, 0.), SLMU(3, 0.), SRTAU(3, 0.), SLTAU(3, 0.);
4934 for (
int p=0;p<3;p++) {
4935 SRE.assign(p, g2*mE/(2.0*MW*
cosb) * sigma1(p));
4936 SLE.assign(p, g2*mE/(2.0*MW*
cosb) * sigma1(p).conjugate());
4937 SRMU.assign(p, g2*mMU/(2.0*MW*
cosb) * sigma1(p));
4938 SLMU.assign(p, g2*mMU/(2.0*MW*
cosb) * sigma1(p).conjugate());
4939 SRTAU.assign(p, g2*mTAU/(2.0*MW*
cosb) * sigma1(p));
4940 SLTAU.assign(p, g2*mTAU/(2.0*MW*
cosb) * sigma1(p).conjugate());
4952 for (
int x=0;x<6;x++) {
4953 for (
int a=0;a<4;a++) {
4954 for (
int b=0;b<4;b++) {
5058 for (
int y=0;y<6;y++) {
5148 gslpp::complex B2HiggsnR = (-0.5*HpengMuEEENR0*SLE(0)/(mh*mh)-0.5*HpengMuEEENR1*SLE(1)/(mH*mH)-0.5*HpengMuEEENR2*SLE(2)/(mA*mA))/(4.0*pi*alph);
5149 gslpp::complex B2HiggsnL = (-0.5*HpengMuEEENL0*SRE(0)/(mh*mh)-0.5*HpengMuEEENL1*SRE(1)/(mH*mH)-0.5*HpengMuEEENL2*SRE(2)/(mA*mA))/(4.0*pi*alph);
5150 gslpp::complex B3HiggsnR = (HpengMuEEENR0*SRE(0)/(mh*mh)+HpengMuEEENR1*SRE(1)/(mH*mH)+HpengMuEEENR2*SRE(2)/(mA*mA))/(4.0*pi*alph);
5151 gslpp::complex B3HiggsnL = (HpengMuEEENL0*SLE(0)/(mh*mh)+HpengMuEEENL1*SLE(1)/(mH*mH)+HpengMuEEENL2*SLE(2)/(mA*mA))/(4.0*pi*alph);
5160 for (
int x=0;x<3;x++) {
5161 for (
int a=0;a<2;a++) {
5162 for (
int b=0;b<2;b++) {
5266 for (
int y=0;y<3;y++) {
5299 HpengMuEEECR0 = HpengMuEEECR0 - 2.0*piconst*(SRMU(0)/(mE*mE-mMU*mMU)*(-
CRlE(a,x)*
CRlMU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5308 HpengMuEEECL0 = HpengMuEEECL0 - 2.0*piconst*(SLMU(0)/(mE*mE-mMU*mMU)*(-
CLlE(a,x)*
CLlMU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5317 HpengMuEEECR1 = HpengMuEEECR1 - 2.0*piconst*(SRMU(1)/(mE*mE-mMU*mMU)*(-
CRlE(a,x)*
CRlMU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5326 HpengMuEEECL1 = HpengMuEEECL1 - 2.0*piconst*(SLMU(1)/(mE*mE-mMU*mMU)*(-
CLlE(a,x)*
CLlMU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5335 HpengMuEEECR2 = HpengMuEEECR2 - 2.0*piconst*(SRMU(2)/(mE*mE-mMU*mMU)*(-
CRlE(a,x)*
CRlMU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5344 HpengMuEEECL2 = HpengMuEEECL2 - 2.0*piconst*(SLMU(2)/(mE*mE-mMU*mMU)*(-
CLlE(a,x)*
CLlMU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5356 gslpp::complex B2HiggscR = (-0.5*HpengMuEEECR0*SLE(0)/(mh*mh)-0.5*HpengMuEEECR1*SLE(1)/(mH*mH)-0.5*HpengMuEEECR2*SLE(2)/(mA*mA))/(4.0*pi*alph);
5357 gslpp::complex B2HiggscL = (-0.5*HpengMuEEECL0*SRE(0)/(mh*mh)-0.5*HpengMuEEECL1*SRE(1)/(mH*mH)-0.5*HpengMuEEECL2*SRE(2)/(mA*mA))/(4.0*pi*alph);
5358 gslpp::complex B3HiggscR = (HpengMuEEECR0*SRE(0)/(mh*mh)+HpengMuEEECR1*SRE(1)/(mH*mH)+HpengMuEEECR2*SRE(2)/(mA*mA))/(4.0*pi*alph);
5359 gslpp::complex B3HiggscL = (HpengMuEEECL0*SLE(0)/(mh*mh)+HpengMuEEECL1*SLE(1)/(mH*mH)+HpengMuEEECL2*SLE(2)/(mA*mA))/(4.0*pi*alph);
5377 for (
int x=0;x<6;x++) {
5378 for (
int a=0;a<4;a++) {
5379 for (
int b=0;b<4;b++) {
5483 for (
int y=0;y<6;y++) {
5516 HpengTauMUMUMUNR0 = HpengTauMUMUMUNR0 - 2.0*piconst*(SRTAU(0)/(mMU*mMU-mTAU*mTAU)*(-
NRlMU(a,x)*
NRlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MNeig(a)*
MNeig(a),
mym_se_sq(x))
5525 HpengTauMUMUMUNL0 = HpengTauMUMUMUNL0 - 2.0*piconst*(SLTAU(0)/(mMU*mMU-mTAU*mTAU)*(-
NLlMU(a,x)*
NLlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MNeig(a)*
MNeig(a),
mym_se_sq(x))
5534 HpengTauMUMUMUNR1 = HpengTauMUMUMUNR1 - 2.0*piconst*(SRTAU(1)/(mMU*mMU-mTAU*mTAU)*(-
NRlMU(a,x)*
NRlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MNeig(a)*
MNeig(a),
mym_se_sq(x))
5543 HpengTauMUMUMUNL1 = HpengTauMUMUMUNL1 - 2.0*piconst*(SLTAU(1)/(mMU*mMU-mTAU*mTAU)*(-
NLlMU(a,x)*
NLlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MNeig(a)*
MNeig(a),
mym_se_sq(x))
5552 HpengTauMUMUMUNR2 = HpengTauMUMUMUNR2 - 2.0*piconst*(SRTAU(2)/(mMU*mMU-mTAU*mTAU)*(-
NRlMU(a,x)*
NRlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MNeig(a)*
MNeig(a),
mym_se_sq(x))
5561 HpengTauMUMUMUNL2 = HpengTauMUMUMUNL2 - 2.0*piconst*(SLTAU(2)/(mMU*mMU-mTAU*mTAU)*(-
NLlMU(a,x)*
NLlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MNeig(a)*
MNeig(a),
mym_se_sq(x))
5573 gslpp::complex B2HiggsnRtm = (-0.5*HpengTauMUMUMUNR0*SLMU(0)/(mh*mh)-0.5*HpengTauMUMUMUNR1*SLMU(1)/(mH*mH)-0.5*HpengTauMUMUMUNR2*SLMU(2)/(mA*mA))/(4.0*pi*alph);
5574 gslpp::complex B2HiggsnLtm = (-0.5*HpengTauMUMUMUNL0*SRMU(0)/(mh*mh)-0.5*HpengTauMUMUMUNL1*SRMU(1)/(mH*mH)-0.5*HpengTauMUMUMUNL2*SRMU(2)/(mA*mA))/(4.0*pi*alph);
5575 gslpp::complex B3HiggsnRtm = (HpengTauMUMUMUNR0*SRMU(0)/(mh*mh)+HpengTauMUMUMUNR1*SRMU(1)/(mH*mH)+HpengTauMUMUMUNR2*SRMU(2)/(mA*mA))/(4.0*pi*alph);
5576 gslpp::complex B3HiggsnLtm = (HpengTauMUMUMUNL0*SLMU(0)/(mh*mh)+HpengTauMUMUMUNL1*SLMU(1)/(mH*mH)+HpengTauMUMUMUNL2*SLMU(2)/(mA*mA))/(4.0*pi*alph);
5585 for (
int x=0;x<3;x++) {
5586 for (
int a=0;a<2;a++) {
5587 for (
int b=0;b<2;b++) {
5691 for (
int y=0;y<3;y++) {
5724 HpengTauMUMUMUCR0 = HpengTauMUMUMUCR0 - 2.0*piconst*(SRTAU(0)/(mMU*mMU-mTAU*mTAU)*(-
CRlMU(a,x)*
CRlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5733 HpengTauMUMUMUCL0 = HpengTauMUMUMUCL0 - 2.0*piconst*(SLTAU(0)/(mMU*mMU-mTAU*mTAU)*(-
CLlMU(a,x)*
CLlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5742 HpengTauMUMUMUCR1 = HpengTauMUMUMUCR1 - 2.0*piconst*(SRTAU(1)/(mMU*mMU-mTAU*mTAU)*(-
CRlMU(a,x)*
CRlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5751 HpengTauMUMUMUCL1 = HpengTauMUMUMUCL1 - 2.0*piconst*(SLTAU(1)/(mMU*mMU-mTAU*mTAU)*(-
CLlMU(a,x)*
CLlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5760 HpengTauMUMUMUCR2 = HpengTauMUMUMUCR2 - 2.0*piconst*(SRTAU(2)/(mMU*mMU-mTAU*mTAU)*(-
CRlMU(a,x)*
CRlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5769 HpengTauMUMUMUCL2 = HpengTauMUMUMUCL2 - 2.0*piconst*(SLTAU(2)/(mMU*mMU-mTAU*mTAU)*(-
CLlMU(a,x)*
CLlTAU(a,x).conjugate()*mMU*mMU*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
5781 gslpp::complex B2HiggscRtm = (-0.5*HpengTauMUMUMUCR0*SLMU(0)/(mh*mh)-0.5*HpengTauMUMUMUCR1*SLMU(1)/(mH*mH)-0.5*HpengTauMUMUMUCR2*SLMU(2)/(mA*mA))/(4.0*pi*alph);
5782 gslpp::complex B2HiggscLtm = (-0.5*HpengTauMUMUMUCL0*SRMU(0)/(mh*mh)-0.5*HpengTauMUMUMUCL1*SRMU(1)/(mH*mH)-0.5*HpengTauMUMUMUCL2*SRMU(2)/(mA*mA))/(4.0*pi*alph);
5783 gslpp::complex B3HiggscRtm = (HpengTauMUMUMUCR0*SRMU(0)/(mh*mh)+HpengTauMUMUMUCR1*SRMU(1)/(mH*mH)+HpengTauMUMUMUCR2*SRMU(2)/(mA*mA))/(4.0*pi*alph);
5784 gslpp::complex B3HiggscLtm = (HpengTauMUMUMUCL0*SLMU(0)/(mh*mh)+HpengTauMUMUMUCL1*SLMU(1)/(mH*mH)+HpengTauMUMUMUCL2*SLMU(2)/(mA*mA))/(4.0*pi*alph);
5802 for (
int x=0;x<6;x++) {
5803 for (
int a=0;a<4;a++) {
5804 for (
int b=0;b<4;b++) {
5908 for (
int y=0;y<6;y++) {
5998 gslpp::complex B2HiggsnRte = (-0.5*HpengTauEEENR0*SLE(0)/(mh*mh)-0.5*HpengTauEEENR1*SLE(1)/(mH*mH)-0.5*HpengTauEEENR2*SLE(2)/(mA*mA))/(4.0*pi*alph);
5999 gslpp::complex B2HiggsnLte = (-0.5*HpengTauEEENL0*SRE(0)/(mh*mh)-0.5*HpengTauEEENL1*SRE(1)/(mH*mH)-0.5*HpengTauEEENL2*SRE(2)/(mA*mA))/(4.0*pi*alph);
6000 gslpp::complex B3HiggsnRte = (HpengTauEEENR0*SRE(0)/(mh*mh)+HpengTauEEENR1*SRE(1)/(mH*mH)+HpengTauEEENR2*SRE(2)/(mA*mA))/(4.0*pi*alph);
6001 gslpp::complex B3HiggsnLte = (HpengTauEEENL0*SLE(0)/(mh*mh)+HpengTauEEENL1*SLE(1)/(mH*mH)+HpengTauEEENL2*SLE(2)/(mA*mA))/(4.0*pi*alph);
6010 for (
int x=0;x<3;x++) {
6011 for (
int a=0;a<2;a++) {
6012 for (
int b=0;b<2;b++) {
6116 for (
int y=0;y<3;y++) {
6149 HpengTauEEECR0 = HpengTauEEECR0 - 2.0*piconst*(SRTAU(0)/(mE*mE-mTAU*mTAU)*(-
CRlE(a,x)*
CRlTAU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
6158 HpengTauEEECL0 = HpengTauEEECL0 - 2.0*piconst*(SLTAU(0)/(mE*mE-mTAU*mTAU)*(-
CLlE(a,x)*
CLlTAU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
6167 HpengTauEEECR1 = HpengTauEEECR1 - 2.0*piconst*(SRTAU(1)/(mE*mE-mTAU*mTAU)*(-
CRlE(a,x)*
CRlTAU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
6176 HpengTauEEECL1 = HpengTauEEECL1 - 2.0*piconst*(SLTAU(1)/(mE*mE-mTAU*mTAU)*(-
CLlE(a,x)*
CLlTAU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
6185 HpengTauEEECR2 = HpengTauEEECR2 - 2.0*piconst*(SRTAU(2)/(mE*mE-mTAU*mTAU)*(-
CRlE(a,x)*
CRlTAU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
6194 HpengTauEEECL2 = HpengTauEEECL2 - 2.0*piconst*(SLTAU(2)/(mE*mE-mTAU*mTAU)*(-
CLlE(a,x)*
CLlTAU(a,x).conjugate()*mE*mE*
PV.
B0(1.,0.,
MChi(a)*
MChi(a),
mym_sn_sq(x))
6206 gslpp::complex B2HiggscRte = (-0.5*HpengTauEEECR0*SLE(0)/(mh*mh)-0.5*HpengTauEEECR1*SLE(1)/(mH*mH)-0.5*HpengTauEEECR2*SLE(2)/(mA*mA))/(4.0*pi*alph);
6207 gslpp::complex B2HiggscLte = (-0.5*HpengTauEEECL0*SRE(0)/(mh*mh)-0.5*HpengTauEEECL1*SRE(1)/(mH*mH)-0.5*HpengTauEEECL2*SRE(2)/(mA*mA))/(4.0*pi*alph);
6208 gslpp::complex B3HiggscRte = (HpengTauEEECR0*SRE(0)/(mh*mh)+HpengTauEEECR1*SRE(1)/(mH*mH)+HpengTauEEECR2*SRE(2)/(mA*mA))/(4.0*pi*alph);
6209 gslpp::complex B3HiggscLte = (HpengTauEEECL0*SLE(0)/(mh*mh)+HpengTauEEECL1*SLE(1)/(mH*mH)+HpengTauEEECL2*SLE(2)/(mA*mA))/(4.0*pi*alph);
6228 double sw2 =
mySUSY.StandardModel::sW2(MW);
6229 double stw =
sqrt(sw2);
6230 double ctw =
sqrt(1.0 - sw2);
6231 double ttw = stw/ctw;
6237 double cdenn = MW*
cosb;
6239 double g2t = g2/
sqrt(2.0);
6244 for (
int a=0;a<4;a++) {
6245 for (
int x=0;x<6;x++) {
6247 NRlE.assign(a, x, - (g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 0) + (mE/cdenn)*
ON(a, 2)*
myRl(x, 3)));
6248 NRlMU.assign(a, x, -(g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 1) + (mMU/cdenn)*
ON(a, 2)*
myRl(x, 4)));
6252 NLlE.assign(a, x, -(g2t)*((mE/cdenn)*
ON(a, 2)*
myRl(x, 0) + 2.0*
ON(a, 0)*ttw*
myRl(x, 3)));
6253 NLlMU.assign(a, x, -(g2t)*((mMU/cdenn)*
ON(a, 2)*
myRl(x, 1) + 2.0*
ON(a, 0)*ttw*
myRl(x, 4)));
6259 for (
int a=0;a<2;a++) {
6260 for (
int x=0;x<3;x++) {
6265 CLlE.assign(a, x, g2*mE/cdenc*
myU(a, 1).conjugate()*
myRn(x, 0));
6266 CLlMU.assign(a, x, g2*mMU/cdenc*
myU(a, 1).conjugate()*
myRn(x, 1));
6269 for (
int a=0;a<2;a++) {
6270 for (
int x=0;x<6;x++) {
6285 for (
int a=0;a<4;a++) {
6286 for (
int b=0;b<4;b++) {
6287 for (
int x=0;x<6;x++) {
6288 for (
int y=0;y<6;y++) {
6323 for (
int a=0;a<2;a++) {
6324 for (
int b=0;b<2;b++) {
6325 for (
int x=0;x<3;x++) {
6326 for (
int y=0;y<6;y++) {
6365 double piconst = 1.0/(32.0 * pi * pi);
6366 double sw2 =
mySUSY.StandardModel::sW2(MW);
6367 double stw =
sqrt(sw2);
6368 double ctw =
sqrt(1.0 - sw2);
6369 double ttw = stw/ctw;
6375 double cdenn = MW*
cosb;
6377 double g2t = g2/
sqrt(2.0);
6384 for (
int a=0;a<4;a++) {
6385 for (
int x=0;x<6;x++) {
6387 NRlE.assign(a, x, - (g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 0) + (mE/cdenn)*
ON(a, 2)*
myRl(x, 3)));
6388 NRlMU.assign(a, x, -(g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 1) + (mMU/cdenn)*
ON(a, 2)*
myRl(x, 4)));
6389 NRlTAU.assign(a, x, -(g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 2) + (mTAU/cdenn)*
ON(a, 2)*
myRl(x, 5)));
6391 NLlE.assign(a, x, -(g2t)*((mE/cdenn)*
ON(a, 2)*
myRl(x, 0) + 2.0*
ON(a, 0)*ttw*
myRl(x, 3)));
6392 NLlMU.assign(a, x, -(g2t)*((mMU/cdenn)*
ON(a, 2)*
myRl(x, 1) + 2.0*
ON(a, 0)*ttw*
myRl(x, 4)));
6393 NLlTAU.assign(a, x, -(g2t)*((mTAU/cdenn)*
ON(a, 2)*
myRl(x, 2) + 2.0*
ON(a, 0)*ttw*
myRl(x, 5)));
6402 for (
int a=0;a<2;a++) {
6403 for (
int x=0;x<3;x++) {
6409 CLlE.assign(a, x, g2*mE/cdenc*
myU(a, 1).conjugate()*
myRn(x, 0));
6410 CLlMU.assign(a, x, g2*mMU/cdenc*
myU(a, 1).conjugate()*
myRn(x, 1));
6411 CLlTAU.assign(a, x, g2*mTAU/cdenc*
myU(a, 1).conjugate()*
myRn(x, 2));
6417 for (
int a=0;a<4;a++) {
6418 for (
int x=0;x<6;x++) {
6423 for (
int x=0;x<6;x++) {
6424 for (
int a=0;a<4;a++) {
6425 for (
int b=0;b<4;b++) {
6426 if (a != b && std::fabs(
Lepty(a,x)-
Lepty(b,x)) > 0.01 && std::fabs(1.0-
Lepty(a,x)) > 0.01 && std::fabs(1.0-
Lepty(b,x)) > 0.01) {
6432 else if (a != b && std::fabs(
Lepty(a,x)-
Lepty(b,x)) > 0.01 && std::fabs(1.0-
Lepty(a,x)) > 0.01 && std::fabs(1.0-
Lepty(b,x)) <= 0.01) {
6436 else if (a != b && std::fabs(
Lepty(a,x)-
Lepty(b,x)) > 0.01 && std::fabs(1.0-
Lepty(b,x)) > 0.01 && std::fabs(1.0-
Lepty(a,x)) <= 0.01) {
6440 else if ((a == b || std::fabs(
Lepty(a,x)-
Lepty(b,x)) <= 0.01) && std::fabs(1.0-
Lepty(a,x)) > 0.01 && std::fabs(1.0-
Lepty(b,x)) > 0.01) {
6453 for (
int a=0;a<2;a++) {
6454 for (
int x=0;x<3;x++) {
6459 for (
int x=0;x<3;x++) {
6460 for (
int a=0;a<2;a++) {
6461 for (
int b=0;b<2;b++) {
6462 if (a != b && std::fabs(
Leptz(a,x)-
Leptz(b,x)) > 0.01 && std::fabs(1.0-
Leptz(a,x)) > 0.01 && std::fabs(1.0-
Leptz(b,x)) > 0.01) {
6468 else if (a != b && std::fabs(
Leptz(a,x)-
Leptz(b,x)) > 0.01 && std::fabs(1.0-
Leptz(a,x)) > 0.01 && std::fabs(1.0-
Leptz(b,x)) <= 0.01) {
6472 else if (a != b && std::fabs(
Leptz(a,x)-
Leptz(b,x)) > 0.01 && std::fabs(1.0-
Leptz(b,x)) > 0.01 && std::fabs(1.0-
Leptz(a,x)) <= 0.01) {
6476 else if ((a == b || std::fabs(
Leptz(a,x)-
Leptz(b,x)) <= 0.01) && std::fabs(1.0-
Leptz(a,x)) > 0.01 && std::fabs(1.0-
Leptz(b,x)) > 0.01) {
6493 for (
int x=0;x<6;x++) {
6494 for (
int a=0;a<4;a++) {
6495 for (
int b=0;b<4;b++) {
6496 ZpengMuEEENR = ZpengMuEEENR -
NLlE(a,x)*
NLlMU(b,x)*piconst*
6498 ZpengMuEEENL = ZpengMuEEENL +
NRlE(a,x)*
NRlMU(b,x)*piconst*
6505 for (
int x=0;x<3;x++) {
6506 for (
int a=0;a<2;a++) {
6507 for (
int b=0;b<2;b++) {
6508 ZpengMuEEEC = ZpengMuEEEC +
CRlE(a,x)*
CRlMU(b,x)*piconst*
6515 FFunctions.assign(0, ZpengMuEEENR/(MZ*MZ*ctw*ctw) );
6516 FFunctions.assign(1, ZpengMuEEENR*(-0.5+sw2)/(MZ*MZ*sw2*ctw*ctw) );
6517 FFunctions.assign(2, (ZpengMuEEENL + ZpengMuEEEC)/(MZ*MZ*ctw*ctw) );
6518 FFunctions.assign(3, (ZpengMuEEENL + ZpengMuEEEC)*(-0.5+sw2)/(MZ*MZ*sw2*ctw*ctw) );
6525 for (
int x=0;x<6;x++) {
6526 for (
int a=0;a<4;a++) {
6527 for (
int b=0;b<4;b++) {
6528 ZpengTauMuMuMuNR = ZpengTauMuMuMuNR -
NLlMU(a,x)*
NLlTAU(b,x)*piconst*
6530 ZpengTauMuMuMuNL = ZpengTauMuMuMuNL +
NRlMU(a,x)*
NRlTAU(b,x)*piconst*
6537 for (
int x=0;x<3;x++) {
6538 for (
int a=0;a<2;a++) {
6539 for (
int b=0;b<2;b++) {
6540 ZpengTauMuMuMuC = ZpengTauMuMuMuC +
CRlMU(a,x)*
CRlTAU(b,x)*piconst*
6547 FFunctions.assign(0, ZpengTauMuMuMuNR/(MZ*MZ*ctw*ctw) );
6548 FFunctions.assign(1, ZpengTauMuMuMuNR*(-0.5+sw2)/(MZ*MZ*sw2*ctw*ctw) );
6549 FFunctions.assign(2, (ZpengTauMuMuMuNL + ZpengTauMuMuMuC)/(MZ*MZ*ctw*ctw) );
6550 FFunctions.assign(3, (ZpengTauMuMuMuNL + ZpengTauMuMuMuC)*(-0.5+sw2)/(MZ*MZ*sw2*ctw*ctw) );
6557 for (
int x=0;x<6;x++) {
6558 for (
int a=0;a<4;a++) {
6559 for (
int b=0;b<4;b++) {
6560 ZpengTauEEENR = ZpengTauEEENR -
NLlE(a,x)*
NLlTAU(b,x)*piconst*
6562 ZpengTauEEENL = ZpengTauEEENL +
NRlE(a,x)*
NRlTAU(b,x)*piconst*
6569 for (
int x=0;x<3;x++) {
6570 for (
int a=0;a<2;a++) {
6571 for (
int b=0;b<2;b++) {
6572 ZpengTauEEEC = ZpengTauEEEC +
CRlE(a,x)*
CRlTAU(b,x)*piconst*
6579 FFunctions.assign(0, ZpengTauEEENR/(MZ*MZ*ctw*ctw) );
6580 FFunctions.assign(1, ZpengTauEEENR*(-0.5+sw2)/(MZ*MZ*sw2*ctw*ctw) );
6581 FFunctions.assign(2, (ZpengTauEEENL + ZpengTauEEEC)/(MZ*MZ*ctw*ctw) );
6582 FFunctions.assign(3, (ZpengTauEEENL + ZpengTauEEEC)*(-0.5+sw2)/(MZ*MZ*sw2*ctw*ctw) );
6589 for (
int x=0;x<6;x++) {
6590 for (
int a=0;a<4;a++) {
6591 for (
int b=0;b<4;b++) {
6592 ZpengTauMuEENR = ZpengTauMuEENR -
NLlMU(a,x)*
NLlTAU(b,x)*piconst*
6594 ZpengTauMuEENL = ZpengTauMuEENL +
NRlMU(a,x)*
NRlTAU(b,x)*piconst*
6601 for (
int x=0;x<3;x++) {
6602 for (
int a=0;a<2;a++) {
6603 for (
int b=0;b<2;b++) {
6604 ZpengTauMuEEC = ZpengTauMuEEC +
CRlMU(a,x)*
CRlTAU(b,x)*piconst*
6611 FFunctions.assign(0, ZpengTauMuEENR/(MZ*MZ*ctw*ctw) );
6612 FFunctions.assign(1, ZpengTauMuEENR*(-0.5+sw2)/(MZ*MZ*sw2*ctw*ctw) );
6613 FFunctions.assign(2, (ZpengTauMuEENL + ZpengTauMuEEC)/(MZ*MZ*ctw*ctw) );
6614 FFunctions.assign(3, (ZpengTauMuEENL + ZpengTauMuEEC)*(-0.5+sw2)/(MZ*MZ*sw2*ctw*ctw) );
6626 double piconst = 1.0/(32.0 * pi * pi);
6627 double sw2 =
mySUSY.StandardModel::sW2(MW);
6628 double stw =
sqrt(sw2);
6629 double ctw =
sqrt(1.0 - sw2);
6630 double ttw = stw/ctw;
6634 double cdenn = MW*
cosb;
6636 double g2t = g2/
sqrt(2.0);
6641 for (
int a=0;a<4;a++) {
6642 for (
int x=0;x<6;x++) {
6644 NRlMU.assign(a, x, -(g2t)*((-
ON(a, 1) -
ON(a, 0)*ttw)*
myRl(x, 1) + (mMU/cdenn)*
ON(a, 2)*
myRl(x, 4)));
6646 NLlMU.assign(a, x, -(g2t)*((mMU/cdenn)*
ON(a, 2)*
myRl(x, 1) + 2.0*
ON(a, 0)*ttw*
myRl(x, 4)));
6653 for (
int a=0;a<2;a++) {
6654 for (
int x=0;x<3;x++) {
6658 CLlMU.assign(a, x, g2*mMU/cdenc*
myU(a, 1).conjugate()*
myRn(x, 1));
6662 for (
int a=0;a<4;a++) {
6663 for (
int x=0;x<6;x++) {
6668 for (
int a=0;a<2;a++) {
6669 for (
int x=0;x<3;x++) {
6674 for (
int a=0;a<4;a++) {
6675 for (
int x=0;x<6;x++) {
6676 if (fabs(1.0 -
Lepty(a, x)) > 0.01) {
6679 (6.0 *
pow((1.0 -
Lepty(a,x)),4.0)) );
6690 for (
int a=0;a<2;a++) {
6691 for (
int x=0;x<3;x++) {
6692 if(fabs(1.0-
Leptz(a, x)) > 0.01) {
6695 (6.0*
pow((1.0 -
Leptz(a, x)),4.0))) );
6710 for (
int a=0;a<4;a++) {
6711 for (
int x=0;x<6;x++) {
6712 g2ARN = g2ARN -mMU*mMU*piconst*(4.0*
NRlMU(a,x)*
NRlMU(a,x).conjugate()*
Leptf1(a,x)
6714 g2ALN = g2ALN -mMU*mMU*piconst*(4.0*
NLlMU(a,x)*
NLlMU(a,x).conjugate()*
Leptf1(a,x)
6722 for (
int a=0;a<2;a++) {
6723 for (
int x=0;x<3;x++) {
6724 g2ARC = g2ARC +mMU*mMU*piconst*(4.0*
CRlMU(a,x)*
CRlMU(a,x).conjugate()*
Leptf3(a,x)
6726 g2ALC = g2ALC +mMU*mMU*piconst*(4.0*
CLlMU(a,x)*
CLlMU(a,x).conjugate()*
Leptf3(a,x)
6746 __fPS_params ¶ms= *reinterpret_cast<__fPS_params *>(p);
6747 double r = params.
a*
log(x*(1.-x)/params.
a) / (x*(1.-x)-params.
a);
6757 gsl_integration_glfixed_table * w
6758 = gsl_integration_glfixed_table_alloc(100);
6762 F.params = reinterpret_cast<void *>(¶ms);
6764 result = gsl_integration_glfixed (&F, 0, 1, w);
6766 gsl_integration_glfixed_table_free (w);
6774 r=(2.*x-1.)*
fPS(x)-2.*x*(2.+
log(x));
6781 r=x*(2.+
log(x)-
fPS(x))/2.0;
6788 if (std::fabs(a-b) < 0.0005 && std::fabs(a-c) < 0.0005 && std::fabs(b-c) < 0.0005)
6790 r=9./(2.*(a+b+c)*(a+b+c));
6794 if (std::fabs(a-b) < 0.0005)
6796 if (std::fabs(a-c) < 0.0005)
6798 r=(-4.*(4.*b*b*c*c*
log(b*b/(c*c))
6799 +((b+c)*(b+c))*(c*c*
log((4.*c*c)/((b+c)*(b+c)))
6800 + b*b*
log(((b+c)*(b+c))/(4.*b*b)))))
6801 /(
pow(b-c,3)*(b+c)*(3.*b+c)*(b+3.*c));
6805 if (std::fabs(b-c) < 0.0005)
6807 r=(4.*(
pow(a-2.*b+c,2)*(-3.*
pow(a,8) - 46.*
pow(a,7)*c - 42.*
pow(a,6)*c*c + 234.*
pow(a,5)*c*c*c
6808 - 234.*a*a*a*
pow(c,5) + 42.*a*a*
pow(c,6) + 46.*a*
pow(c,7) + 3.*
pow(c,8)
6809 + 8.*a*a*c*c*(11.*
pow(a,4) + 28.*a*a*a*c + 50.*a*a*c*c + 28.*a*c*c*c + 11.*
pow(c,4))
6811 + a*a*(7.*
pow(a,6) + 62.*
pow(a,5)*c + 265.*
pow(a,4)*c*c
6812 + 356.*a*a*a*c*c*c + 265.*a*a*
pow(c,4) + 62.*a*
pow(c,5) + 7.*
pow(c,6))
6813 *
log((4.*a*a)/((a+c)*(a+c)))
6814 + 7.*
pow(a,6)*c*c*
log(((a+c)*(a+c))/(4.*c*c))
6815 + 62.*
pow(a,5)*c*c*c*
log(((a+c)*(a+c))/(4.*c*c))
6816 + 265.*
pow(a,4)*
pow(c,4)*
log(((a+c)*(a+c))/(4.*c*c))
6817 + 356.*a*a*a*
pow(c,5)*
log(((a+c)*(a+c))/(4.*c*c))
6818 + 265.*a*a*
pow(c,6)*
log(((a+c)*(a+c))/(4.*c*c))
6819 + 62.*a*
pow(c,7)*
log(((a+c)*(a+c))/(4.*c*c))
6820 + 7.*
pow(c,8)*
log(((a+c)*(a+c))/(4.*c*c)))
6821 +((a-c)*(a-c))*
pow(3.*a+c,2)*
pow(a+3.*c,2)
6822 *(4.*a*a*c*c*
log((c*c)/(a*a))
6823 +((a+c)*(a+c))*(a*a*
log((4.*a*a)/((a+c)*(a+c))) + c*c*
log(((a+c)*(a+c))/(4.*c*c))))))
6824 /(
pow(a-c,5)*(a+c)*
pow(3.*a+c,3)*
pow(a+3.*c,3));
6828 r=((b*b) - (c*c) - (c*c)*
log((b*b)/(c*c)))/
pow((b*b) - (c*c),2) +
6829 ((a - b)*(-
pow(b,4) +
pow(c,4) + 2*(b*b)*(c*c)*
log((b*b)/(c*c))))/(b*
pow((b*b) - (c*c),3));
6835 if (std::fabs(b-c) < 0.0005)
6837 if (std::fabs(a-c) < 0.0005)
6839 r=(-4*(4*(a*a)*(b*b)*
log((a*a)/(b*b)) + ((a+b)*(a+b))*((b*b)*
log((4*(b*b))/((a+b)*(a+b))) + (a*a)*
log(((a+b)*(a+b))/(4.*(a*a))))))/
6840 (
pow(a - b,3)*(a + b)*(3*a + b)*(a + 3*b));
6844 r=(-(a*a) + (c*c) - (a*a)*
log((c*c)/(a*a)))/
pow((a*a) - (c*c),2) +
6845 ((b - c)*(
pow(a,4) -
pow(c,4) + 2*(a*a)*(c*c)*
log((c*c)/(a*a))))/(c*
pow(-(a*a) + (c*c),3));
6850 if (std::fabs(a-c) < 0.0005)
6852 r=(-(b*b) + (c*c) + (b*b)*
log((b*b)/(c*c)))/
pow((b*b) - (c*c),2) +
6853 ((a - c)*(
pow(b,4) -
pow(c,4) - 2*(b*b)*(c*c)*
log((b*b)/(c*c))))/(c*
pow(-(b*b) + (c*c),3));
6857 r=( a*a*b*b*
log(a*a/(b*b))+b*b*c*c*
log(b*b/(c*c))+c*c*a*a*
log(c*c/(a*a)) )
6858 /( (a*a-b*b)*(b*b-c*c)*(a*a-c*c) );
6875 if (std::fabs(x-1.) < 0.005)
6877 r=1. + (1059.*(x-1.))/1175.;
6881 r=(1. - x) * (151. * x*x - 335. * x + 592.)
6882 + 6. * (21. * x*x*x - 108. * x*x - 93. * x + 50.) *
log(x)
6883 - 54. * x * (x*x - 2. * x - 2.) *
pow(
log(x),2)
6885 r*=4./(141.*
pow(1. - x,4));
6893 if (std::fabs(x-1.) < 0.005)
6895 r=1. - (45.*(x-1.))/122.;
6899 r=8.*(x*x-3.*x+2.)+(11.*x*x-40.*x+5.)*
log(x)
6900 -2.*(x*x-2.*x-2.)*
log(x)*
log(x)
6902 r=r*(-9.)/(122.*
pow(1-x,3));
6910 if (std::fabs(x-1.) < 0.005)
6912 r=1. + (76.*(x-1.))/875.;
6916 r=(1-x)*(-97.*x*x-529.*x+2.)+6.*x*x*(13.*x+81.)*
log(x)
6918 r=r*4./(105.*
pow(1-x,4));
6926 if (std::fabs(x-1.) < 0.005)
6928 r=1. - (111.*(x-1.)*(x-1.))/800.;
6933 r=r*(-9.)/(4.*
pow(1-x,3));
6941 if (std::fabs(x-y) < 0.005)
6943 if (std::fabs(y-1.) < 0.005)
6945 r=(49. - 22.*x - 22.*y + 10.*x*y)/60.;
6949 r=(2. + 3.*x - 6.*x*x + x*x*x + 6.*x*
log(x))/(4.*
pow(x-1.,4)*x)
6950 +(2. + 3.*y - 6.*y*y + y*y*y + 6.*y*
log(y))/(4.*
pow(y-1.,4)*y);
6955 if (std::fabs(x-1.) < 0.005)
6957 r=((x-1.)*(-25. + 48.*y - 36.*y*y + 16.*y*y*y - 3.*y*y*y*y - 12.*
log(y)))/(12.*
pow(y-1.,5))
6958 + (-11. + 18.*y - 9.*y*y + 2.*y*y*y - 6.*
log(y))/(6.*
pow(y-1.,4));
6960 else if (std::fabs(y-1.) < 0.005)
6962 r=((y-1.)*(-25. + 48.*x - 36.*x*x + 16.*x*x*x - 3.*x*x*x*x - 12.*
log(x)))/(12.*
pow(x-1.,5))
6963 + (-11. + 18.*x - 9.*x*x + 2.*x*x*x - 6.*
log(x))/(6.*
pow(x-1.,4));
6968 G3x=(1./(2.*
pow(x-1,3))) *( (x-1)*(x-3)+2.*
log(x) );
6969 G3y=(1./(2.*
pow(y-1,3))) *( (y-1)*(y-3)+2.*
log(y) );
6979 if (std::fabs(x-y) < 0.005)
6981 if (std::fabs(y-1.) < 0.005)
6983 r=(13. - 5.*x - 5.*y + 2.*x*y)/60.;
6987 r=((2. + 3.*x - 6.*x*x + x*x*x + 6.*x*
log(x))/
pow(x-1.,4)
6988 + (2. + 3.*y - 6.*y*y + y*y*y + 6.*y*
log(y))/
pow(y-1.,4))/12.;
6993 if (std::fabs(x-1.) < 0.005)
6995 r=(2. + 3.*y - 6.*y*y + y*y*y + 6.*y*
log(y))/(6.*
pow(y-1.,4))
6996 + ((x-1.)*(3. + 10.*y - 18.*y*y + 6.*y*y*y - y*y*y*y + 12.*y*
log(y)))/(12.*
pow(y-1.,5));
6998 else if (std::fabs(y-1.) < 0.005)
7000 r=(2. + 3.*x - 6.*x*x + x*x*x + 6.*x*
log(x))/(6.*
pow(x-1.,4))
7001 + ((y-1.)*(3. + 10.*x - 18.*x*x + 6.*x*x*x - x*x*x*x + 12.*x*
log(x)))/(12.*
pow(x-1.,5));
7006 G4x=(1./(2.*
pow(x-1,3))) *( (x-1)*(x+1) -2.*x*
log(x) );
7007 G4y=(1./(2.*
pow(y-1,3))) *( (y-1)*(y+1) -2.*y*
log(y) );
7022 double g1atMZ = 0.357456;
7023 double g2atMZ = 0.651721;
7032 double vew =
v/
sqrt(2.);
7059 double msq1L = MsQhat2(0,0).real();
7060 double msq2L = MsQhat2(1,1).real();
7061 double msq3L = MsQhat2(2,2).real();
7062 double msuR = MsUhat2(0,0).real();
7063 double mscR = MsUhat2(1,1).real();
7064 double mstR = MsUhat2(2,2).real();
7065 double msdR = MsDhat2(0,0).real();
7066 double mssR = MsDhat2(1,1).real();
7067 double msbR = MsDhat2(2,2).real();
7068 double mseL = MsLhat2(0,0).real();
7069 double msmuL = MsLhat2(1,1).real();
7070 double mstauL = MsLhat2(2,2).real();
7071 double mseR = MsEhat2(0,0).real();
7072 double msmuR = MsEhat2(1,1).real();
7073 double mstauR = MsEhat2(2,2).real();
7080 b1[0]=41./6., b1[1]=-19./6., b1[2]=-7.;
7081 gi[0]=g1atMZ, gi[1]=g2atMZ, gi[2]=g3atMZ;
7083 res_g=1./gi[k]/gi[k] -b1[k]/(8.*pi*pi)*
log(
sqrt(msmuL)/MZ);
7084 gi[k]=
sqrt(1/res_g);
7086 g1=gi[0], g2=gi[1], g3=gi[2];
7088 double alp = (g1*g1*g2*g2/(g1*g1+g2*g2))/(4.0*pi);
7091 double mzq =
sqrt( 0.5*(g1*g1+g2*g2)*vew*vew );
7092 double mwq =
sqrt( 0.5*(g2*g2)*vew*vew );
7093 double sw2 = g1*g1/(g1*g1+g2*g2);
7097 double msneu2 = msmuL+0.5*mzq*mzq*c2b;
7099 TEhat.assign(1,1, 0.);
7107 Rsmu.assign(0,0, msmuL+mmu*mmu+c2b*mzq*mzq*(-1.0/2.0+sw2) );
7110 Rsmu.assign(1,1, msmuR+mmu*mmu-c2b*mzq*mzq*sw2);
7111 Rsmu.eigensystem(Xm,msmu2);
7119 double gm21loop = 1.06978e-9;
7126 double gm2cor, res2, res3;
7127 double res01, res02;
7141 for(
int i=0; i<2; i++){
7143 lxh1.assign(i, tmp1*(
myU(i,0)*
myV(i,1)*ca +
myU(i,1)*
myV(i,0)*(-sa) ) );
7144 lxh2.assign(i, tmp1*(
myU(i,0)*
myV(i,1)*sa +
myU(i,1)*
myV(i,0)*( ca) ) );
7149 double xps, xh1, xh2;
7150 for(
int i=0; i<2; i++){
7155 res += (lA*lxA(i)).real()*
fPS(xps) + (lh1*lxh1(i)).real()*
fS(xh1)
7156 +(lh2*lxh2(i)).real()*
fS(xh2);
7158 res *= alp*alp * mmu*mmu / (8.*pi*pi*mwq*mwq*sw2);
7166 stauM.assign(0,0, mstauL+mTAU*mTAU + mzq*mzq*c2b*(0.5-(-1.)*sw2) );
7167 stauM.assign(1,1, mstauR+mTAU*mTAU + mzq*mzq*c2b*(-1.)*sw2);
7168 stauM.assign(0,1, mTAU*(a3tau-muH*
tanb));
7169 stauM.assign(1,0, mTAU*(a3tau-muH*
tanb));
7171 sbottomM.assign(0,0, msq3L+mb*mb + mzq*mzq*c2b*(0.5-(-1./3.)*sw2) );
7172 sbottomM.assign(1,1, msbR+mb*mb + mzq*mzq*c2b*(-1./3.)*sw2);
7173 sbottomM.assign(0,1, mb*(a3b-muH*
tanb));
7174 sbottomM.assign(1,0, mb*(a3b-muH*
tanb));
7176 stopM.assign(0,0, msq3L+mt*mt + mzq*mzq*c2b*(0.5-(2./3.)*sw2) );
7177 stopM.assign(1,1, mstR+mt*mt + mzq*mzq*c2b*(2./3.)*sw2);
7178 stopM.assign(0,1, mt*(a3t-muH/
tanb));
7179 stopM.assign(1,0, mt*(a3t-muH/
tanb));
7181 stauM.eigensystem(Ustau,mstau2);
7182 sbottomM.eigensystem(Usbottom,msbottom2);
7183 stopM.eigensystem(Ustop,mstop2);
7185 Ustau=Ustau.hconjugate();
7186 Usbottom=Usbottom.hconjugate();
7187 Ustop=Ustop.hconjugate();
7197 for(
int i=0; i<2; i++){
7199 rr=2.*mTAU/(
cosb*mstau2(i)) * Ustau(i,0).
conjugate() * Ustau(i,1);
7200 lstauh1.assign(i,rr*(-muH*ca + a3tau*(-sa)) );
7201 lstauh2.assign(i,rr*(-muH*sa + a3tau*(ca)) );
7203 rr=2.*mb/(
cosb*msbottom2(i)) * Usbottom(i,0).
conjugate() * Usbottom(i,1);
7204 lsbottomh1.assign(i,rr*(-muH*ca + a3b*(-sa)) );
7205 lsbottomh2.assign(i,rr*(-muH*sa + a3b*(ca)) );
7207 rr=2.*mt/(
cosb*mstop2(i)) * Ustop(i,0).conjugate() * Ustop(i,1);
7208 lstoph1.assign(i,rr*(-muH*sa + a3t*(ca)) );
7209 lstoph2.assign(i,rr*(-muH*(-ca) + a3t*(sa)) );
7218 for(
int i=0;i<2;i++){
7219 xx1=mstau2(i)/(mh*mh);
7220 xx2=mstau2(i)/(mhh*mhh);
7221 res1 += (lh1*lstauh1(i)).real()*
fft(xx1)
7222 +(lh2*lstauh2(i)).real()*
fft(xx2);
7226 xx1=msbottom2(i)/(mh*mh);
7227 xx2=msbottom2(i)/(mhh*mhh);
7228 res1 += (lh1*lsbottomh1(i)).real()*qe2*
fft(xx1)
7229 +(lh2*lsbottomh2(i)).real()*qe2*
fft(xx2);
7233 xx1=mstop2(i)/(mh*mh);
7234 xx2=mstop2(i)/(mhh*mhh);
7235 res1 += (lh1*lstoph1(i)).real()*qe2*
fft(xx1)
7236 +(lh2*lstoph2(i)).real()*qe2*
fft(xx2);
7239 res1 *= alp*alp * mmu*mmu / (8.*pi*pi*mwq*mwq*sw2);
7251 gm2cor=0, res2=0, res3=0;
7262 double x0, x1a, x2a, xL, xR;
7265 x0=
sqrt(std::fabs(msneu2));
7266 xL=
sqrt( msmuL-mzq*mzq*(sw2-0.5) );
7267 xR=
sqrt( msmuR +mzq*mzq*sw2 );
7269 tmp2=M2.
abs2()+muH.
abs2() +2.*mwq*mwq;
7270 tmp3= tmp2*tmp2 -4.*M2.
abs2()*muH.
abs2();
7271 x1a=
sqrt( 0.5*( tmp2-
sqrt(tmp3) ) );
7272 x2a=
sqrt( 0.5*( tmp2+
sqrt(tmp3) ) );
7274 dmu=-muH.
real()*
tanb*g2*g2*M2.
real()/(16.*pi*pi)*(
It(x1a,x2a,x0) + 0.5*
It(x1a,x2a,xL))
7278 gm2cor=gm21loop/(1+dmu);
7287 res2=gm2cor*alp/(4.*pi)*16.*
log(mmu/
sqrt(msmuL));
7296 if (sub_leading==1){
7301 double ymu=mmu/(vew*
cosb);
7303 ckR.assign(0, ymu*
myU(0,1));
7304 ckR.assign(1, ymu*
myU(1,1));
7307 ckL.assign(0, -g2*
myV(0,0));
7308 ckL.assign(1, -g2*
myV(1,0));
7311 for(
int i=0; i<2; i++){
7313 amch=amch -(47.*mmu/(12.*msneu2))*( (ckL(i).abs()*ckL(i).abs() + ckR(i).abs()*ckR(i).abs())*
F3C(xk) );
7314 amch=amch -(122.*
MChi(i)/(9.*msneu2))*
F4C(xk)*( (ckL(i)*ckR(i)).real() );
7316 amch=(1./(1+dmu))*amch*mmu/(16.*pi*pi)*alp/(4.*pi);
7323 for(
int i=0; i<4; i++){
7324 nR1.assign(i,
sqrt(2.)*g1*
myN(i,0)*Xm(0,1) + ymu*
myN(i,2)*Xm(0,0));
7325 nR2.assign(i,
sqrt(2.)*g1*
myN(i,0)*Xm(1,1) + ymu*
myN(i,2)*Xm(1,0));
7326 nL1.assign(i, (1./
sqrt(2.))*(g2*
myN(i,1) + g1*
myN(i,0))*Xm(0,0).conjugate() -ymu*
myN(i,2)*Xm(0,1).conjugate());
7327 nL2.assign(i, (1./
sqrt(2.))*(g2*
myN(i,1) + g1*
myN(i,0))*Xm(1,0).conjugate() -ymu*
myN(i,2)*Xm(1,1).conjugate());
7333 for(
int i=0; i<4; i++){
7336 tmp4=nL1(i).abs2() + nR1(i).abs2();
7337 tmp5=nL2(i).abs2() + nR2(i).abs2();
7338 r1=35.*mmu/(72.*msmu2(0))*
F3N(xi1)*tmp4 - 16.*
MChi0(i)/(9.*msmu2(0))*
F4N(xi1)*( (nL1(i)*nR1(i)).real() );
7339 r2=35.*mmu/(72.*msmu2(1))*
F3N(xi2)*tmp5 - 16.*
MChi0(i)/(9.*msmu2(1))*
F4N(xi2)*( (nL2(i)*nR2(i)).real() );
7343 amne*=(1./(1.+dmu))*mmu/(16.*pi*pi)*alp/(4.*pi);
7347 res2=res2+amch+amne;
7360 double x1,y1,x2,y2,x3,y3,x4,y4,x5,y5;
7361 double awhn=mmu*mmu*muH.
real()*
tanb/(1+dmu);
7368 x1=M2.
abs2()/msneu2;
7369 y1=muH.
abs2()/msneu2;
7370 awhn *= g2*g2*M2.
real()/(8.*pi*pi*msneu2*msneu2)*
Fa(x1,y1);
7374 y2=muH.
abs2()/msmuL;
7375 awhl *= -g2*g2*M2.
real()/(16.*pi*pi*
pow(msmuL,2))*
Fb(x2,y2);
7379 y3=muH.
abs2()/msmuL;
7380 abhl *= g1*g1*M1.
real()/(16.*pi*pi*
pow(msmuL,2))*
Fb(x3,y3);
7384 y4=muH.
abs2()/msmuR;
7385 abhr *= -g1*g1*M1.
real()/(8.*pi*pi*
pow(msmuR,2))*
Fb(x4,y4);
7388 x5=msmuL/(M1.
abs2());
7389 y5=msmuR/(M1.
abs2());
7390 ablr *= g1*g1/(8.*pi*pi*
pow(M1.
real(),3))*
Fb(x5,y5);
7392 double dg2, dg1, dh, dwh, dbh, dtb;
7393 double msusy=
sqrt(msmuL);
7395 dg1=g1*g1/(16.*pi*pi)*(4./3.)
7402 dg2=g2*g2/(16.*pi*pi)*(4./3.)
7407 double yb, ytau, yt;
7408 double as_mt, delta_mt;
7413 as_mt=as_mt - (-7.)*
log(mt/
sqrt(msmuL))/(8.*pi*pi);
7414 as_mt=(1./as_mt)/(4.*pi);
7416 delta_mt=-4./3.*(as_mt/pi)-9.1*
pow( (as_mt/pi), 2)-80.*
pow((as_mt/pi),3);
7419 yt=mt/(vew*
sinb)*(1+delta_mt);
7420 ytau=mTAU/(vew*
cosb);
7423 yt=yt*(1+as_mt/(8.*pi)*(-4./3.));
7424 yb=yb*(1+as_mt/(8.*pi)*(-4./3.));
7426 dh=0.5/(16.*pi*pi)*(3.*yt*yt*
log(
sqrt(mstR)/msusy) +3.*yb*yb*
log(
sqrt(msbR)/msusy)
7427 +3.*(yt*yt+yb*yb)*
log(
sqrt(msq3L)/msusy)
7428 +ytau*ytau*
log(
sqrt(mstauR)/msusy) +ytau*ytau*
log(
sqrt(mstauL)/msusy));
7430 dwh=yt*yt/(16.*pi*pi)*(-6.*
log(
sqrt(msq3L)/msusy));
7432 dbh=yt*yt/(16.*pi*pi)*( 2.*
log(
sqrt(msq3L)/msusy)-8.*
log(
sqrt(mstR)/msusy) );
7436 dtb=1./(16.*pi*pi)*( 3.*yb*yb -3.*yt*yt +ytau*ytau)*
log(
Q_S/msusy);
7439 res3= awhn*(dg2 + dh + dwh + dtb)
7440 + awhl*(dg2 + dh + dwh + dtb)
7441 + abhl*(dg1 + dh + dbh + dtb)
7442 + abhr*(dg1 + dh + dbh + dtb)
7452 res01=gm2cor+res2+res3-gm21loop;
7477 double sw2 =
mySUSY.StandardModel::sW2(MW);
7478 double cw2 = 1.0 - sw2;
7484 + (1.0-1.0/(4.0*sw2))/(cw2*MZ*MZ)*
FFunctions(li_to_lj)(1)
7487 + (1.0-1.0/(4.0*sw2))/(cw2*MZ*MZ)*
FFunctions(li_to_lj)(0)
7497 double sw2 =
mySUSY.StandardModel::sW2(MW);
7498 double cw2 = 1.0 - sw2;
7503 C10.assign(0, 1.0/(4.0*sw2*cw2*MZ*MZ)*
FFunctions(li_to_lj)(1)
7505 C10.assign(1, 1.0/(4.0*sw2*cw2*MZ*MZ)*
FFunctions(li_to_lj)(0)
7546 vmcDLij = StandardModelMatching::CMDLij(li_lj);
7557 std::stringstream out;
7559 throw std::runtime_error(
"SUSYMatching::CMDLij(): order " + out.str() +
" not implemented.\nFor lepton flavour violating observables only Leading Order (LO) necessary.");
7562 vmcDLij.push_back(
mcDLij);
7569 vmcDLi3j = StandardModelMatching::CMDLi3j(li_lj);
7613 std::stringstream out;
7615 throw std::runtime_error(
"SUSYMatching::CMDLi3j(): order " + out.str() +
" not implemented.\nFor lepton flavour violating observables only Leading Order (LO) necessary.");
7626 double sw2 =
mySUSY.StandardModel::sW2(MW);
7628 vmcmueconv = StandardModelMatching::CMmueconv();
7647 std::stringstream out;
7649 throw std::runtime_error(
"SUSYMatching::CMmueconv(): order " + out.str() +
" not implemented.\nFor lepton flavour violating observables only Leading Order (LO) necessary.");
7659 vmcgminus2mu = StandardModelMatching::CMgminus2mu();
7675 std::stringstream out;
7677 throw std::runtime_error(
"SUSYMatching::CMgminus2mu(): order " + out.str() +
" not implemented.\nFor lepton flavour violating observables only Leading Order (LO) necessary.");
7681 return(vmcgminus2mu);