SoftSUSY is hosted by Hepforge, IPPP Durham
SOFTSUSY  4.0
lowe.h
Go to the documentation of this file.
1 
12 #ifndef LOWE_H
13 #define LOWE_H
14 
15 #include <iostream>
16 #include <fstream>
17 #include <sstream>
18 #include <cstring>
19 #include "def.h"
20 #include "utils.h"
21 #include "linalg.h"
22 #include "rge.h"
23 
25 namespace softsusy {
26 const double MUP = 2.2e-3;
27 const double MDOWN = 4.7e-3;
28 const double MSTRANGE = 0.096;
29 const double MCHARM = 1.27;
30 const double MBOTTOM = 4.18;
31 const double MTOP = 165.0;
32 const double MELECTRON = 5.109989461e-4;
34 const double MMUON = 1.0565835745e-1;
35 const double MTAU = 1.77686;
36 const double ALPHASMZ = 0.1181;
37 const double ALPHAMZ = 1.0 / 127.950;
38 
39 const double PMTOP = 173.21;
40 const double PMBOTTOM = 4.9;
41 const double THETA12CKM = 0.229206;
43 const double THETA13CKM = 0.003960;
44 const double THETA23CKM = 0.042223;
45 
47 typedef enum {mUp=1, mCharm, mTop, mDown, mStrange, mBottom, mElectron,
48  mMuon, mTau} mass;
50 typedef enum {ALPHA=1, ALPHAS} leGauge;
51 
53 DoubleVector gaugeDerivs(double, const DoubleVector &);
54 
56  class QedQcd: public RGE, public Approx {
57  private:
58  DoubleVector a;
59  DoubleVector mf;
60  double mtPole, mbPole;
61  double mbMb;
62  double mtauPole;
63 
64 public:
65  QedQcd();
66  QedQcd(const QedQcd &);
67  const QedQcd& operator=(const QedQcd & m);
68  virtual ~QedQcd() {};
69 
70  void setPoleMt(double mt) { mtPole = mt; };
71  void setPoleMb(double mb) { mbPole = mb; };
72  void setPoleMtau(double mtau) { mtauPole = mtau; };
73  void setMbMb(double mb) { mbMb = mb; };
74  void setMass(mass mno, double m) { mf(mno) = m; };
77  void setAlpha(leGauge ai, double ap) { a(ai) = ap; };
79  void set(const DoubleVector &);
80 
82  double displayPoleMt() const { return mtPole; };
84  double displayPoleMtau() const { return mtauPole; };
86  double displayPoleMb() const { return mbPole; };
88  const DoubleVector & displayMass() const { return mf; };
90  double displayMass(mass mno) const { return mf.display(mno); };
92  double displayAlpha(leGauge ai) const { return a.display(ai); };
94  const DoubleVector display() const;
96  double displayMbMb() const { return mbMb; }
97 
98  int flavours(double) const;
99 
100  double qedBeta() const;
101  double qcdBeta() const;
102  void massBeta(DoubleVector &) const;
103  DoubleVector beta() const;
105 
107  void runGauge(double start, double end);
109  double extractPoleMb(double asMb);
111  double extractRunningMb(double asMb);
113  void calcRunningMb();
116  void calcPoleMb();
117 
119  void toMt();
121  void toMz();
126  // alpha1 is in the GUT normalisation. sinth = sin^2 thetaW(Q) in MSbar
127  // scheme
128  DoubleVector getGaugeMu(const double m2, const
129  double sinth) const;
130 };
131 
133 ostream & operator <<(ostream &, const QedQcd &);
135 istream & operator >>(istream &left, QedQcd &m);
136 
140 void readIn(QedQcd & oneset, const char fname[80]);
143 double getRunMt(double poleMt, double asmt);
146 double getAsmt(double mtop, double alphasMz);
148 double getRunMtFromMz(double poleMt, double asMZ);
149 
150 inline QedQcd::QedQcd(const QedQcd &m)
151  : RGE(), Approx(m.displayApprox()), a(m.a), mf(m.mf), mtPole(m.mtPole),
152  mbPole(m.mbPole), mbMb(m.mbMb),
153  mtauPole(m.mtauPole) {
154  setPars(11);
155  setMu(m.displayMu());
156 }
157 
159 void massFermions(const QedQcd & r, DoubleMatrix & mDon,
160  DoubleMatrix & mUpq, DoubleMatrix & mEle);
161 
162 } // namespace softsusy
163 
164 #endif
165 
void setAlpha(leGauge ai, double ap)
sets QED or QCD structure constant
Definition: lowe.h:77
global variable declaration
Definition: def.cpp:13
void readIn(QedQcd &mset, const char fname[80])
Definition: lowe.cpp:443
DoubleVector is of variable length, and contains double precision.
Definition: linalg.h:35
leGauge
order of gauge couplings stored in QedQcd
Definition: lowe.h:50
void toMz()
Evolves object to MZ.
Definition: lowe.cpp:379
const double MTOP
Definition: lowe.h:31
switches (options) and parameters such as default fermion masses, required accuracy etc ...
DoubleVector and DoubleMatrix classes of doubles and operations between them, complexified copies als...
void calcRunningMb()
calculates running bottom mass given alpha_s(Mb)^{MSbar} from pole m_b
Definition: lowe.cpp:314
void setPoleMb(double mb)
set pole top mass
Definition: lowe.h:71
mass
used to give order of quark masses stored
Definition: lowe.h:47
void setMbMb(double mb)
set pole tau mass
Definition: lowe.h:73
double displayPoleMtau() const
Display pole tau mass.
Definition: lowe.h:84
void setMass(mass mno, double m)
sets a running quark mass
Definition: lowe.h:75
const double MBOTTOM
default running quark mass from PDG
Definition: lowe.h:30
DoubleVector beta() const
Beta functions of both beta-functions and all MSbar masses.
Definition: lowe.cpp:237
const double MUP
default running quark mass from PDG at 2 GeV
Definition: lowe.h:26
RGE objects consisting of energy scale and parameters and loops (order in perturbation theory) and th...
double displayMass(mass mno) const
Returns a single running mass.
Definition: lowe.h:90
void setPars(int i)
Set number of parameters in RGE object.
Definition: rge.h:65
double getAsmt(double mtop, double alphasMz)
Definition: lowe.cpp:498
Quark and lepton masses and gauge couplings in QEDxQCD effective theory.
Definition: lowe.h:56
double getRunMtFromMz(double poleMt, double asMZ)
Given pole mass and alphaS(MZ), returns running top mass – one loop qcd.
Definition: lowe.cpp:486
A few handy bits and pieces - little mathematical functions and the like.
void calcPoleMb()
Definition: lowe.cpp:340
double display(int i) const
display one element
Definition: linalg.h:90
double qedBeta() const
returns number of active flavours
Definition: lowe.cpp:151
const double ALPHAMZ
default running alpha(MZ) from PDG
Definition: lowe.h:37
const double PMTOP
PDG 2014.
Definition: lowe.h:39
istream & operator>>(istream &left, flavourPhysical &s)
Formatted input of physical parameters.
Definition: flavoursoft.cpp:1220
const double THETA12CKM
default central values of CKM matrix elements from PDG 2006 in radians
Definition: lowe.h:42
double displayPoleMb() const
Returns bottom "pole" mass.
Definition: lowe.h:86
QedQcd()
Initialises with default values defined in lowe.h.
Definition: lowe.cpp:17
double extractPoleMb(double asMb)
calculates pole bottom mass given alpha_s(Mb)^{MSbar} from running b mass
Definition: lowe.cpp:290
ostream & operator<<(ostream &left, const FlavourMssmSoftsusy &m)
Formatted output.
Definition: flavoursoft.cpp:217
double getRunMt(double poleMt, double asmt)
Definition: lowe.cpp:492
Describes a set of parameters and RGEs in general.
Definition: rge.h:50
void massBeta(DoubleVector &) const
Definition: lowe.cpp:204
double extractRunningMb(double asMb)
Done at pole mb: extracts running mb(polemb)
Definition: lowe.cpp:263
Matrix from 1..rows, 1..cols of double values.
Definition: linalg.h:211
const double MELECTRON
default pole lepton mass from PDG
Definition: lowe.h:33
void setPoleMtau(double mtau)
set pole bottom mass
Definition: lowe.h:72
void runGauge(double start, double end)
Does not run the masses, just gauge couplings from start to end.
Definition: lowe.cpp:248
double displayMu() const
Return renomalisation scale.
Definition: rge.h:68
DoubleVector gaugeDerivs(double x, const DoubleVector &y)
Returns beta functions of alpha, alpha_s only.
Definition: lowe.cpp:474
const double MCHARM
default running quark mass from PDG at 2 GeV
Definition: lowe.h:29
const double MMUON
default pole lepton mass from PDG
Definition: lowe.h:34
double qcdBeta() const
QCD beta function.
Definition: lowe.cpp:175
double displayPoleMt() const
Display pole top mass.
Definition: lowe.h:82
const double PMBOTTOM
Definition: lowe.h:40
Definition: rge.h:24
const double MSTRANGE
default running quark mass from PDG at 2 GeV
Definition: lowe.h:28
void massFermions(const QedQcd &r, DoubleMatrix &mDon, DoubleMatrix &mUpq, DoubleMatrix &mEle)
Returns diagonal fermion mass matrices given input object r.
Definition: lowe.cpp:504
void setMu(double e)
Sets renormalisation scale to e.
Definition: rge.h:60
double displayAlpha(leGauge ai) const
Returns a single gauge structure constant.
Definition: lowe.h:92
const DoubleVector & displayMass() const
Returns a vector of running fermion masses.
Definition: lowe.h:88
const double THETA23CKM
From Vcb/Vtb in global CKM fit, PDG.
Definition: lowe.h:44
const double MTAU
default pole lepton mass from PDG
Definition: lowe.h:35
DoubleVector getGaugeMu(const double m2, const double sinth) const
Definition: lowe.cpp:403
const double THETA13CKM
From Vub in global CKM fit, PDG.
Definition: lowe.h:43
void toMt()
Evolves object to running top mass.
Definition: lowe.cpp:357
double displayMbMb() const
Returns mb(mb) MSbar.
Definition: lowe.h:96
const QedQcd & operator=(const QedQcd &m)
Sets two objects equal.
Definition: lowe.cpp:32
const double MDOWN
default running quark mass from PDG at 2 GeV
Definition: lowe.h:27
const DoubleVector display() const
Obgligatory: returns vector of all running parameters.
Definition: lowe.cpp:54
const double ALPHASMZ
default running value from PDG
Definition: lowe.h:36