a Code for the Combination of Indirect and Direct Constraints on High Energy Physics Models Logo
Charm_Kpnunu.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 "Charm_Kpnunu.h"
9 #include "StandardModel.h"
10 
11 Charm_Kpnunu::Charm_Kpnunu(const StandardModel& model_i)
12 : model(model_i),
13 cp(3, 0.), dcp(3, 0.), c_p(3, 0.), cpmuW0(3, 0.), cpmuW1(3, 0.), cpmuW2(3, 0.),
14 cb(2, 0.), dcb(2, 0.), c_b(2, 0.), cbmuW0(2, 0.), cbmuW1(2, 0.), cbmuW2(2, 0.),
15 U4p(3, 0.), U5p(3, 0.), J5p1(3, 0.), J4p1(3, 0.), J5p2(3, 0.), J4p2(3, 0.), dc_p(3, 0.),
16 U4b(2, 0.), U5b(2, 0.), J5b1(2, 0.), J4b1(2, 0.), J5b2(2, 0.), J4b2(2, 0.), dc_b(2, 0.)
17 {
18 
19  CP = 0.;
20  CBe = 0.;
21  CBt = 0.;
22  etab = model.Als(model.getMuw()) / model.Als(model.getMub());
23  etacb = model.Als(model.getMub()) / model.Als(model.getMuc());
24  // WARNING: we need to put mc(mc) for xc, line right below to be checked!
25  mc = model.Mrun(model.getMuc(), model.getQuarks(QCD::CHARM).getMass_scale(),
26  model.getQuarks(QCD::CHARM).getMass(), FULLNNLO);
27  etac = model.Als(model.getMuc()) / model.Als(mc);
28  kc = pow(etac, 24. / 25.);
29  xi1c = 15212. / 1875. * (etac - 1.) / etac;
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  xc = kc * (1. + model.Als(mc) / 4. / M_PI * xi1c + (model.Als(mc) / 4. / M_PI)*(model.Als(mc) / 4. / M_PI) * xi2c);
33 }
34 
35 Charm_Kpnunu::~Charm_Kpnunu()
36 {}
37 
38 gslpp::vector<double> Charm_Kpnunu::Cp(orders order)
39 {
40 
41  double x = model.getMatching().x_t(model.getMuw());
42  double L = log(model.getMuw() * model.getMuw() / model.Mw_tree() / model.Mw_tree());
43 
44  switch (order) {
45  case(LO):
46  cp(0) = 4.;
47  cp(1) = 4.;
48  cp(2) = 0.;
49 
50  return (cp);
51 
52  case(NLO):
53  cp(0) = 4. / 3. * (11. + 6. * L);
54  cp(1) = -8. / 3. * (11. + 6. * L);
55  cp(2) = 8. * (2. + L);
56 
57  return (cp);
58 
59  case(NNLO):
60  cp(0) = 4. * (-(135677. - 124095.) / 3600. + 58. / 18. * M_PI * M_PI - 0.5 * (2. / 3.)*
61  (112. / 9. + 32. * x + (20. / 3. + 16. * x) * log(x) - (8. + 16. * x) *
62  sqrt(4. * x - 1.) * gsl_sf_clausen(2. * asin(1. / (2. * sqrt(x)))))
63  +(5. / 36. * 238. * L) + 58. / 6. * L * L);
64  cp(1) = 4. * (-(135677. + 124095.) / 3600. - 44. / 18. * M_PI * M_PI + 0.5 * (4. / 3.)*
65  (112. / 9. + 32. * x + (20. / 3. + 16. * x) * log(x) - (8. + 16. * x) *
66  sqrt(4. * x - 1.) * gsl_sf_clausen(2. * asin(1. / (2. * sqrt(x)))))
67  - (5. / 36. * 260. * L) - 44. / 6. * L * L);
68  cp(2) = 4. * model.getCF()*(33. + 4. * M_PI * M_PI + 34. * L + 12. * L * L);
69 
70  return (cp);
71 
72  default:
73  std::stringstream out;
74  out << order;
75  throw std::runtime_error("Charm_Kpnunu::Cp() order" + out.str() + " not implemented");
76  }
77 }
78 
79 gslpp::matrix<double> Charm_Kpnunu::RGevolP(int nf, orders order)
80 {
81 
82  gslpp::matrix<double> evo(3, 3, 0.);
83 
84  switch (order) {
85  case (LO):
86  switch (nf) {
87  case(5):
88  evo(0, 0) = pow(etab, 6. / 23.);
89  evo(1, 1) = pow(etab, -12. / 23.);
90  evo(2, 0) = 12. / 5. * (pow(etab, 6. / 23.) - pow(etab, 1. / 23.));
91  evo(2, 1) = 6. / 13. * (pow(etab, -12. / 23.) - pow(etab, 1. / 23.));
92  evo(2, 2) = pow(etab, 1. / 23.);
93 
94  return (evo);
95 
96  case(4):
97  evo(0, 0) = pow(etacb, 6. / 25.);
98  evo(1, 1) = pow(etacb, -12. / 25.);
99  evo(2, 0) = 12. / 7. * (pow(etacb, 6. / 25.) - pow(etacb, -1. / 25.));
100  evo(2, 1) = 6. / 11. * (pow(etacb, -12. / 25.) - pow(etacb, -1. / 25.));
101  evo(2, 2) = pow(etacb, -1. / 25.);
102 
103  return (evo);
104 
105  default:
106  std::stringstream out;
107  out << nf;
108  throw std::runtime_error("Charm_Kpnunu::RgevolP() " + out.str() + " wrong number of falvours");
109  }
110  case(NLO):
111  switch (nf) {
112  case(5):
113  evo(0, 0) = 5165. / 3174.;
114  evo(1, 1) = -2267. / 1587.;
115  evo(2, 0) = -15857. / 1587.;
116  evo(2, 1) = 15305. / 3174.;
117  evo(2, 2) = -14924. / 1587.;
118 
119  return (evo);
120 
121  case(4):
122  evo(0, 0) = 6719. / 3750.;
123  evo(1, 1) = -3569. / 1875.;
124  evo(2, 0) = -15931. / 1875.;
125  evo(2, 1) = 5427. / 1250.;
126  evo(2, 2) = -15212. / 1875.;
127 
128  return (evo);
129 
130  default:
131  std::stringstream out;
132  out << nf;
133  throw std::runtime_error("Charm_Kpnunu::RgevolP() " + out.str() + " wrong number of flavours");
134  }
135  case(NNLO):
136  switch (nf) {
137  case(5):
138  evo(0, 0) = -7.35665;
139  evo(1, 1) = -54.9107;
140  evo(2, 0) = 17.7699;
141  evo(2, 1) = -1.7514;
142  evo(2, 2) = 18.3025;
143 
144  return (evo);
145 
146  case(4):
147  evo(0, 0) = -10.2451;
148  evo(1, 1) = -50.3422;
149  evo(2, 0) = 8.0325;
150  evo(2, 1) = -0.3657;
151  evo(2, 2) = 4.91177;
152 
153  return (evo);
154 
155  default:
156  std::stringstream out;
157  out << nf;
158  throw std::runtime_error("Charm_Kpnunu::RgevolP() " + out.str() + " wrong number of flavours");
159  }
160  default:
161  std::stringstream out;
162  out << order;
163  throw std::runtime_error("Charm_Kpnunu::RgevolP() " + out.str() + " wrong order assignment ");
164  }
165 }
166 
167 gslpp::vector<double> Charm_Kpnunu::ThresholdCp(orders order)
168 {
169 
170  double mub = model.getMub();
171  double Mb = model.Mrun(model.getMub(), model.getQuarks(QCD::BOTTOM).getMass_scale(),
172  model.getQuarks(QCD::BOTTOM).getMass(), FULLNNLO);
173  double L = log(mub * mub / Mb / Mb);
174 
175  switch (order) {
176  case(LO):
177  dcp = 0.;
178 
179  return (dcp);
180 
181  case(NLO):
182  dcp = 0.;
183 
184  return (dcp);
185 
186  case(NNLO):
187  dcp(0) = -pow(etab, 6. / 23.)*(2. / 3. * L * ((631. + 9699.) / 6348. * (1 - etab) + etab * Cp(NLO)(0))
188  - 2. / 3. * (59. / 36. + 1. / 3. * L + L * L));
189  dcp(1) = -pow(etab, 8. / 23.)*(2. / 3. * L * ((631. - 9699.) / 6348. * (1 - etab) + etab * Cp(NLO)(1))
190  + 4. / 3. * (59. / 36. + 1. / 3. * L + L * L));
191  dcp(2) = (-2. / 3.) * L * ((284704. / 2645. + 694522. / 20631. * etab) * pow(etab, 1. / 23.)
192  -(1033492. / 7935. + 8264. / 529. * etab) * pow(etab, 6. / 23.) + (3058. / 1587. + 18136. / 6877. * etab)
193  * pow(etab, -12. / 23.) + etab * (pow(etab, 1. / 23.) * Cp(NLO)(2) + 48. / 5. * (pow(etab, 6. / 23.)
194  - pow(etab, 1. / 23.)) * Cp(NLO)(0) + 24. / 13. * (pow(etab, -12. / 23.) - pow(etab, 1. / 23.)) * Cp(NLO)(1)));
195 
196  return (dcp);
197 
198  default:
199  std::stringstream out;
200  out << order;
201  throw std::runtime_error("Charm_Kpnunu::MatchingCp() " + out.str() + " wrong order assignment ");
202  }
203 }
204 
205 gslpp::vector<double> Charm_Kpnunu::C_p(orders order)
206 {
207 
208  cpmuW0 = Cp(LO);
209  cpmuW1 = Cp(NLO);
210  cpmuW2 = Cp(NNLO);
211 
212  U4p = RGevolP(4, LO);
213  U5p = RGevolP(5, LO);
214  J4p1 = RGevolP(4, NLO);
215  J5p1 = RGevolP(5, NLO);
216  J4p2 = RGevolP(4, NNLO);
217  J5p2 = RGevolP(5, NLO);
218 
219  for (int i = 0; i < 3; i++) {
220  dc_p(i, i) = ThresholdCp(NNLO)(i);
221  }
222 
223  switch (order) {
224  case(LO):
225  c_p = U4p * U5p*cpmuW0;
226 
227  return (c_p);
228 
229  case(NLO):
230  c_p = J4p1 * U4p * U5p * cpmuW0 + etacb * U4p * (J5p1 - J4p1) * U5p * cpmuW0 +
231  etab * etacb * U4p * U5p * (cpmuW1 - J5p1 * cpmuW0);
232 
233  return (c_p);
234 
235  case(NNLO):
236  //threshold term need to be checked (not really clear in (86) of hep-ph/0603079v2!)
237  c_p = J4p2 * U4p * U5p * cpmuW0 + etacb * J4p1 * U4p * (J5p1 - J4p1) * U5p * cpmuW0 +
238  etab * etacb * J4p1 * U4p * U5p * (cpmuW1 - J5p1 * cpmuW0) + etacb * etacb * U4p *
239  (J5p2 - J4p2 - J4p1 * (J5p1 - J4p1) - dc_p) * U5p * cpmuW0 + etab * etacb * etacb *
240  U4p * (J5p1 - J4p1) * U5p * (cpmuW1 - J5p1 * cpmuW0) + etab * etab * etacb * etacb * U4p * U5p *
241  (cpmuW2 - J5p1 * cpmuW1 - (J5p2 - J5p1 * J5p1) * cpmuW0);
242 
243  return (c_p);
244 
245  default:
246  std::stringstream out;
247  out << order;
248  throw std::runtime_error("Charm_Kpnunu::C_p() order" + out.str() + " not implemented");
249  }
250 
251 }
252 
253 double Charm_Kpnunu::C_P(orders order)
254 {
255 
256  double L = log(model.getMuc() * model.getMuc() / mc / mc);
257  double C_P0 = C_p(LO)(2);
258  double rhoP1p = 4 * (1. - L) + 4. * log(kc);
259  double rhoP1m = -2. * (1. - L) - 2 * log(kc);
260  double C_P1 = C_p(NLO)(2) + 4. * (C_p(LO)(0) * rhoP1p + C_p(LO)(1) * rhoP1m) + xi1c * C_p(LO)(2);
261  double rhoP2p = 11. + 20. * L - 12. * L * L - 20. * log(kc) - 12. * log(kc) * log(kc) +
262  24. * log(kc) * L + 4. * xi1c;
263  double rhoP2m = -7. + 13. * L + 12. * L * L - 12. * log(kc) + 12. * log(kc) * log(kc)
264  - 24. * log(kc) * L - 2. * xi1c;
265  double C_P2 = C_p(NNLO)(2) + 4. * (C_p(NLO)(0) * rhoP1p + C_p(NLO)(1) * rhoP1p + C_p(LO)(0) * rhoP2p
266  + C_p(LO)(1) * rhoP2m) + xi1c * (C_p(NLO)(2) + 4. * (C_p(LO)(0) * rhoP1p + C_p(LO)(1)
267  * rhoP1m)) + xi2c * C_p(LO)(2);
268 
269  switch (order) {
270  case(LO):
271  CP = kc * xc / 32. * (4. * M_PI / model.Als(model.getMuc()) * C_P0);
272 
273  return (CP);
274 
275  case(NLO):
276  CP = kc * xc / 32. * (4. * M_PI / model.Als(model.getMuc()) * C_P0 + C_P1);
277 
278  return (CP);
279 
280  case(NNLO):
281  CP = kc * xc / 32. * (4. * M_PI / model.Als(model.getMuc()) * C_P0 + C_P1 + model.Als(model.getMuc())
282  / 4. / M_PI * C_P2);
283 
284  return (CP);
285 
286  default:
287  std::stringstream out;
288  out << order;
289  throw std::runtime_error("Charm_Kpnunu::C_P() order" + out.str() + " not implemented");
290  }
291 }
292 
293 gslpp::vector<double> Charm_Kpnunu::Cb(orders order)
294 {
295 
296  double L = log(model.getMuw() * model.getMuw() / model.Mw_tree() / model.Mw_tree());
297 
298  switch (order) {
299  case(LO):
300  cb(0) = 4.;
301  cb(1) = 0.;
302 
303  return (cb);
304 
305  case(NLO):
306  cb(0) = 0.;
307  cb(1) = -4. * (9. + 4. * L);
308 
309  return (cb);
310 
311  case(NNLO):
312  cb(0) = 0.;
313  cb(1) = -8. * model.getCF()*(20. + 2. * M_PI * M_PI + 25. * L + 6. * L * L);
314 
315  return (cb);
316 
317  default:
318  std::stringstream out;
319  out << order;
320  throw std::runtime_error("Charm_Kpnunu::Cb() order" + out.str() + " not implemented");
321  }
322 }
323 
324 gslpp::matrix<double> Charm_Kpnunu::RGevolB(int nf, orders order)
325 {
326 
327  gslpp::matrix<double> evo(2, 2, 0.);
328 
329  switch (order) {
330  case (LO):
331  switch (nf) {
332  case(5):
333  evo(0, 0) = 1.;
334  evo(0, 1) = 0.;
335  evo(1, 0) = 12. * (1. - pow(etab, 1. / 23.));
336  evo(1, 1) = pow(etab, 1. / 23.);
337 
338  return (evo);
339 
340  case(4):
341  evo(0, 0) = 1.;
342  evo(0, 1) = 0.;
343  evo(1, 0) = -12. * (1. - pow(etacb, -1. / 25.));
344  evo(1, 1) = pow(etacb, -1. / 25.);
345 
346  return (evo);
347 
348  default:
349  std::stringstream out;
350  out << nf;
351  throw std::runtime_error("Charm_Kpnunu::RgevolB() " + out.str() + " wrong number of falvours");
352  }
353  case(NLO):
354  switch (nf) {
355  case(5):
356  evo(0, 0) = 0.;
357  evo(0, 1) = 0.;
358  evo(1, 0) = 2402. / 1587.;
359  evo(1, 1) = -14924. / 1587.;
360 
361  return (evo);
362 
363  case(4):
364  evo(0, 0) = 0.;
365  evo(0, 1) = 0.;
366  evo(1, 0) = 581. / 1875.;
367  evo(1, 1) = -15212. / 1875.;
368 
369  return (evo);
370 
371  default:
372  std::stringstream out;
373  out << nf;
374  throw std::runtime_error("Charm_Kpnunu::RgevolB() " + out.str() + " wrong number of flavours");
375  }
376  case(NNLO):
377  switch (nf) {
378  case(5):
379  evo(0, 0) = 0.;
380  evo(0, 1) = 0.;
381  evo(1, 0) = 1296371522. / 39457581. - 34624. / 1081. * gsl_sf_zeta_int(3);
382  evo(1, 1) = -177621017. / 7555707. + 800. / 23. * gsl_sf_zeta_int(3);
383 
384  return (evo);
385 
386  case(4):
387  evo(0, 0) = 0.;
388  evo(0, 1) = 0.;
389  evo(1, 0) = 684990354. / 19140625. - 6976. / 245. * gsl_sf_zeta_int(3);
390  evo(1, 1) = -272751559. / 10546875. + 128. / 5. * gsl_sf_zeta_int(3);
391 
392  return (evo);
393 
394  default:
395  std::stringstream out;
396  out << nf;
397  throw std::runtime_error("Charm_Kpnunu::RgevolB() " + out.str() + " wrong number of flavours");
398  }
399  default:
400  std::stringstream out;
401  out << order;
402  throw std::runtime_error("Charm_Kpnunu::RgevolB() " + out.str() + " wrong order assignment ");
403  }
404 }
405 
406 gslpp::vector<double> Charm_Kpnunu::ThresholdCb(orders order)
407 {
408 
409  double mub = model.getMub();
410  double Mb = model.Mrun(model.getMub(), model.getQuarks(QCD::BOTTOM).getMass_scale(),
411  model.getQuarks(QCD::BOTTOM).getMass(), FULLNNLO);
412  double L = log(mub * mub / Mb / Mb);
413 
414  switch (order) {
415  case(LO):
416  dcb = 0.;
417 
418  return (dcb);
419 
420  case(NLO):
421  dcb = 0.;
422 
423  return (dcb);
424 
425  case(NNLO):
426  dcb(0) = 0.;
427  dcb(1) = -2. / 3. * L * ((238784. / 529. - 9608. / 1587 * etab) * pow(etab, 1. / 23.) - 1336. / 3. +
428  pow(etab, 24. / 23.) * Cb(NLO)(1));
429 
430  return (dcb);
431 
432  default:
433  std::stringstream out;
434  out << order;
435  throw std::runtime_error("Charm_Kpnunu::MatchingCp() " + out.str() + " wrong order assignment ");
436  }
437 }
438 
439 gslpp::vector<double> Charm_Kpnunu::C_b(orders order)
440 {
441 
442  cbmuW0 = Cb(LO);
443  cbmuW1 = Cb(NLO);
444  cbmuW2 = Cb(NNLO);
445 
446  U4b = RGevolB(4, LO);
447  U5b = RGevolB(5, LO);
448  J4b1 = RGevolB(4, NLO);
449  J5b1 = RGevolB(5, NLO);
450  J4b2 = RGevolB(4, NNLO);
451  J5b2 = RGevolB(5, NLO);
452 
453  for (int i = 0; i < 3; i++) {
454  dc_b(i, i) = ThresholdCb(NNLO)(i);
455  }
456 
457  switch (order) {
458  case(LO):
459  c_b = U4b * U5b*cbmuW0;
460 
461  return (c_b);
462 
463  case(NLO):
464  c_b = J4b1 * U4b * U5b * cbmuW0 + etacb * U4b * (J5b1 - J4b1) * U5b * cbmuW0 +
465  etab * etacb * U4b * U5b * (cbmuW1 - J5b1 * cbmuW0);
466 
467  return (c_b);
468 
469  case(NNLO):
470  //threshold term need to be checked (not really clear in (86) of hep-ph/0603079v2!)
471  c_b = J4b2 * U4b * U5b * cbmuW0 + etacb * J4b1 * U4b * (J5b1 - J4b1) * U5b * cbmuW0 +
472  etab * etacb * J4b1 * U4b * U5b * (cbmuW1 - J5b1 * cbmuW0) + etacb * etacb * U4b *
473  (J5b2 - J4b2 - J4b1 * (J5b1 - J4b1) - dc_b) * U5b * cbmuW0 + etab * etacb * etacb *
474  U4b * (J5b1 - J4b1) * U5b * (cbmuW1 - J5b1 * cbmuW0) + etab * etab * etacb * etacb * U4b * U5b *
475  (cbmuW2 - J5b1 * cbmuW1 - (J5b2 - J5b1 * J5b1) * cbmuW0);
476 
477  return (c_b);
478 
479  default:
480  std::stringstream out;
481  out << order;
482  throw std::runtime_error("Charm_Kpnunu::C_p() order" + out.str() + " not implemented");
483  }
484 
485 }
486 
487 double Charm_Kpnunu::C_Be(orders order)
488 {
489 
490  double L = log(model.getMuc() * model.getMuc() / mc / mc);
491  double rhoB1_e = 5. + 4. * L - 4. * log(kc);
492  double rhoB2_e = -2. * model.getCF()*(9. - L - 6. * L * L) - 8. / 3. * log(kc) + 16. * log(kc) * log(kc)
493  - 32. * log(kc) * L - 4. * xi1c;
494  double C_B0e = C_b(LO)(1);
495  double C_B1e = C_b(NLO)(1) + 4. * rhoB1_e + xi1c * C_b(LO)(1);
496  double C_B2e = C_b(NNLO)(1) + 4. * rhoB2_e + xi1c * C_b(NLO)(1) + xi2c * C_b(LO)(1);
497 
498  switch (order) {
499  case(LO):
500  CBe = kc * xc / 32. * (4. * M_PI / model.Als(model.getMuc()) * C_B0e);
501 
502  return (CBe);
503 
504  case(NLO):
505  CBe = kc * xc / 32. * (4. * M_PI / model.Als(model.getMuc()) * C_B0e + C_B1e);
506 
507  return (CBe);
508 
509  case(NNLO):
510  CBe = kc * xc / 32. * (4. * M_PI / model.Als(model.getMuc()) * C_B0e + C_B1e + model.Als(model.getMuc()) / 4. / M_PI * C_B2e);
511 
512  return (CBe);
513 
514  default:
515  std::stringstream out;
516  out << order;
517  throw std::runtime_error("Charm_Kpnunu::C_Be() order" + out.str() + " not implemented");
518  }
519 }
520 
521 double Charm_Kpnunu::C_Bt(orders order)
522 {
523 
524  double L = log(model.getMuc() * model.getMuc() / mc / mc);
525  double x_t = model.getLeptons(model.TAU).getMass() * model.getLeptons(model.TAU).getMass() / mc / mc;
526 
527  double C_B0t = C_b(LO)(1);
528  double rhoB1_t = 5. + 4. * L - 4. * log(kc) + 4. / (x_t - kc)*(kc * log(kc)
529  - x_t * (1. - kc) / (1. - x_t) * log(x_t));
530  double C_B1t = C_b(NLO)(1) + 4. * rhoB1_t + xi1c * C_b(LO)(1);
531  double rhoB2_t = -2. * model.getCF()*((9. + 7. * x_t) / (1. - x_t) + x_t * (3. + 13. * x_t) /
532  (1. - x_t) / (1. - x_t) - 12. * x_t / (1. - x_t) * gsl_sf_dilog(1. - x_t)
533  - ((1. - 13. * x_t) / (1. - x_t) - 12. * x_t * x_t / (1. - x_t) / (1. - x_t) * log(x_t))
534  * L - 6. * L * L) + 32. / (x_t - kc)*(4. * x_t * (1 - kc) / 3. / (1. - x_t) - x_t
535  * ((x_t * (13. - 29. * x_t)) + kc * (3. + 29. * x_t * x_t) - kc * kc * (3. + 13. * x_t))
536  / 12. / (x_t - kc) / (1. - x_t) / (1. - x_t) * log(x_t) + kc * (17. * x_t - kc) / 12.
537  / (x_t - kc) * log(kc) + x_t * x_t / (x_t - kc) * log(x_t) * log(kc) - (x_t * x_t +
538  2. * x_t * kc - kc * kc) / 2. / (x_t - kc) * log(kc) * log(kc) - x_t * (x_t - kc) / 2.
539  / (1. - x_t) * gsl_sf_dilog(1. - x_t) - x_t * gsl_sf_dilog(1. - x_t / kc))
540  + 32. / (x_t - kc) / (1. - x_t)*(x_t * (1. - kc) + kc * (1. - x_t)*(2. * x_t - kc)
541  / (x_t - kc) * log(x_t) + x_t * x_t * (1. - kc)*(1. - 2. * x_t + kc) / (x_t - kc)
542  / (1. - x_t) * log(kc)) * L + 4. * kc / (x_t - kc)*(1. - x_t / (x_t - kc) * log(x_t)
543  + x_t / (x_t - kc) * log(kc)) * xi1c;
544  double C_B2t = C_b(NNLO)(1) + 4. * rhoB2_t + xi1c * C_b(NLO)(1) + xi2c * C_b(LO)(1);
545 
546  switch (order) {
547  case(LO):
548  CBt = kc * xc / 16. * (4. * M_PI / model.Als(model.getMuc()) * C_B0t);
549 
550  return (CBt);
551 
552  case(NLO):
553  CBt = kc * xc / 16. * (4. * M_PI / model.Als(model.getMuc()) * C_B0t + C_B1t);
554 
555  return (CBt);
556 
557  case(NNLO):
558  CBt = kc * xc / 16. * (4. * M_PI / model.Als(model.getMuc()) * C_B0t + C_B1t
559  + model.Als(model.getMuc()) / 4. / M_PI * C_B2t);
560 
561  return (CBt);
562 
563  default:
564  std::stringstream out;
565  out << order;
566  throw std::runtime_error("Charm_Kpnunu::C_Bt() order" + out.str() + " not implemented");
567  }
568 }
569 
570 double Charm_Kpnunu::P_C(orders order)
571 {
572  double Xe = C_P(order) + C_Be(order);
573  double Xt = C_P(order) + C_Bt(order);
574  double lambda4 = model.getCKM().getLambda() * model.getCKM().getLambda() * model.getCKM().getLambda() * model.getCKM().getLambda();
575  double pc = 1. / lambda4 * (2. / 3. * Xe + 1. / 3. * Xt);
576 
577  return (pc);
578 }
579 
580 double Charm_Kpnunu::C_TOT(orders order, orders_qed order_qed)
581 {
582 
583  double xt = model.getMatching().x_t(model.getMut());
584  double Muw = model.getMuw();
585  double Ale = model.getAle();
586  double a = 1. / model.getMatching().mt2omh2(Muw);
587  gslpp::complex lambdat = model.getCKM().computelamt();
588  gslpp::complex lambdac = model.getCKM().computelamc();
589  double lambda = model.getCKM().getLambda();
590  double lambda5 = model.getCKM().getLambda() * model.getCKM().getLambda() * model.getCKM().getLambda()
591  * model.getCKM().getLambda() * model.getCKM().getLambda();
592  double IBT = model.getOptionalParameter("DeltaP_cu");
593  double X = 0.;
594 
595  if ((order == NNLO) && (order_qed == NLO_QED11)) {
596  X = lambdat.real() / lambda5 * (model.getMatching().X0t(xt) + model.Als(model.getMut(), FULLNLO)
597  / 4. / M_PI * model.getMatching().X1t(xt) + Ale / 4. / M_PI * model.getMatching().Xewt(xt, a, Muw));
598  return (X + (lambdac.real() / lambda)*(P_C(NNLO) + IBT));
599  }
600 
601  if ((order == NLO) && (order_qed == LO_QED)) {
602  X = lambdat.real() / lambda5 * (model.getMatching().X0t(xt) / lambda5 + model.Als(model.getMut(), FULLNLO) / 4. / M_PI * model.getMatching().X1t(xt));
603  return (X + (lambdac.real() / lambda)*(P_C(NLO) + IBT));
604  }
605 
606  if ((order == LO) && (order_qed == LO_QED)) {
607  X = lambdat.real() * model.getMatching().X0t(xt) / lambda5;
608  return (X + (lambdac.real() / lambda)*(P_C(LO) + IBT));
609  }
610  else {
611  std::stringstream out;
612  out << order_qed;
613  throw std::runtime_error("BRKppnunu::BRKppnunu(): order_qed " + out.str() + "not implemented\n");
614  out << order;
615  throw std::runtime_error("BRKppnunu::BRKppnunu(): order " + out.str() + "not implemented");
616  }
617 }
QCD::TAU
Definition: QCD.h:316
Charm_Kpnunu::C_Be
double C_Be(orders order)
Definition: Charm_Kpnunu.cpp:487
QCD::BOTTOM
Definition: QCD.h:329
Charm_Kpnunu::xc
double xc
Definition: Charm_Kpnunu.h:150
Charm_Kpnunu::ThresholdCp
gslpp::vector< double > ThresholdCp(orders order)
Definition: Charm_Kpnunu.cpp:167
Charm_Kpnunu::dc_b
gslpp::matrix< double > dc_b
Definition: Charm_Kpnunu.h:148
gslpp::matrix< double >
A class for constructing and defining operations on real matrices.
Definition: gslpp_matrix_double.h:48
Charm_Kpnunu::U5p
gslpp::matrix< double > U5p
Definition: Charm_Kpnunu.h:148
NLO_QED11
Definition: OrderScheme.h:51
Charm_Kpnunu::cp
gslpp::vector< double > cp
Definition: Charm_Kpnunu.h:146
Charm_Kpnunu::kc
double kc
Definition: Charm_Kpnunu.h:150
Charm_Kpnunu::dc_p
gslpp::matrix< double > dc_p
Definition: Charm_Kpnunu.h:148
Charm_Kpnunu::P_C
double P_C(orders order)
Definition: Charm_Kpnunu.cpp:570
Charm_Kpnunu::J5b1
gslpp::matrix< double > J5b1
Definition: Charm_Kpnunu.h:148
Charm_Kpnunu::dcb
gslpp::vector< double > dcb
Definition: Charm_Kpnunu.h:146
Charm_Kpnunu::C_b
gslpp::vector< double > C_b(orders order)
Definition: Charm_Kpnunu.cpp:439
Charm_Kpnunu.h
LO
Definition: OrderScheme.h:33
Charm_Kpnunu::mc
double mc
Definition: Charm_Kpnunu.h:150
StandardModel.h
Charm_Kpnunu::etab
double etab
Definition: Charm_Kpnunu.h:150
QCD::CHARM
Definition: QCD.h:326
Charm_Kpnunu::RGevolP
gslpp::matrix< double > RGevolP(int nf, orders order)
Definition: Charm_Kpnunu.cpp:79
Charm_Kpnunu::CBt
double CBt
Definition: Charm_Kpnunu.h:150
gslpp::complex
A class for defining operations on and functions of complex numbers.
Definition: gslpp_complex.h:35
gslpp::log
complex log(const complex &z)
Definition: gslpp_complex.cpp:342
Charm_Kpnunu::cpmuW2
gslpp::vector< double > cpmuW2
Definition: Charm_Kpnunu.h:146
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
Charm_Kpnunu::ThresholdCb
gslpp::vector< double > ThresholdCb(orders order)
Definition: Charm_Kpnunu.cpp:406
Charm_Kpnunu::CP
double CP
Definition: Charm_Kpnunu.h:150
Charm_Kpnunu::c_b
gslpp::vector< double > c_b
Definition: Charm_Kpnunu.h:146
StandardModel
A model class for the Standard Model.
Definition: StandardModel.h:477
Charm_Kpnunu::xi1c
double xi1c
Definition: Charm_Kpnunu.h:150
Charm_Kpnunu::CBe
double CBe
Definition: Charm_Kpnunu.h:150
Charm_Kpnunu::C_Bt
double C_Bt(orders order)
Definition: Charm_Kpnunu.cpp:521
Charm_Kpnunu::J5p2
gslpp::matrix< double > J5p2
Definition: Charm_Kpnunu.h:148
Charm_Kpnunu::cbmuW0
gslpp::vector< double > cbmuW0
Definition: Charm_Kpnunu.h:146
Charm_Kpnunu::C_TOT
double C_TOT(orders order, orders_qed order_qed)
Definition: Charm_Kpnunu.cpp:580
QCD::getCF
double getCF() const
A get method to access the Casimir factor of QCD.
Definition: QCD.h:597
CKM::getLambda
double getLambda() const
A member for returning the value of the Wolfenstein parameter .
Definition: CKM.h:78
Charm_Kpnunu::U4b
gslpp::matrix< double > U4b
Definition: Charm_Kpnunu.h:148
Charm_Kpnunu::Charm_Kpnunu
Charm_Kpnunu(const StandardModel &model_i)
constructor
Definition: Charm_Kpnunu.cpp:11
Charm_Kpnunu::~Charm_Kpnunu
~Charm_Kpnunu()
destructor
Definition: Charm_Kpnunu.cpp:35
lambda4
An observable class for the quartic Higgs potential coupling .
Definition: THDMquantities.h:451
Particle::getMass
const double & getMass() const
A get method to access the particle mass.
Definition: Particle.h:61
Charm_Kpnunu::model
const StandardModel & model
Definition: Charm_Kpnunu.h:145
Charm_Kpnunu::U4p
gslpp::matrix< double > U4p
Definition: Charm_Kpnunu.h:148
gslpp::pow
complex pow(const complex &z1, const complex &z2)
Definition: gslpp_complex.cpp:395
Charm_Kpnunu::RGevolB
gslpp::matrix< double > RGevolB(int nf, orders order)
Definition: Charm_Kpnunu.cpp:324
gslpp::sqrt
complex sqrt(const complex &z)
Definition: gslpp_complex.cpp:385
Charm_Kpnunu::etac
double etac
Definition: Charm_Kpnunu.h:150
QCD::getMut
double getMut() const
A get method to access the threshold between six- and five-flavour theory in GeV.
Definition: QCD.h:561
Charm_Kpnunu::J4b2
gslpp::matrix< double > J4b2
Definition: Charm_Kpnunu.h:148
LO_QED
Definition: OrderScheme.h:50
QCD::getMuc
double getMuc() const
A get method to access the threshold between four- and three-flavour theory in GeV.
Definition: QCD.h:579
Charm_Kpnunu::C_p
gslpp::vector< double > C_p(orders order)
Definition: Charm_Kpnunu.cpp:205
NNLO
Definition: OrderScheme.h:35
StandardModel::Als
double Als(double mu, orders order=FULLNLO, bool qed_flag=false, bool Nf_thr=true) const
The running QCD coupling in the scheme including QED corrections.
Definition: StandardModel.cpp:602
QCD::getQuarks
Particle getQuarks(const QCD::quark q) const
A get method to access a quark as an object of the type Particle.
Definition: QCD.h:534
Charm_Kpnunu::cpmuW1
gslpp::vector< double > cpmuW1
Definition: Charm_Kpnunu.h:146
orders_qed
orders_qed
An enum type for orders in electroweak.
Definition: OrderScheme.h:47
Charm_Kpnunu::J4p2
gslpp::matrix< double > J4p2
Definition: Charm_Kpnunu.h:148
gslpp::vector< double >
A class for constructing and defining operations on real vectors.
Definition: gslpp_vector_double.h:33
Charm_Kpnunu::Cb
gslpp::vector< double > Cb(orders order)
Definition: Charm_Kpnunu.cpp:293
Charm_Kpnunu::cbmuW1
gslpp::vector< double > cbmuW1
Definition: Charm_Kpnunu.h:146
orders
orders
An enum type for orders in QCD.
Definition: OrderScheme.h:31
StandardModel::Mw_tree
virtual double Mw_tree() const
The tree-level mass of the boson, .
Definition: StandardModel.cpp:951
Charm_Kpnunu::c_p
gslpp::vector< double > c_p
Definition: Charm_Kpnunu.h:146
Charm_Kpnunu::U5b
gslpp::matrix< double > U5b
Definition: Charm_Kpnunu.h:148
QCD::getMub
double getMub() const
A get method to access the threshold between five- and four-flavour theory in GeV.
Definition: QCD.h:570
Charm_Kpnunu::etacb
double etacb
Definition: Charm_Kpnunu.h:150
StandardModel::getMatching
virtual StandardModelMatching & getMatching() const
A get method to access the member reference of type StandardModelMatching.
Definition: StandardModel.h:952
QCD::Mrun
double Mrun(const double mu, const double m, const orders order=FULLNNLO) const
Computes a running quark mass from .
Definition: QCD.cpp:1064
Charm_Kpnunu::J4b1
gslpp::matrix< double > J4b1
Definition: Charm_Kpnunu.h:148
Charm_Kpnunu::cpmuW0
gslpp::vector< double > cpmuW0
Definition: Charm_Kpnunu.h:146
lambda5
An observable class for the quartic Higgs potential coupling .
Definition: THDMquantities.h:474
gslpp::complex::real
const double & real() const
Definition: gslpp_complex.cpp:53
QCD::getOptionalParameter
double getOptionalParameter(std::string name) const
A method to get parameters that are specific to only one set of observables.
Definition: QCD.h:448
CKM::computelamt
gslpp::complex computelamt() const
The product of the CKM elements .
Definition: CKM.cpp:104
NLO
Definition: OrderScheme.h:34
Charm_Kpnunu::Cp
gslpp::vector< double > Cp(orders order)
Definition: Charm_Kpnunu.cpp:38
Charm_Kpnunu::J5b2
gslpp::matrix< double > J5b2
Definition: Charm_Kpnunu.h:148
Charm_Kpnunu::J5p1
gslpp::matrix< double > J5p1
Definition: Charm_Kpnunu.h:148
Charm_Kpnunu::dcp
gslpp::vector< double > dcp
Definition: Charm_Kpnunu.h:146
CKM::computelamc
gslpp::complex computelamc() const
The product of the CKM elements .
Definition: CKM.cpp:109
FULLNNLO
Definition: OrderScheme.h:38
Charm_Kpnunu::J4p1
gslpp::matrix< double > J4p1
Definition: Charm_Kpnunu.h:148
FULLNLO
Definition: OrderScheme.h:37
StandardModel::getAle
double getAle() const
A get method to retrieve the fine-structure constant .
Definition: StandardModel.h:748
StandardModel::getCKM
CKM getCKM() const
A get method to retrieve the member object of type CKM.
Definition: StandardModel.h:879
StandardModel::getMuw
double getMuw() const
A get method to retrieve the matching scale around the weak scale.
Definition: StandardModel.h:938
Charm_Kpnunu::cb
gslpp::vector< double > cb
Definition: Charm_Kpnunu.h:146
Charm_Kpnunu::xi2c
double xi2c
Definition: Charm_Kpnunu.h:150
Charm_Kpnunu::cbmuW2
gslpp::vector< double > cbmuW2
Definition: Charm_Kpnunu.h:146
Charm_Kpnunu::C_P
double C_P(orders order)
Definition: Charm_Kpnunu.cpp:253
StandardModel::getLeptons
Particle getLeptons(const QCD::lepton p) const
A get method to retrieve the member object of a lepton.
Definition: StandardModel.h:712