13 CW0b(2,0.),CW1b(2,0.),CW2b(2,0.),CWeb(2,0.),CWesb(2,0.),
14 U0_4b(2,0.),U0_5b(2,0.),J1_4b(2,0.),J2_4b(2,0.),J1_5b(2,0.),J2_5b(2,0.),
15 R0_4b(2,0.),R0_5b(2,0.),R1_4b(2,0.),R1_5b(2,0.),
16 CW0p(3,0.),CW1p(3,0.),CW2p(3,0.),CWep(3,0.),CWesp(3,0.),
17 U0_4p(3,0.),U0_5p(3,0.),J1_4p(3,0.),J2_4p(3,0.),J1_5p(3,0.),J2_5p(3,0.),
18 R0_4p(3,0.),R0_5p(3,0.),R1_4p(3,0.),R1_5p(3,0.)
30 xi2c=966966391. / 10546875. - 231404944. / 3515625. * (1. /
etac) - 272751559. / 10546875. *
31 (1. /
etac)*(1. /
etac) - 128. / 5. * (1. - (1. /
etac)*(1. /
etac)) * gsl_sf_zeta_int(3);;
32 xice= 8. / 3. / (11. - 2. / 3. * 3. ) * (
etac - 1.);
33 xices=(32. / 9. / (11. - 2. / 3. * 3. ) - (-8. / 9. * (1. + 2./4. ) ) * 8. / (11. - 2. / 3. * 3. ) / (11. - 2. / 3. * 3. ) - (102. - 38. / 3. * 3. ) * 8. / 3. / (11. - 2. / 3. * 3. ) / (11. - 2. / 3. * 3. )) *
34 log(
etac) + 8. / 3. / (11. - 2. / 3. * 3. ) * (8. * (102. - 38. / 3. * 3. ) / (11. - 2. / 3. * 3. ) / (11. - 2. / 3. * 3. ) - (404. / 3. - 40. / 9. * 3.) / (11. - 2. / 3. * 3. )) * (1. - 1. /
etac) * (1. -
etac);
618 switch (contribution)
630 double gamma_nuB = -8. ;
631 double gamma_M = 8. ;
632 double beta = 11. - 2. / 3. * nf ;
633 double gamma_nu = 2. * ( gamma_M - beta ) ;
635 gslpp::matrix<double> gamma_T(2, 0.);
638 gamma_T(1,0) = gamma_nuB;
639 gamma_T(1,1) = gamma_nu;
649 double gamma_M = 404. / 3. - 40. / 9. * nf ;
650 double beta = 102. - 38. / 3. * nf ;
651 double gamma_nu = 2. * ( gamma_M - beta ) ;
653 gslpp::matrix<double> gamma_T(2, 0.);
656 gamma_T(1,0) = gamma_nuB;
657 gamma_T(1,1) = gamma_nu;
666 double gamma_nuB = 2. *
model.
getCF() * (69. / 3. - 458./3.*3. - (48. / 3. - 96. * 3.) * gsl_sf_zeta_int(3) + 38. / 3. * nf) ;
667 double gamma_M = 2498. - (4432./27. + 320.*gsl_sf_zeta_int(3)/3.) * nf - (280. * nf * nf )/81. ;
668 double beta = 2857./2. - 5033.*nf/18. + 325.*nf*nf/54. ;
669 double gamma_nu = 2. * ( gamma_M - beta ) ;
671 gslpp::matrix<double> gamma_T(2, 0.);
674 gamma_T(1,0) = gamma_nuB;
675 gamma_T(1,1) = gamma_nu;
683 std::stringstream out;
685 throw std::runtime_error(
"Charm_Kpnunu::ADM: order " + out.str() +
" not implemented\n");
691 double gamma_W = -4. ;
692 double gamma_nuB = 0. ;
693 double gamma_M = 8. / 3. ;
695 double gamma_nu = 2. * ( gamma_M - beta ) ;
697 gslpp::matrix<double> gamma_T(2, 0.);
699 gamma_T(0,0) = 2. * gamma_W;
700 gamma_T(1,0) = gamma_nuB;
701 gamma_T(1,1) = gamma_nu;
709 double gamma_W = 4. ;
710 double gamma_nuB = -316. / 9. ;
711 double gamma_M = 32. / 9. ;
712 double beta = -8. / 9. * (nfu + nfd/4. ) ;
713 double gamma_nu = 2. * ( gamma_M - beta ) ;
715 gslpp::matrix<double> gamma_T(2, 0.);
717 gamma_T(0,0) = 2. * gamma_W;
718 gamma_T(1,0) = gamma_nuB;
719 gamma_T(1,1) = gamma_nu;
726 std::stringstream out, out2;
728 out2 << contribution;
729 throw std::runtime_error(
"Charm_Kpnunu::ADM: order_qed " + out.str() +
" of contribution " + out2.str() +
" not implemented\n");
742 gslpp::matrix<double> gamma_pm(2,0.);
743 gamma_pm(0,0)= +6.*(1.-1./3.);
746 gamma_pm(1,1)= -6.*(1.+1./3.);
748 double gamma_pnu = -0.5 * ( -4. * ( 1. + 3. ) ) ;
749 double gamma_mnu = -0.5 * ( -4. * ( 1. - 3. ) ) ;
750 double gamma_M = 8. ;
751 double beta = 11. - 2. / 3. * nf ;
752 double gamma_nu = 2. * ( gamma_M - beta ) ;
754 gslpp::matrix<double> gamma_T(3, 0.);
756 gamma_T(0,0) = gamma_pm(0,0);
757 gamma_T(1,0) = gamma_pm(0,1);
758 gamma_T(0,1) = gamma_pm(1,0);
759 gamma_T(1,1) = gamma_pm(1,1);
761 gamma_T(2,0) = gamma_pnu;
762 gamma_T(2,1) = gamma_mnu;
764 gamma_T(2,2) = gamma_nu;
772 gslpp::matrix<double> gamma_pm(2,0.);
773 gamma_pm(0,0)= (-21./2. + 2.*nf/3.) * (1. - 1./3.);
776 gamma_pm(1,1)= (-21./2. - 2.*nf/3.) * (1. + 1./3.);
778 double gamma_pnu = -0.5 * ( 16. * (2. - 11.) ) ;
779 double gamma_mnu = -0.5 * ( 16. * (2. + 11.) ) ;
780 double gamma_M = 404. / 3. - 40. / 9. * nf ;
781 double beta = 102. - 38. / 3. * nf ;
782 double gamma_nu = 2. * ( gamma_M - beta ) ;
784 gslpp::matrix<double> gamma_T(3, 0.);
786 gamma_T(0,0) = gamma_pm(0,0);
787 gamma_T(1,0) = gamma_pm(0,1);
788 gamma_T(0,1) = gamma_pm(1,0);
789 gamma_T(1,1) = gamma_pm(1,1);
791 gamma_T(2,0) = gamma_pnu;
792 gamma_T(2,1) = gamma_mnu;
794 gamma_T(2,2) = gamma_nu;
802 gslpp::matrix<double> gamma_pm(2,0.);
803 gamma_pm(0,0)= 1./300. * (349049. + 201485.) - 1./1350. * (115577. - 9795.)*nf - 130./27. * (1. - 1./3.)*nf*nf - (672. + 80. * (1. - 1./3.)*nf ) * gsl_sf_zeta_int(3);
806 gamma_pm(1,1)= 1./300. * (349049. - 201485.) - 1./1350. * (115577. + 9795.)*nf + 130./27. * (1. + 1./3.)*nf*nf + (672. + 80. * (1. + 1./3.)*nf ) * gsl_sf_zeta_int(3);
808 double gamma_pnu = -0.5 * ( -2./225. * (45124. + 484917.) + 32.*(13. + 15.)*gsl_sf_zeta_int(3) + 144.*nf ) ;
809 double gamma_mnu = -0.5 * ( -2./225. * (45124. - 484917.) + 32.*(13. - 15.)*gsl_sf_zeta_int(3) - 144.*nf );
810 double gamma_M = 2498. - (4432./27. + 320.*gsl_sf_zeta_int(3)/3.) * nf - (280. * nf * nf )/81. ;
811 double beta = 2857./2. - 5033.*nf/18. + 325.*nf*nf/54. ;
812 double gamma_nu = 2. * ( gamma_M - beta ) ; ;
814 gslpp::matrix<double> gamma_T(3, 0.);
816 gamma_T(0,0) = gamma_pm(0,0);
817 gamma_T(1,0) = gamma_pm(0,1);
818 gamma_T(0,1) = gamma_pm(1,0);
819 gamma_T(1,1) = gamma_pm(1,1);
821 gamma_T(2,0) = gamma_pnu;
822 gamma_T(2,1) = gamma_mnu;
824 gamma_T(2,2) = gamma_nu;
831 std::stringstream out;
833 throw std::runtime_error(
"Charm_Kpnunu::ADM: order " + out.str() +
" not implemented\n");
839 gslpp::matrix<double> gamma_pm(2,0.);
840 gamma_pm(0,0)= -8./3.;
843 gamma_pm(1,1)= -8./3.;
845 double gamma_pnu = 0. ;
846 double gamma_mnu = 0. ;
847 double gamma_M = 8. / 3. ;
849 double gamma_nu = 2. * ( gamma_M - beta ) ;
851 gslpp::matrix<double> gamma_T(3, 0.);
853 gamma_T(0,0) = gamma_pm(0,0);
854 gamma_T(1,0) = gamma_pm(0,1);
855 gamma_T(0,1) = gamma_pm(1,0);
856 gamma_T(1,1) = gamma_pm(1,1);
858 gamma_T(2,0) = gamma_pnu;
859 gamma_T(2,1) = gamma_mnu;
861 gamma_T(2,2) = gamma_nu;
869 gslpp::matrix<double> gamma_12(2,0.);
870 gslpp::matrix<complex> trash(2,0.);
871 gslpp::vector<complex> ev_gamma12(2,0.);
873 gamma_12.assign(0,0,194./9.);
874 gamma_12.assign(0,1,-2./3.);
875 gamma_12.assign(1,0,25./3.);
876 gamma_12.assign(1,1,-49./9.) ;
878 gamma_12.eigensystem(trash,ev_gamma12);
880 double gamma_pnu = 0. ;
881 double gamma_mnu = 0. ;
882 double gamma_M = 32. / 9. ;
883 double beta = -8. / 9. * (nfu + nfd/4. ) ;
884 double gamma_nu = 2. * ( gamma_M - beta ) ;
886 gslpp::matrix<double> gamma_T(3, 0.);
888 gamma_T(0,0) = ev_gamma12(0).real();
889 gamma_T(1,1) = ev_gamma12(1).real();
891 gamma_T(2,0) = gamma_pnu;
892 gamma_T(2,1) = gamma_mnu;
894 gamma_T(2,2) = gamma_nu;
901 std::stringstream out, out2;
903 out2 << contribution;
904 throw std::runtime_error(
"Charm_Kpnunu::ADM: order_qed " + out.str() +
" of contribution " + out2.str() +
" not implemented\n");
910 std::stringstream out ;
912 throw std::runtime_error(
"Charm_Kpnunu::ADM: contribution " + out.str() +
" not implemented\n");
921 switch (contribution)
934 gslpp::vector<double> CWin(2, 0.);
949 gslpp::vector<double> CWin(2, 0.);
950 double Cnu= -4. * (9. + 4. * LW) ;
964 gslpp::vector<double> CWin(2, 0.);
965 double Cnu= -8. *
model.
getCF()*(20. + 2. * M_PI * M_PI + 25. * LW + 6. * LW * LW) ;
976 std::stringstream out, out2;
978 out2 << contribution;
979 throw std::runtime_error(
"Charm_Kpnunu::CWin_muw: order " + out.str() +
" of contribution " + out2.str() +
" not implemented\n");
986 gslpp::vector<double> CWin(2, 0.);
990 CWin(0)= 4. * CW * CW ;
1001 gslpp::vector<double> CWin(2, 0.);
1002 double CW= -11. / 3. - 2. * LZ ;
1005 CWin(0)= 4. * CW * CW ;
1012 std::stringstream out, out2;
1014 out2 << contribution;
1015 throw std::runtime_error(
"Charm_Kpnunu::Cwin_qed: order_qed " + out.str() +
" of contribution " + out2.str() +
" not implemented\n");
1033 gslpp::vector<double> CWin(3, 0.);
1051 gslpp::vector<double> CWin(3, 0.);
1053 double Cp= 0.5*(1. - 1./3.) * (11. + 6. * LW) ;
1054 double Cm= -0.5*(1. + 1./3.) * (11. + 6. * LW) ;
1055 double Cnu= 8. * (2. + LW) ;
1072 gslpp::vector<double> CWin(3, 0.);
1074 double Cp= (-(135677. - 124095.) / 3600. + 58. / 18. * M_PI * M_PI - 0.5 * (2. / 3.)*
1075 (112. / 9. + 32. * x + (20. / 3. + 16. * x) * log(x) - (8. + 16. * x) *
1076 sqrt(4. * x - 1.) * gsl_sf_clausen(2. * asin(1. / (2. * sqrt(x)))))
1077 +(5. / 36. * 238. * LW) + 58. / 6. * LW * LW) ;
1078 double Cm= (-(135677. + 124095.) / 3600. - 44. / 18. * M_PI * M_PI + 0.5 * (4. / 3.)*
1079 (112. / 9. + 32. * x + (20. / 3. + 16. * x) * log(x) - (8. + 16. * x) *
1080 sqrt(4. * x - 1.) * gsl_sf_clausen(2. * asin(1. / (2. * sqrt(x)))))
1081 - (5. / 36. * 260. * LW) - 44. / 6. * LW * LW);
1082 double Cnu= 4. *
model.
getCF() * (33. + 4. * M_PI * M_PI + 34. * LW + 12. * LW * LW) ;
1095 std::stringstream out, out2;
1097 out2 << contribution;
1098 throw std::runtime_error(
"Charm_Kpnunu::CWin_muw: order " + out.str() +
" of contribution " + out2.str() +
" not implemented\n");
1105 gslpp::vector<double> CWin(3, 0.);
1111 CWin(0)= 4. * Cp * CA ;
1112 CWin(1)= 4. * Cm * CA ;
1127 double mt2 = mt * mt;
1128 double MW2 = MW * MW;
1129 double MZ2 = MZ * MZ;
1130 double MH2 = MH * MH;
1131 double MH4 = MH2 * MH2;
1133 double sw4 = sw2 * sw2;
1134 double cw2 = 1. - sw2;
1136 gslpp::vector<double> CWin(3, 0.);
1137 double Cp= -22. / 9. - 4. / 3. * LZ ;
1138 double Cm= -22. / 9. - 4. / 3. * LZ ;
1139 double CA= 3. * mt2 / 4. / sw2 / MW2 + (11. * sw2 - 6.) / 4. / sw2 / cw2 -
1140 3. / 4. * (MW2 - cw2 * MH2) / (MH2 - MW2) / sw4 * log(MW2 / MZ2) +
1141 3. * MH4 / 4. / (MH2 - MW2) / (MW2 - cw2 * MH2) * log(MH2 / MZ2);
1144 CWin(0)= 4. * Cp * CA ;
1145 CWin(1)= 4. * Cm * CA ;
1152 std::stringstream out, out2;
1154 out2 << contribution;
1155 throw std::runtime_error(
"Charm_Kpnunu::Cwin_qed: order_qed " + out.str() +
" of contribution " + out2.str() +
" not implemented\n");
1161 std::stringstream out ;
1162 out << contribution;
1163 throw std::runtime_error(
"Charm_Kpnunu::Cwin_qed: contribution " + out.str() +
" not implemented\n");
1175 switch (contribution)
1186 std::stringstream out ;
1187 out << contribution;
1188 throw std::runtime_error(
"Charm_Kpnunu::RGevol_J : contribution " + out.str() +
" not implemented\n");
1193 gslpp::matrix<gslpp::complex> V(dimension,0.) , V_inv(dimension,0.) ;
1194 gslpp::matrix<double> G1(dimension,0.) ;
1195 gslpp::matrix<double> S1(dimension,0.) ;
1196 gslpp::matrix<double> J1(dimension,0.) ;
1197 gslpp::matrix<double> ADM_0(dimension,0.) , ADM_1(dimension,0.) , ADM_2(dimension,0.);
1200 gslpp::vector<gslpp::complex> e(dimension,0.) ;
1202 double beta_0 = 11. - 2. / 3. * nf ;
1203 double beta_1 = 102. - 38. / 3. * nf ;
1204 double beta_2 = 2857. / 2. - 5033. / 18. * nf + 325. / 54. * nf * nf ;
1211 ADM_0.eigensystem(V,e);
1215 G1 = V_inv.real() * ADM_1 * V.real() ;
1218 for(
unsigned int i = 0; i < G1.size_i(); i++){
1219 for (
unsigned int j = 0; j < G1.size_j(); j++) {
1220 S1.assign(i , j, beta_1 / 2. / beta_0 / beta_0 * e(i).real() * (i==j) - G1(i,j) / (2. * beta_0 + e(i).real() - e(j).real() ));
1224 J1 = V.real() * S1 * V_inv.real() ;
1227 switch (contribution){
1233 std::stringstream out ;
1234 out << contribution;
1235 throw std::runtime_error(
"Charm_Kpnunu::RGevol_J : contribution " + out.str() +
" not implemented\n");
1237 }
else if (order==
NLO){
1239 }
else if (order==
NNLO){
1240 gslpp::matrix<double> S2(dimension,0.) , G2(dimension,0.) , J2(dimension,0.);
1241 G2 = V_inv.real() * ADM_2 * V.real() ;
1242 for(
unsigned int i = 0; i < G1.size_i(); i++){
1243 for (
unsigned int j = 0; j < G1.size_j(); j++) {
1245 for (
unsigned int k = 0; k < G1.size_j(); k++){
1246 term += ( 1. + e(i).real() / 2. / beta_0 - e(k).real() / 2. / beta_0 ) / ( 2. + e(i).real() / 2. / beta_0 - e(j).real() / 2. / beta_0 ) * ( S1(i,k) * S1(k,j) - beta_1 / beta_0 * S1(i,j) * (j==k) );
1248 S2.assign(i , j , beta_2 / 4. / beta_0 / beta_0 * e(i).real() * (i==j) + term - G2(i,j) / (4. * beta_0 + e(i).real() - e(j).real() ));
1251 J2 = V.real() * S2 * V_inv.real() ;
1255 std::stringstream out ;
1257 throw std::runtime_error(
"Charm_Kpnunu::RGevol_J : order " + out.str() +
" not implemented\n");
1266 switch (contribution)
1277 std::stringstream out ;
1278 out << contribution;
1279 throw std::runtime_error(
"Charm_Kpnunu::RGevol_K : contribution " + out.str() +
" not implemented\n");
1284 gslpp::matrix<gslpp::complex> V(dimension,0.) , V_inv(dimension,0.) ;
1285 gslpp::matrix<double> M_0(dimension,0.) , K_0(dimension,0.) ;
1286 gslpp::matrix<double> ADM_0(dimension,0.) , ADM_1(dimension,0.) , ADM_2(dimension,0.), ADM_e(dimension,0.), ADM_es(dimension,0.);
1287 gslpp::matrix<double> M_1(dimension,0.), S_1(dimension,0.) , K1_1(dimension,0.) , K2_1(dimension,0.) , K3_1(dimension,0.), G_1(dimension,0.) ;
1290 gslpp::vector<gslpp::complex> e(dimension,0.) ;
1292 double beta_0 = 11. - 2. / 3. * nf ;
1293 double beta_1 = 102. - 38. / 3. * nf ;
1303 ADM_0.eigensystem(V,e);
1307 M_0 = V_inv.real() * ADM_e * V.real();
1308 M_1 = V_inv.real() * ( (ADM_es - beta_1 / beta_0 * ADM_e) + (ADM_e *
RGevol_J(
NLO,nf,contribution) -
RGevol_J(
NLO,nf,contribution) * ADM_e) ) * V.real();
1309 G_1 = V_inv.real() * ADM_1 * V.real();
1311 for(
unsigned int i = 0; i < G_1.size_i(); i++){
1312 for (
unsigned int j = 0; j < G_1.size_j(); j++) {
1313 S_1.assign(i , j, beta_1 / 2. / beta_0 / beta_0 * e(i).real() * (i==j) - G_1(i,j) / (2. * beta_0 + e(i).real() - e(j).real() ));
1318 for(
unsigned int i = 0; i < M_0.size_i(); i++){
1319 for (
unsigned int j = 0; j < M_0.size_j(); j++) {
1320 if (e(i).real() != e(j).real() + 2. * beta_0){
1321 K_0.assign(i , j, M_0(i,j) / (e(i).real() / 2. / beta_0 - e(j).real() / 2. / beta_0 - 1.) * ( pow(
etacb,e(j).real() /2. / beta_0) - pow(
etacb,e(i).real() /2. / beta_0) /
etacb) );
1323 K_0.assign(i , j, M_0(i,j) * pow(
etacb, e(j).real() / 2. / beta_0) * log( 1. /
etacb) );
1329 for(
unsigned int i = 0; i < M_0.size_i(); i++){
1330 for (
unsigned int j = 0; j < M_0.size_j(); j++) {
1331 if (e(i).real() != e(j).real() + 2. * beta_0){
1332 K_0.assign(i , j, M_0(i,j) / (e(i).real() / 2. / beta_0 - e(j).real() / 2. / beta_0 - 1.) * ( pow(
etab,e(j).real() /2. / beta_0) /
etacb - pow(
etab,e(i).real() /2. / beta_0) /
etacb /
etab) );
1334 K_0.assign(i , j, M_0(i,j) * pow(
etab, e(j).real() / 2. / beta_0) * log( 1. /
etab) /
etacb );
1340 std::stringstream out ;
1342 throw std::runtime_error(
"Charm_Kpnunu::RGevol_K : nf " + out.str() +
" not implemented\n");
1352 return ((-2. * M_PI / beta_0) * (V.real() * K_0 * V_inv.real()) );
1355 for(
unsigned int i = 0; i < M_1.size_i(); i++){
1356 for (
unsigned int j = 0; j < M_1.size_j(); j++) {
1358 K1_1.assign(i , j, M_1(i,j) / (e(i).real() / 2. / beta_0 - e(j).real() / 2. / beta_0) * ( pow(
etacb,e(j).real() /2. / beta_0) - pow(
etacb,e(i).real() /2. / beta_0) ) );
1360 K1_1.assign(i , j, M_1(i,j) * pow(
etacb, e(i).real() / 2. / beta_0) * log( 1. /
etacb) );
1364 K2_1 = -
etacb * K_0 * S_1;
1367 return ( ( -0.5 / beta_0) * V.real() * (K1_1 + K2_1 + K3_1) * V_inv.real());
1370 std::stringstream out ;
1372 throw std::runtime_error(
"Charm_Kpnunu::RGevol_K : order " + out.str() +
" not implemented\n");
1378 return ((-2. * M_PI / beta_0) * (V.real() * K_0 * V_inv.real()) );
1381 for(
unsigned int i = 0; i < M_1.size_i(); i++){
1382 for (
unsigned int j = 0; j < M_1.size_j(); j++) {
1384 K1_1.assign(i , j, M_1(i,j) / (e(i).real() / 2. / beta_0 - e(j).real() / 2. / beta_0) * ( pow(
etab,e(j).real() /2. / beta_0) - pow(
etab,e(i).real() /2. / beta_0) ) );
1386 K1_1.assign(i , j, M_1(i,j) * pow(
etab, e(i).real() / 2. / beta_0) * log( 1. /
etab) );
1391 K3_1 =
etacb * S_1 * K_0;
1393 return ( ( -0.5 / beta_0) * V.real() * (K1_1 + K2_1 + K3_1) * V_inv.real());
1396 std::stringstream out ;
1398 throw std::runtime_error(
"Charm_Kpnunu::RGevol_K : order " + out.str() +
" not implemented\n");
1403 std::stringstream out ;
1405 throw std::runtime_error(
"Charm_Kpnunu::RGevol_K : nf " + out.str() +
" not implemented\n");