SoftSUSY is hosted by Hepforge, IPPP Durham
SOFTSUSY  4.1
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 "utils.h"
20 #include "linalg.h"
21 #include "rge.h"
22 
24 namespace softsusy {
25 const double MUP = 2.2e-3;
26 const double MDOWN = 4.7e-3;
27 const double MSTRANGE = 0.096;
28 const double MCHARM = 1.27;
29 const double MBOTTOM = 4.18;
30 const double MTOP = 165.0;
31 const double MELECTRON = 5.109989461e-4;
33 const double MMUON = 1.0565835745e-1;
34 const double MTAU = 1.77686;
35 const double ALPHASMZ = 0.1181;
36 const double ALPHAMZ = 1.0 / 127.950;
37 
38 const double PMTOP = 173.21;
39 const double PMBOTTOM = 4.9;
40 const double THETA12CKM = 0.229206;
42 const double THETA13CKM = 0.003960;
43 const double THETA23CKM = 0.042223;
44 
46 typedef enum {mUp=1, mCharm, mTop, mDown, mStrange, mBottom, mElectron,
47  mMuon, mTau} mass;
49 typedef enum {ALPHA=1, ALPHAS} leGauge;
50 
52 DoubleVector gaugeDerivs(double, const DoubleVector &);
53 
55  class QedQcd: public RGE, public Approx {
56  private:
57  DoubleVector a;
58  DoubleVector mf;
59  double mtPole, mbPole;
60  double mbMb;
61  double mtauPole;
62 
63 public:
64  QedQcd();
65  QedQcd(const QedQcd &);
66  const QedQcd& operator=(const QedQcd & m);
67  virtual ~QedQcd() {};
68 
69  void setPoleMt(double mt) { mtPole = mt; };
70  void setPoleMb(double mb) { mbPole = mb; };
71  void setPoleMtau(double mtau) { mtauPole = mtau; };
72  void setMbMb(double mb) { mbMb = mb; };
73  void setMass(mass mno, double m) { mf(mno) = m; };
76  void setAlpha(leGauge ai, double ap) { a(ai) = ap; };
78  void set(const DoubleVector &);
79 
81  double displayPoleMt() const { return mtPole; };
83  double displayPoleMtau() const { return mtauPole; };
85  double displayPoleMb() const { return mbPole; };
87  const DoubleVector & displayMass() const { return mf; };
89  double displayMass(mass mno) const { return mf.display(mno); };
91  double displayAlpha(leGauge ai) const { return a.display(ai); };
93  const DoubleVector display() const;
95  double displayMbMb() const { return mbMb; }
96 
97  int flavours(double) const;
98 
99  double qedBeta() const;
100  double qcdBeta() const;
101  void massBeta(DoubleVector &) const;
102  DoubleVector beta() const;
104 
106  void runGauge(double start, double end);
108  double extractPoleMb(double asMb);
110  double extractRunningMb(double asMb);
112  void calcRunningMb();
115  void calcPoleMb();
116 
118  void toMt();
120  void toMz();
125  // alpha1 is in the GUT normalisation. sinth = sin^2 thetaW(Q) in MSbar
126  // scheme
127  DoubleVector getGaugeMu(const double m2, const
128  double sinth) const;
129 };
130 
132 ostream & operator <<(ostream &, const QedQcd &);
134 istream & operator >>(istream &left, QedQcd &m);
135 
138 double getRunMt(double poleMt, double asmt);
141 double getAsmt(double mtop, double alphasMz);
143 double getRunMtFromMz(double poleMt, double asMZ);
144 
145 inline QedQcd::QedQcd(const QedQcd &m)
146  : RGE(), Approx(m.displayApprox()), a(m.a), mf(m.mf), mtPole(m.mtPole),
147  mbPole(m.mbPole), mbMb(m.mbMb),
148  mtauPole(m.mtauPole) {
149  setPars(11);
150  setMu(m.displayMu());
151 }
152 
154 void massFermions(const QedQcd & r, DoubleMatrix & mDon,
155  DoubleMatrix & mUpq, DoubleMatrix & mEle);
156 
157 } // namespace softsusy
158 
159 #endif
160 
void setAlpha(leGauge ai, double ap)
sets QED or QCD structure constant
Definition: lowe.h:76
global variable declaration
Definition: def.cpp:13
DoubleVector is of variable length, and contains double precision.
Definition: linalg.h:35
leGauge
order of gauge couplings stored in QedQcd
Definition: lowe.h:49
void toMz()
Evolves object to MZ.
Definition: lowe.cpp:379
const double MTOP
Definition: lowe.h:30
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:70
mass
used to give order of quark masses stored
Definition: lowe.h:46
void setMbMb(double mb)
set pole tau mass
Definition: lowe.h:72
double displayPoleMtau() const
Display pole tau mass.
Definition: lowe.h:83
void setMass(mass mno, double m)
sets a running quark mass
Definition: lowe.h:74
const double MBOTTOM
default running quark mass from PDG
Definition: lowe.h:29
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:25
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:89
void setPars(int i)
Set number of parameters in RGE object.
Definition: rge.h:64
double getAsmt(double mtop, double alphasMz)
Definition: lowe.cpp:460
Quark and lepton masses and gauge couplings in QEDxQCD effective theory.
Definition: lowe.h:55
double getRunMtFromMz(double poleMt, double asMZ)
Given pole mass and alphaS(MZ), returns running top mass – one loop qcd.
Definition: lowe.cpp:448
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:36
const double PMTOP
PDG 2014.
Definition: lowe.h:38
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:41
double displayPoleMb() const
Returns bottom "pole" mass.
Definition: lowe.h:85
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:454
Describes a set of parameters and RGEs in general.
Definition: rge.h:49
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:32
void setPoleMtau(double mtau)
set pole bottom mass
Definition: lowe.h:71
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:67
DoubleVector gaugeDerivs(double x, const DoubleVector &y)
Returns beta functions of alpha, alpha_s only.
Definition: lowe.cpp:436
const double MCHARM
default running quark mass from PDG at 2 GeV
Definition: lowe.h:28
const double MMUON
default pole lepton mass from PDG
Definition: lowe.h:33
double qcdBeta() const
QCD beta function.
Definition: lowe.cpp:175
double displayPoleMt() const
Display pole top mass.
Definition: lowe.h:81
const double PMBOTTOM
Definition: lowe.h:39
Definition: rge.h:23
const double MSTRANGE
default running quark mass from PDG at 2 GeV
Definition: lowe.h:27
void massFermions(const QedQcd &r, DoubleMatrix &mDon, DoubleMatrix &mUpq, DoubleMatrix &mEle)
Returns diagonal fermion mass matrices given input object r.
Definition: lowe.cpp:466
void setMu(double e)
Sets renormalisation scale to e.
Definition: rge.h:59
double displayAlpha(leGauge ai) const
Returns a single gauge structure constant.
Definition: lowe.h:91
const DoubleVector & displayMass() const
Returns a vector of running fermion masses.
Definition: lowe.h:87
const double THETA23CKM
From Vcb/Vtb in global CKM fit, PDG.
Definition: lowe.h:43
const double MTAU
default pole lepton mass from PDG
Definition: lowe.h:34
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:42
void toMt()
Evolves object to running top mass.
Definition: lowe.cpp:357
double displayMbMb() const
Returns mb(mb) MSbar.
Definition: lowe.h:95
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:26
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:35