a Code for the Combination of Indirect and Direct Constraints on High Energy Physics Models Logo
CPenguinBox.cpp
Go to the documentation of this file.
1/*
2 * Copyright (C) 2012 HEPfit Collaboration
3 *
4 *
5 * For the licensing terms see doc/COPYING.
6 */
7
8#include "CPenguinBox.h"
9#include "StandardModel.h"
10#include "std_make_vector.h"
11
12CPenguinBox::CPenguinBox(const StandardModel& model_i)
13: model(model_i),
14up(3, 0.), j5p1(3, 0.), j4p1(3, 0.), j5p2(3, 0.), j4p2(1, 0.),
15ub(2, 0.), j5b1(2, 0.), j4b1(2, 0.), j5b2(2, 0.), j4b2(2, 0.),
16up4(3, 0.), up5(3, 0.), ub4(2, 0.), ub5(2, 0.),
17p(3, 0.), b(2, 0.),
18p0(3, 0.), p1(3, 0.), p2(3, 0.), b0(2, 0.), b1(2, 0.), b2(2, 0.),
19r1(3, 0.), r2(3, 0.), CPM0(3, 0.), CPM1(3, 0.), CPM2(3, 0.),
20CBM0(3, 0.), CBM1(3, 0.), CBM2(3, 0.)
21{
22
23 j5p1(0, 0) = 5165. / 3174.;
24 j5p1(1, 1) = -2267. / 1587.;
25 j5p1(2, 2) = -14924. / 1587.;
26 j5p1(2, 0) = -15857. / 1587.;
27 j5p1(2, 1) = 15305. / 3174.;
28
29 j5p2(0, 0) = -7.35665;
30 j5p2(1, 1) = -54.9107;
31 j5p2(2, 2) = 18.3025;
32 j5p2(2, 0) = 17.7699;
33 j5p2(2, 1) = -1.7514;
34
35 j4p1(0, 0) = 6719. / 3750.;
36 j4p1(1, 1) = -3569. / 1875.;
37 j4p1(2, 2) = -15212. / 1875.;
38 j4p1(2, 0) = -15931. / 1875.;
39 j4p1(2, 1) = 5427. / 1250.;
40
41 j4p2(0, 0) = -10.2451;
42 j4p2(1, 1) = -50.3422;
43 j4p2(2, 2) = 4.91177;
44 j4p2(2, 0) = 8.0325;
45 j4p2(2, 1) = -0.3657;
46
47 j5b1(1, 0) = 2402. / 1587.;
48 j5b1(1, 1) = -14924. / 1587.;
49
50 j5b2(1, 0) = 1296371522. / 39457581. - 34624. / 1081. * gsl_sf_zeta(3);
51 j5b2(1, 1) = -177621017. / 7555707. + 800. / 23. * gsl_sf_zeta(3);
52
53 j5b1(1, 0) = 581. / 1875.;
54 j5b1(1, 1) = -15212. / 1875.;
55
56 j5b2(1, 0) = 684990354. / 19140625. - 6976. / 245. * gsl_sf_zeta(3);
57 j5b2(1, 1) = -272751559. / 10546875. + 128. / 5. * gsl_sf_zeta(3);
58}
59
60CPenguinBox::~CPenguinBox()
61{}
62
63gslpp::vector<double> CPenguinBox::WMatchP(orders order)
64{
65 double x = model.getMatching().x_t(model.getMuw());
66
67 double l = log(model.getMuw() * model.getMuw() / model.Mw_tree() / model.Mw_tree());
68
69 switch (order) {
70 case(LO):
71 p(0) = 4.;
72 p(1) = 4.;
73 p(2) = 0.;
74
75 return (p);
76 case(NLO):
77 p(0) = 4. * 0.5 * 2. / 3. * (11. + 6. * l);
78 p(1) = -4. * 0.5 * 4. / 3. * (11. + 6. * l);
79 p(2) = 8. * (2. + l);
80
81 return (p);
82 case(NNLO):
83 p(0) = 4. * (-(135677. - 124095.) / 3600. + 58. / 18. * M_PI * M_PI - 0.5 * (2. / 3.)*
84 (112. / 9. + 34. * x + (20. / 3. + 16. * x) * log(x) - (8. + 16. * x) *
85 sqrt(4. * x - 1.) * gsl_sf_clausen(2. * asin(1. / (2. * sqrt(x))))))
86 + (5. / 36. * 238. * l) + 58. / 6. * l*l;
87 p(1) = -4. * ((135677. + 124095.) / 3600. - 44. / 18. * M_PI * M_PI + 0.5 * (4. / 3.)*
88 (112. / 9. + 34. * x + (20. / 3. + 16. * x) * log(x) - (8. + 16. * x) *
89 sqrt(4. * x - 1.) * gsl_sf_clausen(2. * asin(1. / (2. * sqrt(x))))))
90 + (-5. / 36. * 260. * l) - 44. / 6. * l*l;
91 p(2) = 4. * model.getCF()*(33. + 4. * M_PI * M_PI + 34. * l + 12. * l * l);
92
93 return (p);
94 default:
95 std::stringstream out;
96 out << order;
97 throw std::runtime_error("CPenguinBox::WMatchB() order" + out.str() + " not implemented");
98 }
99}
100
101gslpp::matrix<double> CPenguinBox::RGevolP(int nf)
102{
103
104 switch (nf) {
105 case(5):
106 eta = model.Als(model.getMuw()) / model.Als(model.getMub());
107
108 up(0, 0) = pow(eta, 6. / 23.);
109 up(1, 1) = pow(eta, -12. / 23.);
110 up(2, 2) = pow(eta, 1. / 23.);
111 up(2, 0) = 12. / 5. * (pow(eta, 6. / 23.) - pow(eta, 1. / 23.));
112 up(2, 1) = 6. / 13. * (pow(eta, -12. / 23.) - pow(eta, 1. / 23.));
113
114 return (up);
115 case(4):
116 eta = model.Als(model.getMub()) / model.Als(model.getMuc());
117
118 up(0, 0) = pow(eta, 6. / 25.);
119 up(1, 1) = pow(eta, -12. / 25.);
120 up(2, 2) = pow(eta, -1. / 25.);
121 up(2, 0) = 12. / 7. * (pow(eta, 6. / 25.) - pow(eta, -1. / 25.));
122 up(2, 1) = 6. / 11. * (pow(eta, -12. / 25.) - pow(eta, -1. / 25.));
123
124 return (up);
125 default:
126 std::stringstream out;
127 out << nf;
128 throw std::runtime_error("CPenguinBox::RgevolP() " + out.str() + " wrong number of falvours");
129 }
130}
131
132double CPenguinBox::BmatchP(orders order)
133{
134
135 double mub = model.getMub();
136 double Mb = model.Mrun(model.getMub(), model.getQuarks(QCD::BOTTOM).getMass_scale(),
137 model.getQuarks(QCD::BOTTOM).getMass(), QCD::BOTTOM, FULLNNLO);
138 eta = model.Als(model.getMuw()) / mub;
139
140 switch (order) {
141 case(NNLO):
142 return ((-2. / 3.) * log(mub * mub / Mb / Mb)
143 *((284704. / 2645. + 694522. / 20631. * eta) * pow(eta, 1. / 23.) - (1033492. / 7935.
144 + 8264. / 529. * eta) * pow(eta, 6. / 23.) + (3058. / 1587.
145 + 18136. / 6877. * eta) * pow(eta, 12. / 23.) + eta * (pow(eta, 1. / 23.) * WMatchP(NLO)(2)
146 + 48. / 5. * (pow(eta, 6. / 23.) - pow(eta, 1. / 23.)) * WMatchP(NLO)(0)
147 + 24. / 13. * (pow(eta, -12. / 23.) - pow(eta, 1. / 23.)) * WMatchP(NLO)(1))));
148 case(NLO):
149 return (0.);
150 case(LO):
151 return (0.);
152 default:
153 std::stringstream out;
154 out << order;
155 throw std::runtime_error("CPenguinBox::BmatchB() order " + out.str() + " not implemented");
156 }
157}
158
159gslpp::vector<double> CPenguinBox::WMatchB(orders order)
160{
161 double l = log(model.getMuw() * model.getMuw() / model.Mw_tree() / model.Mw_tree());
162
163 switch (order) {
164 case(LO):
165 b(0) = 0.;
166 b(1) = 4.;
167
168 return (b);
169 case(NLO):
170 b(0) = -4. * (9. + 4. * l);
171 b(1) = 4.;
172
173 return (b);
174 case(NNLO):
175 b(0) = -8. * model.getCF()*(20. + 2. * M_PI * M_PI + 25. * l + 6. * l * l);
176 b(1) = 4.;
177
178 return (b);
179 default:
180 std::stringstream out;
181 out << order;
182 throw std::runtime_error("CPenguinBox::WMatchB() order" + out.str() + " not implemented");
183 }
184}
185
186gslpp::matrix<double> CPenguinBox::RGevolB(int nf)
187{
188
189 switch (nf) {
190 case(5):
191 eta = model.Als(model.getMuw()) / model.Als(model.getMub());
192
193 ub(0, 0) = 1.;
194 ub(1, 1) = pow(eta, 1. / 23.);
195 ub(1, 0) = 12. * (1. - pow(eta, 1. / 23.));
196
197 return (ub);
198 case(4):
199 eta = model.Als(model.getMub()) / model.Als(model.getMuc());
200
201 ub(0, 0) = 1.;
202 ub(1, 1) = pow(eta, -1. / 25.);
203 ub(1, 0) = -12. * (1. - pow(eta, -1. / 25.));
204
205 return (ub);
206 default:
207 std::stringstream out;
208 out << nf;
209 throw std::runtime_error("CPenguinBox::RgevolP() " + out.str() + " wrong number of flavours");
210 }
211}
212
213double CPenguinBox::BmatchB(orders order)
214{
215
216 double mub = model.getMub();
217 double Mb = model.Mrun(model.getMub(), model.getQuarks(QCD::BOTTOM).getMass_scale(),
218 model.getQuarks(QCD::BOTTOM).getMass(), QCD::BOTTOM, FULLNNLO);
219 eta = model.Als(model.getMuw()) / mub;
220
221 switch (order) {
222 case(NNLO):
223 return ( -2. / 3. * log(mub * mub / Mb / Mb)*(
224 (238784. / 529. - 9608. / 1587. * eta) * pow(eta, 1. / 23.) - 1336. / 3.
225 + pow(eta, 24. / 23.) * WMatchB(NLO)(1)));
226 case(NLO):
227 return (0.);
228 case(LO):
229 return (0.);
230 default:
231 std::stringstream out;
232 out << order;
233 throw std::runtime_error("CPenguinBox::BmatchB() order " + out.str() + " not implemented");
234 }
235}
236
237double CPenguinBox::Cmatch(orders order)
238{
239
240 up5 = RGevolP(5);
241 up4 = RGevolP(4);
242 ub5 = RGevolB(5);
243 ub4 = RGevolB(4);
244
245 p0 = WMatchP(LO);
246 b0 = WMatchB(LO);
247
248 double Muc = model.getMuc();
249 double Muw = model.getMuw();
250 double mc = model.Mrun(model.getMuc(), model.getQuarks(QCD::CHARM).getMass_scale(),
251 model.getQuarks(QCD::CHARM).getMass(), QCD::CHARM, FULLNNLO);
252 double l = log(Muc * Muc / mc / mc);
253 double A = model.Als(Muc, FULLNNLO);
254
255 //double etac = model.Als(Muc)/model.Als(mc);
256 double etab = model.Als(Muw) / model.Als(model.getMub());
257 double etacb = model.Als(model.getMub()) / model.Als(Muc);
258 double xtau = pow(1.7769 / Muw, 2.);
259 //double kc = pow(etac, 24./25.);
260
261 r1(0) = 4. * (1. - l);
262 r1(1) = -2. * (1. - l);
263
264 r2(0) = 11. + 20. * l - 12. * l*l;
265 r2(0) = -7. - 12. * l + 12. * l*l;
266
267 r1tau = 5. + 4. * xtau / (1 - xtau) * log(xtau) + 4. * l;
268 r2tau = -2. * model.getCF()*((9. + 7. * xtau) / (1. - xtau) + (3. * xtau + 13. * xtau * xtau) / (1. - xtau) / (1. - xtau)
269 - 12. * xtau * gsl_sf_dilog(1. - xtau) / (1. - xtau)
270 -((1. - 13. * xtau) / (1. - xtau) - 12. * xtau * xtau / (1. - xtau) / (1. - xtau) * log(xtau)) * l - 6. * l * l);
271
272 r1e = 5. + 4. * l;
273 r2e = -2. * model.getCF()*(9. - l - 6. * l * l);
274
275 /*
276 * this commented lines gives the matching in the approximation given in
277 * hep-ph/060379
278 *
279 *
280 rho1(0) = (4.*(1-2.*log(Muc/mc))
281 + 4.*log(kc));
282 rho1(1) = (-2.*(1-2.*log(Muc/mc))
283 - 2.*log(kc));
284
285 rho2(0) = 11.-40.*log(Muc/mc)
286 + 24.*Muc/mc
287 *Muc/mc
288 -20.*log(kc)-12.*log(kc)*log(kc)
289 +48.*log(kc)*log(Muc/mc)
290 +4.*15212./1875.*(1.-1/etac);
291
292 rho2(1) = -7.-24.*log(Muc/mc)
293 - 24.*Muc/mc
294 *Muc/mc
295 +12.*log(kc)+12.*log(kc)*log(kc)
296 -48.*log(kc)*log(Muc/mc)
297 -2.*15212./1875.*(1.-1/etac);
298
299 rho1tau = 5. + 4.*xtau/(1-xtau)*log(xtau)
300 + 16.*log(Muc/mc)
301 + 4./(xtau-kc)*(kc*log(kc)-xtau*(1.-kc)/(1.-xtau)*log(xtau));
302
303 rho2tau = -2.*model.getCF()*((9.+7.*xtau)/(1.-xtau) + (3.*xtau+13.*xtau*xtau)/(1.-xtau)/(1.-xtau)
304 -12.*xtau*gsl_sf_dilog(1.-xtau)/(1.-xtau)
305 -((1.-13.*xtau)/(1.-xtau) - 12.*xtau*xtau/(1.-xtau)/(1.-xtau)*log(xtau))*2.*log(Muc/mc)
306 -6.*log(Muc*Muc/mc/mc)*log(Muc*Muc/mc/mc))
307 + 32./(xtau-kc)*(4.*xtau*(1.-kc)/3.*(1.-xtau) -
308 - xtau*(xtau*13.-29.*xtau*xtau+kc*(3.+29*xtau*xtau)-kc*kc*(3.+13.*xtau))/
309 (12.*(xtau-kc)*(1.-xtau)*(1.-xtau))*log(xtau)
310 + kc*(17.*xtau-kc)/12./(xtau-kc)*log(kc)
311 + xtau*xtau/(xtau-kc)*log(xtau)+log(kc)
312 - (xtau*xtau+2.*xtau*kc-kc*kc)/(2.*(xtau-kc)*log(kc)*log(kc))
313 - (xtau*(xtau-kc))/2./(1.-xtau)*gsl_sf_dilog(1.-xtau)
314 - xtau*gsl_sf_dilog(1.-xtau/kc))
315 + 32./(xtau-kc)/(1.-xtau)*log(Muc*Muc/mc/mc)*(
316 xtau*(1.-kc)
317 + kc*(1.-xtau)*(2.*xtau-kc)/(xtau-kc)*log(xtau)
318 + xtau*xtau*(1.-kc)*(1.-2.*xtau+kc)/(xtau-kc)/(1.-xtau)*log(kc))
319 +4.*kc/(xtau-kc)*15212./1875.*(1.-1./etac)*(
320 1.-xtau/(xtau-kc)*log(xtau)+xtau/(xtau-kc)*log(kc));
321
322 rho1e = 5.+4.*log(Muc+Muc/mc/mc)-4.*log(kc);
323
324 rho2e = -2.*model.getCF()*(9.-log(Muc/Muc/mc/mc)-6.*log(Muc/Muc/mc/mc)*log(Muc/Muc/mc/mc))
325 -8./3.*log(kc)+16.*log(kc)*log(kc) - 32.*log(kc)*log(Muc/Muc/mc/mc)
326 -4.*15212./1875.*(1.-1./etac);
327 */
328
329 switch (order) {
330 case(NNLO):
331 p2 = WMatchP(NNLO);
332 p1 = WMatchP(NLO);
333
334 b2 = WMatchB(NNLO);
335 b1 = WMatchB(NLO);
336
337 CPM0 = up4 * up5 * p0;
338 CBM0 = ub4 * ub5 * b0;
339
340 CPM1 = j4p1 * up4 * up5 * p0
341 + etacb * up4 * (j5p1 - j4p1) * up5 * p0
342 + etab * etacb * up4 * up5 * (p1 - j5p1 * p0);
343 CBM1 = j4b1 * ub4 * ub5 * b0
344 + etacb * ub4 * (j5b1 - j4b1) * ub5 * b0
345 + etab * etacb * ub4 * ub5 * (b1 - j5b1 * b0);
346
347 CPM2 = j4p2 * up4 * up5 * p0
348 + etacb * j4p1 * up4 * (j5p1 - j4p1) * up5 * p0
349 + etab * etacb * j4p1 * up4 * up5 * (p1 - j5p1 * p0)
350 + etacb * etacb * up4 * (j5p2 - j4p2 - j4p1 * (j5p1 - j4p1) - BmatchP(NNLO)) * up5 * p0
351 + etab * etacb * etacb * up4 * (j5p1 - j4p1) * up5 * (p1 - j5p1 * p0)
352 + etab * etab * etacb * etacb * up4 * up5 * (p2 - j5p1 * p1 - (j5p2 - j5p1 * j5p1) * p0);
353 CBM2 = j4b2 * ub4 * ub5 * b0
354 + etacb * j4b1 * ub4 * (j5b1 - j4b1) * ub5 * b0
355 + etab * etacb * j4b1 * ub4 * ub5 * (b1 - j5b1 * b0)
356 + etacb * etacb * ub4 * (j5b2 - j4b2 - j4b1 * (j5b1 - j4b1) - BmatchB(NNLO)) * ub5 * b0
357 + etab * etacb * etacb * ub4 * (j5b1 - j4b1) * ub5 * (b1 - j5b1 * b0)
358 + etab * etab * etacb * etacb * ub4 * ub5 * (b2 - j5b1 * b1 - (j5b2 - j5b1 * j5b1) * b0);
359
360 cpm = pow(mc / Muw, 2.) / 32. * 4. * M_PI / A * (
361 CPM0(2) +
362
363 A / 4. / M_PI * (CPM1(2) + 4. * CPM0(0) * r1(0) + 4. * CPM0(1) * r1(1)) +
364
365 A / 4. / M_PI * A / 4. / M_PI * (CPM2(2) + 4. * CPM0(0) * r2(0) + 4. * CPM1(0) * r1(0)
366 + 4. * CPM0(1) * r2(1) + 4. * CPM1(1) * r1(1))
367 );
368
369 cbmt = pow(mc / Muw, 2.) / 16. * 4. * M_PI / A * (
370 CBM0(1) +
371
372 A / 4. / M_PI * (CBM1(1) + 4. * CBM0(0) * r1tau) +
373
374 A / 4. / M_PI * A / 4. / M_PI * (CBM2(1) + 4. * CBM0(0) * r2tau + 4. * CBM1(0) * r1tau)
375 );
376
377 cbme = pow(mc / Muw, 2.) / 16. * 4. * M_PI / A * (
378 CBM0(1) +
379
380 A / 4. / M_PI * (CBM1(1) + 4. * CBM0(0) * r1e) +
381
382 A / 4. / M_PI * A / 4. / M_PI * (CBM2(1) + 4. * CBM0(0) * r2e + 4. * CBM1(0) * r1e));
383
384 /*
385 cpm = kc * pow(mc/Muw, 2.) / 32. *(
386
387 4.*M_PI/model.Als(model.getMuc(), 4, FULLNNLO) *CPM0(2) +
388
389 CPM1(2) +
390 4.*CPM0(0)*r1(0) +
391 4.*CPM1(0)*r1(1) +
392 15212./1875.*(1.-1/etac)*CPM0(2) +
393
394 model.Als(Muc)/4./M_PI*(
395 CPM2(2) + 4.*CPM1(0)*r1(0) + 4*CPM0(0)*r2(0) + 4.*CPM1(1)*r1(1) + 4*CPM0(1)*r2(1) +
396 15212./1872.*(1-1./kc)*CPM1(2) +
397 4.*15212./1872.*(1.-1./kc)*(CPM0(0)*r1(0) + CPM0(1)*r1(1)) +
398 CPM0(2)*(966966391./10546875. -231404944./3515625./etac-272751559./10546875./etac/etac -
399 128./5.*(1.-1./etac/etac)*gsl_sf_zeta(3))));
400
401 cbmt = kc * pow(mc/Muw, 2.) / 16. * (
402
403 4.*M_PI/model.Als(Muc)*CBM0(1)
404
405 + CBM1(1) + 4.*r1tau*CBM1(0) + 15212./1875.*(1.-1./etac)*CBM0(1)
406
407 + model.Als(Muc)/4./M_PI * (CBM2(1) + 4.*r2tau*CBM0(0)
408 + 15212./1875.*(1.-1./etac)*CBM1(1)
409 + ( 966966391./10546875. - 231404944./3515625./etac - 272751559./10546875./etac/etac
410 -128./5.*(1.-1/etac/etac)*gsl_sf_zeta(3))*CBM0(1) ) );
411
412 cbme = kc * pow(model.Mrun(Muc, mc)/Muw, 2.) / 16. * (
413
414 4.*M_PI/model.Als(Muc)*CBM0(1)
415
416 + CBM1(1) + 4.*r1e*CBM1(0) + 15212./1875.*(1.-1./etac)*CBM0(1)
417
418 + model.Als(Muc)/4./M_PI * (CBM2(1) + 4.*r2e*CBM0(0)
419 + 15212./1875.*(1.-1./etac)*CBM1(1)
420 + ( 966966391./10546875. - 231404944./3515625./etac - 272751559./10546875./etac/etac
421 -128./5.*(1.-1/etac/etac)*gsl_sf_zeta(3))*CBM0(1) ) );
422 */
423 break;
424 case(NLO):
425 p1 = WMatchP(NLO);
426 b1 = WMatchB(NLO);
427
428 CPM0 = up4 * up5 * p0;
429 CBM0 = ub4 * ub5 * b0;
430
431 CPM1 = j4p1 * up4 * up5 * p0
432 + etacb * up4 * (j5p1 - j4p1) * up5 * p0
433 + etab * etacb * up4 * up5 * (p1 - j5p1 * p0);
434 CBM1 = j4b1 * ub4 * ub5 * b0
435 + etacb * ub4 * (j5b1 - j4b1) * ub5 * b0
436 + etab * etacb * ub4 * ub5 * (b1 - j5b1 * b0);
437
438 cpm = pow(mc / Muw, 2.) / 32. * 4. * M_PI / A * (
439
440 CPM0(2) +
441
442 A / 4. / M_PI * (CPM1(2) + 4. * CPM0(0) * r1(0) + 4. * CPM0(1) * r1(1)));
443
444 cbmt = pow(mc / Muw, 2.) / 16. * 4. * M_PI / A * (
445 CBM0(1) +
446
447 A / 4. / M_PI * (CBM1(1) + 4. * CBM0(0) * r1tau));
448
449 cbme = pow(mc / Muw, 2.) / 16. * 4. * M_PI / A * (
450 CBM0(1) +
451
452 A / 4. / M_PI * (CBM1(1) + 4. * CBM0(0) * r1e));
453 /*
454 cpm = kc * pow(mc/Muw, 2.) / 32. *(
455
456 4.*M_PI/model.Als(Muc) *CPM0(2) +
457
458 CPM1(2) + 4.*CPM0(0)*(4.*(1-2.*log(Muc/
459 mc/Muw)) + 4.*log(kc)) +
460 4.*CPM1(0)*(-2.*(1-2.*log(Muc/
461 mc/Muw)) - 2.*log(kc)) +
462 15212./1875.*(1.-1/etac)*CPM0(2));
463
464 cbmt = kc * pow(model.Mrun(Muc, mc)/Muw, 2.) / 16. * (
465
466 4.*M_PI/model.Als(Muc)*CBM0(1)
467
468 + CBM1(1) + 4.*r1tau*CBM1(0) + 15212./1875.*(1.-1./etac)*CBM0(1) );
469
470 cbme = kc * pow(model.Mrun(Muc, mc)/Muw, 2.) / 16. * (
471
472 4.*M_PI/model.Als(Muc)*CBM0(1)
473
474 + CBM1(1) + 4.*r1e*CBM1(0) + 15212./1875.*(1.-1./etac)*CBM0(1) );
475 */
476 break;
477 case(LO):
478 CPM0 = up4 * up5 * p0;
479
480 cpm = pow(mc / Muw, 2.) / 32. * 4. * M_PI / A * CPM0(2);
481
482 cbmt = pow(mc / Muw, 2.) / 16. * 4. * M_PI / A * CBM0(1);
483
484 cbme = pow(mc / Muw, 2.) / 16. * 4. * M_PI / A * CBM0(1);
485 break;
486 default:
487 std::stringstream out;
488 out << order;
489 throw std::runtime_error("CPenguinBox::Cmatch() order " + out.str() + " not implemented");
490 }
491
492 return (3. * cpm + 2. / 3. * cbme + 1. / 3. * cbmt);
493}
494
495double CPenguinBox::CT_tot(orders order, orders_qed order_qed)
496{
497
498 double IBT = model.getOptionalParameter("DeltaP_cu");
499 double xt = model.getMatching().x_t(model.getMuw());
500 double a = model.getMatching().mt2omh2(model.getMuw());
501
502 switch (order) {
503 case(NNLO):
504 switch (order_qed) {
505 case(NLO_QED11):
506 return (model.getCKM().computelamc().real() / model.getCKM().getLambda()*(Cmatch(NNLO) + IBT) +
507 model.getCKM().computelamt().real() / model.getCKM().getLambda()*(model.Als(model.getMuw()) / 4. / M_PI *
508 model.getMatching().X1t(xt) + model.getMatching().X1t(xt) + model.getAle() / 4. / M_PI *
509 model.getMatching().Xewt(xt, a, model.getMuc())));
510 case(LO_QED):
511 return (model.getCKM().computelamc().real() / model.getCKM().getLambda()*(Cmatch(NNLO) + IBT) +
512 model.getCKM().computelamt().real() / model.getCKM().getLambda()*(model.Als(model.getMuw()) / 4. / M_PI *
513 model.getMatching().X1t(xt) + model.getMatching().X1t(xt)));
514 default:
515 std::stringstream out;
516 out << order_qed;
517 throw std::runtime_error("CPenguinBox::CT_tot() order_qed " + out.str() + " not implemented");
518 }
519 case(NLO):
520 switch (order_qed) {
521 case(NLO_QED11):
522 return (model.getCKM().computelamc().real() / model.getCKM().getLambda()*(Cmatch(NLO) + IBT) +
523 model.getCKM().computelamt().real() / model.getCKM().getLambda()*(model.Als(model.getMuw()) / 4. / M_PI *
524 model.getMatching().X1t(xt) + model.getMatching().X1t(xt) + model.getAle() / 4. / M_PI *
525 model.getMatching().Xewt(xt, a, model.getMuc())));
526 case(LO_QED):
527 return (model.getCKM().computelamc().real() / model.getCKM().getLambda()*(Cmatch(NLO) + IBT) +
528 model.getCKM().computelamt().real() / model.getCKM().getLambda()*(model.Als(model.getMuw()) / 4. / M_PI *
529 model.getMatching().X1t(xt) + model.getMatching().X1t(xt)));
530 default:
531 std::stringstream out;
532 out << order_qed;
533 throw std::runtime_error("CPenguinBox::CT_tot() order_qed " + out.str() + " not implemented");
534 }
535 case(LO):
536 switch (order_qed) {
537 case(NLO_QED11):
538 return (model.getCKM().computelamc().real() / model.getCKM().getLambda()*(Cmatch(LO) + IBT) +
539 model.getCKM().computelamt().real() / model.getCKM().getLambda()*(model.Als(model.getMuw()) / 4. / M_PI *
540 model.getMatching().X1t(xt) + model.getAle() / 4. / M_PI *
541 model.getMatching().Xewt(xt, a, model.getMuc())));
542 case(LO_QED):
543 return (model.getCKM().computelamc().real() / model.getCKM().getLambda()*(Cmatch(NLO) + IBT) +
544 model.getCKM().computelamt().real() / model.getCKM().getLambda()*(model.Als(model.getMuw()) / 4. / M_PI *
545 model.getMatching().X1t(xt)));
546 default:
547 std::stringstream out;
548 out << order_qed;
549 throw std::runtime_error("CPenguinBox::CT_tot() order_qed " + out.str() + " not implemented");
550 }
551 default:
552 std::stringstream out;
553 out << order;
554 throw std::runtime_error("CPenguinBox::CT_tot() order " + out.str() + " not implemented");
555 }
556}
NNLO
@ NNLO
Definition: OrderScheme.h:36
LO
@ LO
Definition: OrderScheme.h:34
NLO
@ NLO
Definition: OrderScheme.h:35
FULLNNLO
@ FULLNNLO
Definition: OrderScheme.h:39
NLO_QED11
@ NLO_QED11
Definition: OrderScheme.h:59
LO_QED
@ LO_QED
Definition: OrderScheme.h:58
CKM::computelamc
const gslpp::complex computelamc() const
The product of the CKM elements .
Definition: CKM.cpp:147
CKM::computelamt
const gslpp::complex computelamt() const
The product of the CKM elements .
Definition: CKM.cpp:142
CKM::getLambda
const double getLambda() const
A member for returning the value of the Wolfenstein parameter .
Definition: CKM.h:88
CPenguinBox::CBM1
gslpp::vector< double > CBM1
Definition: CPenguinBox.h:108
CPenguinBox::p0
gslpp::vector< double > p0
Definition: CPenguinBox.h:102
CPenguinBox::model
const StandardModel & model
Definition: CPenguinBox.h:93
CPenguinBox::~CPenguinBox
~CPenguinBox()
destructor
Definition: CPenguinBox.cpp:60
CPenguinBox::r1e
double r1e
Definition: CPenguinBox.h:110
CPenguinBox::b
gslpp::vector< double > b
Definition: CPenguinBox.h:100
CPenguinBox::j4p2
gslpp::matrix< double > j4p2
Definition: CPenguinBox.h:95
CPenguinBox::p1
gslpp::vector< double > p1
Definition: CPenguinBox.h:102
CPenguinBox::r2tau
double r2tau
Definition: CPenguinBox.h:110
CPenguinBox::up5
gslpp::matrix< double > up5
Definition: CPenguinBox.h:98
CPenguinBox::RGevolP
gslpp::matrix< double > RGevolP(int nf)
Definition: CPenguinBox.cpp:101
CPenguinBox::CPenguinBox
CPenguinBox(const StandardModel &model_i)
constructor
Definition: CPenguinBox.cpp:12
CPenguinBox::RGevolB
gslpp::matrix< double > RGevolB(int nf)
Definition: CPenguinBox.cpp:186
CPenguinBox::BmatchB
double BmatchB(orders order)
Definition: CPenguinBox.cpp:213
CPenguinBox::r2e
double r2e
Definition: CPenguinBox.h:110
CPenguinBox::BmatchP
double BmatchP(orders order)
Definition: CPenguinBox.cpp:132
CPenguinBox::WMatchB
gslpp::vector< double > WMatchB(orders order)
Definition: CPenguinBox.cpp:159
CPenguinBox::CBM2
gslpp::vector< double > CBM2
Definition: CPenguinBox.h:108
CPenguinBox::cbme
double cbme
Definition: CPenguinBox.h:114
CPenguinBox::j4b2
gslpp::matrix< double > j4b2
Definition: CPenguinBox.h:96
CPenguinBox::j5b1
gslpp::matrix< double > j5b1
Definition: CPenguinBox.h:96
CPenguinBox::p
gslpp::vector< double > p
Definition: CPenguinBox.h:100
CPenguinBox::b0
gslpp::vector< double > b0
Definition: CPenguinBox.h:103
CPenguinBox::r2
gslpp::vector< double > r2
Definition: CPenguinBox.h:105
CPenguinBox::up4
gslpp::matrix< double > up4
Definition: CPenguinBox.h:98
CPenguinBox::up
gslpp::matrix< double > up
Definition: CPenguinBox.h:95
CPenguinBox::p2
gslpp::vector< double > p2
Definition: CPenguinBox.h:102
CPenguinBox::cpm
double cpm
Definition: CPenguinBox.h:114
CPenguinBox::r1
gslpp::vector< double > r1
Definition: CPenguinBox.h:105
CPenguinBox::CT_tot
double CT_tot(orders order, orders_qed order_qed)
NNLO QCD correction from hep-ph/0603079, two loop EW corrections from hep-ph/1009....
Definition: CPenguinBox.cpp:495
CPenguinBox::Cmatch
double Cmatch(orders order)
Definition: CPenguinBox.cpp:237
CPenguinBox::eta
double eta
Definition: CPenguinBox.h:112
CPenguinBox::j5p1
gslpp::matrix< double > j5p1
Definition: CPenguinBox.h:95
CPenguinBox::ub5
gslpp::matrix< double > ub5
Definition: CPenguinBox.h:98
CPenguinBox::j5b2
gslpp::matrix< double > j5b2
Definition: CPenguinBox.h:96
CPenguinBox::CPM0
gslpp::vector< double > CPM0
Definition: CPenguinBox.h:107
CPenguinBox::ub4
gslpp::matrix< double > ub4
Definition: CPenguinBox.h:98
CPenguinBox::j4p1
gslpp::matrix< double > j4p1
Definition: CPenguinBox.h:95
CPenguinBox::cbmt
double cbmt
Definition: CPenguinBox.h:114
CPenguinBox::j5p2
gslpp::matrix< double > j5p2
Definition: CPenguinBox.h:95
CPenguinBox::b2
gslpp::vector< double > b2
Definition: CPenguinBox.h:103
CPenguinBox::CBM0
gslpp::vector< double > CBM0
Definition: CPenguinBox.h:108
CPenguinBox::r1tau
double r1tau
Definition: CPenguinBox.h:110
CPenguinBox::b1
gslpp::vector< double > b1
Definition: CPenguinBox.h:103
CPenguinBox::CPM1
gslpp::vector< double > CPM1
Definition: CPenguinBox.h:107
CPenguinBox::CPM2
gslpp::vector< double > CPM2
Definition: CPenguinBox.h:107
CPenguinBox::ub
gslpp::matrix< double > ub
Definition: CPenguinBox.h:96
CPenguinBox::WMatchP
gslpp::vector< double > WMatchP(orders order)
Definition: CPenguinBox.cpp:63
CPenguinBox::j4b1
gslpp::matrix< double > j4b1
Definition: CPenguinBox.h:96
Particle::getMass_scale
double getMass_scale() const
A get method to access the scale at which the particle mass is defined.
Definition: Particle.h:133
Particle::getMass
const double & getMass() const
A get method to access the particle mass.
Definition: Particle.h:61
QCD::getMuc
const double getMuc() const
A get method to access the threshold between four- and three-flavour theory in GeV.
Definition: QCD.h:591
QCD::getCF
const double getCF() const
A get method to access the Casimir factor of QCD.
Definition: QCD.h:618
QCD::getOptionalParameter
const double getOptionalParameter(std::string name) const
A method to get parameters that are specific to only one set of observables.
Definition: QCD.h:450
QCD::Mrun
const double Mrun(const double mu, const double m, const quark q, const orders order=FULLNNLO) const
Computes a running quark mass from .
Definition: QCD.cpp:1353
QCD::BOTTOM
@ BOTTOM
Definition: QCD.h:329
QCD::CHARM
@ CHARM
Definition: QCD.h:326
QCD::getQuarks
const Particle & getQuarks(const QCD::quark q) const
A get method to access a quark as an object of the type Particle.
Definition: QCD.h:536
QCD::getMub
const double getMub() const
A get method to access the threshold between five- and four-flavour theory in GeV.
Definition: QCD.h:582
StandardModel
A model class for the Standard Model.
Definition: StandardModel.h:509
StandardModel::getMuw
const double getMuw() const
A get method to retrieve the matching scale around the weak scale.
Definition: StandardModel.h:951
StandardModel::getCKM
const CKM & getCKM() const
A get method to retrieve the member object of type CKM.
Definition: StandardModel.h:920
StandardModel::getMatching
virtual StandardModelMatching & getMatching() const
A get method to access the member reference of type StandardModelMatching.
Definition: StandardModel.h:965
StandardModel::Als
const double Als(const double mu, const orders order, const bool Nf_thr, const bool qed_flag) const
The running QCD coupling in the scheme including QED corrections.
Definition: StandardModel.h:1088
StandardModel::Mw_tree
const double Mw_tree() const
The tree-level mass of the boson, .
Definition: StandardModel/src/StandardModel.cpp:997
StandardModel::getAle
const double getAle() const
A get method to retrieve the fine-structure constant .
Definition: StandardModel.h:789
orders
orders
An enum type for orders in QCD.
Definition: OrderScheme.h:33
orders_qed
orders_qed
An enum type for orders in electroweak.
Definition: OrderScheme.h:56