qcdevol
Module¶
Introduction¶
This package provides tools for evolving QCD couplings, masses and other renormalization parameters. Typical usage is:
>>> import qcdevol as qcd
>>> al = qcd.Alpha_s(alpha0=0.2128, mu0=5, nf=4) # coupling
>>> mc = qcd.M_msb(m0=0.9851, mu0=3, alpha=al) # c-quark mass
>>> al(91.2), mc(91.2) # run to mu=91.2 GeV
(0.11269900754817504, 0.6333892033863693)
Uncertainties can be introduced into parameters using the gvar
module, and these are propagated through the evolution. QED effects
can be included. There are also tools for manipulating and evolving
perturbation series (in alpha_s).
QCD Coupling¶
The central component of the qcdevol
package is the class that implements
the QCD coupling constant:
- class qcdevol.Alpha_s(alpha0, mu0, nf, scheme='msb', alpha_qed=None, beta_msb=None, maxalpha=1.0, rtol=None, order=None)¶
Running coupling in QCD.
Typical usage is:
>>> import qcdevol as qcd >>> al = qcd.Alpha_s(alpha0=0.2128, mu0=5, nf=4) >>> al(5) # coupling at 5 GeV 0.2128 >>> al(91.2) # coupling at 91.2 GeV 0.11269900754817504 >>> al([5., 91.2]) array([0.2128 , 0.11269901])
- Parameters:
alpha0 – Value of coupling at scale
mu0
. This can be a number, or a variable of typegvar.GVar
if the value has an uncertainty (see the documentation forgvar
).mu0 – Scale at which
self(mu0) == alpha0
.nf (int) – Number of quark flavors included in vacuum polarization.
scheme (str) – Definition/scheme used for the coupling. By default this is either
'msb'
for the MSbar scheme or'v'
for the scheme defined in terms of the static quark potential; but new schemes can be added toqcdevol.SCHEMES
. Default is'msb'
.alpha_qed – QED coupling constant. If specified, the leading-order QED correction is added to the beta function. The running of the QED coupling is ignored in leading order;
alpha_qed
should be chosen to be appropriate for the range of scales of interest. Ignored ifalpha_qed=None
(default) or if parameterbeta_msb
is specified.beta_msb (array) – Coefficients of the MSbar beta function. If specified these are used in place of what is returned by
qcdevol.BETA_MSB()
. Ignored ifbeta_msb=None
(default).maxalpha (float) – Maximum value for alpha allowed in intermediate calculations, which are halted if it is exceeded. Default value is
maxalpha=1.
. Setmaxalpha=None
to have the maximum coupling calculated from the beta function.rtol (float) – Temporarily replaces the default
rtol
parameter; seeqcdevol.setdefaults()
.order (int) – Temporarily replaces the default
order
parameter; seeqcdevol.setdefaults()
.
qcdevol.Alpha_s
objects have the following methods:- __call__(mu, scheme=None)¶
Return the coupling at scale
mu
, which can be an array.If
scheme
is notNone
, the coupling is evaluated in schemescheme
, otherwise the scheme specified when creating the coupling is used. (Default isNone
)
- del_quark(m, mu=None, zeta=None)¶
Create a new coupling by removing a heavy quark from the vac. polarization.
Creates a new coupling with
self.nf-1
flavors by removing a quark whoseself.nf
-flavor mass ism
at scalemu
. Ifmu
is omitted (orNone
), it is set equal tom
.The default
zeta
parameters (fromqcdevol.ZETA2_G_MSB(self.nf, m, mu)
) can be replaced by specifying arrayzeta
. Ignored ifzeta=None
(default).
- add_quark(m, mu=None, zeta=None)¶
Create a new coupling by adding a heavy quark to the vac. polarization.
Creates a new coupling with
self.nf+1
flavors by adding a quark whose(self.nf+1)
-flavor mass ism
at scalemu
. Ifmu
is omitted (orNone
), it is set equal tom
.The default
zeta
parameters (fromqcdevol.ZETA2_G_MSB(self.nf+1, m, mu)
) can be replaced by specifying arrayzeta
. Ignored ifzeta=None
(default).
- clone(scheme=None, mu0=None)¶
Create clone of
self
.This method creates couplings that are equivalent to
self
but potentially in a different scheme or with a different normalization scalemu0
.- Parameters:
scheme (str or None) – Scheme used to define new coupling. Ignored if
scheme=None
(default).mu0 – New normalization scale for coupling. Ignored if
mu0=None
(default).
- Returns:
New object of type
Alpha_s
equivalent toself
but possibly in a different scheme (or with a differentmu0
).
- ps(mu, mu0, order=None)¶
Power series for
self(mu)
in terms ofself(mu0)
.Returns a power series
ps
of orderorder
and typegvar.powerseries.PowerSeries
such that such thatps(self(mu0))
is equal toself(mu)
through the order specified.- Parameters:
mu – The expansion gives
self(mu)
. Ifmu
is a tuple,(mu, scheme)
, the expansion is forself(*mu)
.mu0 – The expansion is in powers of
self(mu0)
. Ifmu0
is a tuple,(mu0,scheme0)
, the expansion is in powers ofself(*mu0)
.order – The order to which the expansion is carried out. Leaving
order=None
(default) setsorder
equal to the default set byqcdevol.setdefaults()
.
- exact(mu, scheme=None)¶
Same as
Alpha_s.__call__
except it integrates differential evolution equation directly. It is much slower, and is only included for testing purposes.
Operator Z Factors¶
Multiplying a local operator by its Z factor removes its scale dependence. The following class represents such Z factors.
- class qcdevol.OPZ(z0, mu0, alpha, gamma=None, gammat=None, order=None, rtol=None)¶
Z factor for operator with anomalous dimenion gamma:
mu**2 dlog(Z(mu))/dmu**2 = - gammat - alpha(mu) * (gamma[0] + gamma[1] * alpha(mu) + ....)
Typical usage is:
>>> import qcdevol as qcd >>> al = qcd.Alpha_s(alpha0=0.2128, mu0=5., nf=4) >>> Z = qcd.OPZ(z0=4., mu0=4., alpha=al, gamma=[0.318, 0.370, 0.321]) >>> Z(4.) # Z at mu = 4 GeV 4.0 >>> Z(1000.) # Z at mu = 1000 GeV 2.2907519504498515 >>> z([1., 4., 1000.]) array([5.74602581, 4. , 2.38456863])
- Parameters:
z0 – Value of Z factor at scale
mu0
. This can be a number, or a variable of typegvar.GVar
if the value has an uncertainty (see the documentation for thegvar
module).z0
can also be an array, representing multiple Z factors.mu0 – Scale at which
self(mu0)=z0
.alpha – QCD coupling (type
qcdevol.Alpha_s
) governing the Z factor’s evolution.gamma – An array specifying the Z factor’s anomalous dimension.
gammat – Tree-level gamma. Set to zero if
gammat=None
(default).order (int) – Temporarily replaces the default
order
parameter; seeqcdevol.setdefaults()
.rtol (float) – Temporarily replaces the default
rtol
parameter; seeqcdevol.setdefaults()
.
qcdevol.OPZ
objects have the following methods:- __call__(mu)¶
Return the Z factor at scale
mu
, which may be an array.If
self.z0
andmu
are both arrays, they must have the same shape as each element ofmu
is matched with the corresponding element ofself.z0
.
- clone(mu0=None)¶
Create clone of
self
.If
mu0
is notNone
, parametermu0
is reset in the clone (with the corresponding value forz0
).
- exact(mu)¶
Same as
OPZ.__call__
except it integrates the differential evolution equation directly. It is much slower, and is only included for testing purposes. (Particularly inefficient ifself.z0
and`mu
are both arrays.)
Running Quark Masses¶
MSbar quark masses are represented by objects of type qcdevol.M_msb
,
which are derived from qcdevol.OPZ
.
- class qcdevol.M_msb(m0, mu0, alpha, Q=None, gamma_msb=None, order=None, rtol=None)¶
MSbar quark masses.
Typical usage:
>>> import qcdevol as qcd >>> al = qcd.Alpha_s(alpha0=0.2128, mu0=5., nf=4) >>> mc = qcd.M_msb(m0=0.9851, mu0=3., alpha=al) >>> mc(3) # mc at 3 GeV 0.9851 >>> mc(1000.) # mc at 1000 GeV 0.5377724915576455 >>> mc([1., 3., 1000.]) array([1.43226913, 0.9851 , 0.53777249])
- Parameters:
m0 – Value of mass at scale mu0. This can be a number, or a variable of type
gvar.GVar
if the value has an uncertainty (see the documentation for thegvar
module).mu0 – Scale at which quark mass equals
m0
.alpha – QCD coupling (type
qcdevol.Alpha_s
) governing the mass’s evolution. Must be in the'msb'
scheme.gamma_msb (array) – Coefficients of the mass’s MSbar anomalous dimension. If specified, these are used in place of the values used by
qcdevol.GAMMA_MSB()
. Ignored ifgamma_msb=None
(default).Q – QED charge of quark in units of the proton charge. Ignored if
None
(default). When specified,alpha
must also include QED corrections.order (int) – Temporarily replaces the default
order
parameter; seeqcdevol.setdefaults()
.rtol (float) – Temporarily replaces the default
rtol
parameter; seeqcdevol.setdefaults()
.
qcdevol.M_msb
objects have the following methods:- __call__(mu)¶
Return mass’s value at scale
mu
, which may be an array.mu
is typically a number but can also be a string containing an arithmetic expression involving parameterm
. For example,>>> mb = qcd.M_msb(m0=4.513, mu0=3., alpha=al) >>> mb('m') # mu = mb(mu) 4.231289159789313 >>> mb('3*m') # mu = 3*mb(mu) 3.674028791373564
calculates
mb(mu)
first formu=mb(mu)
and then formu=mb(3*mu)
.
- del_quark(m=None, mu=None, zeta=None)¶
Create a new mass by removing a heavy quark from the vac. polarization.
Creates a new mass with
nf-1
flavors by removing a quark whosenf
-flavor mass ism
at scalemu
, wherenf
equalsself.alpha.nf
. If mu is omitted (or None), it is set equal to m.If
m
is unspecified orm=None
, the quark mass is specified by the running mass itself – that is, it removes itself from the vacuum polarization.The default zeta parameters (from qcdevol.ZETA_M_MSB(self.nf, m, mu)) can be replaced by specifying array
zeta
. Ignored ifzeta=None
(default).
- add_quark(m=None, mu=None, zeta=None)¶
Create a new coupling by adding a heavy quark to the vac. polarization.
Creates a new mass with
nf+1
flavors by adding a quark whose(nf+1)
-flavor mass is m at scale mu, wherenf
equalsself.alpha.nf
. If mu is omitted (or None), it is set equal to m.If
m
is unspecified orm=None
, the quark mass is specified by the running mass itself – that is, it adds itself to the vacuum polarization.The default zeta parameters (from
qcdevol.ZETA_M_MSB(nf+1, m, mu)
) can be replaced by specifying arrayzeta
. Ignored ifzeta=None
(default).
- clone(mu0=None)¶
Create clone of
self
.If
mu0
is notNone
, parametermu0
is reset in the clone (with the corresponding value form0
).
Functions¶
The following function is used to evolve and manipulate perturbation series:
- qcdevol.evol_ps(pseries, mu, mu0, nf, gamma=None, beta_msb=None, order=None)¶
Evolve power series in
alpha
from scalemu0
to scalemu
.Given a power series in coupling
alpha(mu0)
S(mu0) = sum_n=0..order c_n(mu0) * alpha(mu0)**n
this function calculates the coefficients
c_n(mu)
for the the expansion ofS(mu)
in powers ofalpha(mu)
whered/dt ln(S(mu)) = - alpha(mu) * sum_n=0..N gamma[n] * alpha(mu)**n
and
t=ln(mu**2)
. The new coefficients depend only on the coefficients of the original power series, the ratiomu/mu0
, the numbernf
of quark flavors included in the coupling and the coefficientsgamma[n]
of the anomalous dimension. The initial expansionS(mu0)
is specified bypseries
which is an object of typegvar.powerseries.PowerSeries
(for more information, see the documentation for thegvar
module). The final power seriesS(mu)
is returned as an object of the same type.Typical usage is
>>> import qcdevol as qcd >>> import gvar as gv >>> Z10 = gv.powerseries.PowerSeries([4.], order=5) >>> print(Z10.c) # Z(10) coefficients [4. 0. 0. 0. 0. 0.] >>> Z12 = qcd.evol_ps(Z10, mu=12., mu0=10., nf=4, gamma=[0.318, 0.370, 0.321]) >>> print(Z12.c) # Z(12) coefficients [ 4. -0.46382604 -0.56885921 -0.56718909 -0.18133736 0.34568772]
where
Z12
is the expansion atmu=12
of the quantity whose expansion atmu0=10
is given byZ10
and whose anomalous dimension is specified bygamma
.- Parameters:
pseries – Object of type
gvar.powerseries.PowerSeries
representingS(mu0)
, a series in powers ofalpha(mu0)
with coefficientspseries.c[n]
.mu – Results returned for the series
S(mu)
in terms ofalpha(mu)
. Ifmu
is tuple(mu,scheme)
,S(mu)
is a series in powers ofalpha(mu, scheme)
. The MSbar scheme is assumed ifscheme
is not specified.mu0 – Parameter
pseries
gives the expansion ofS(mu0)
in powers ofalpha(mu0)
. Ifmu0
is tuple(mu0,scheme0)
.S(mu0)
is a series in powers ofalpha(mu0,scheme0)
. The MSbar scheme is assumed ifscheme0
is not specified.nf (int) – Number of quark flavors included in the coupling.
gamma (array) – Coefficients
gamma[n]
of the anomalous dimension.beta_msb (array or None) – If not
None
, replaces the default MSbar beta function (fromqcdevol.BETA_MSB()
). Ignored ifNone
(default).order – Order to which the expansion of
S(mu)
is carried. Iforder=None
,order
is set equal topseries.order
. Larger values should be chosen fororder
if the ratiomu/mu0
is large, to capture contributions involving powers of large logarithms in high orders.
- Returns:
Power series for
S(mu)
in powers ofalpha(mu)
, an object of typegvar.powerseries.PowerSeries
.
The default parameters for the beta function, mass anomalous dimension, and matching formulas (adding and deleting quarks) are obtained from the following functions:
- qcdevol.BETA_MSB(nf, terr=None, alpha_qed=None)¶
Compute
b[n]
for thenf
-flavor MSbar couplingal
:mu**2 dal(mu)/dmu**2 = -b[0]*al**2 -b[1]*al**3 -b[2]*al**4 -b[3]*al**5 -b[4]*al**6
If
terr
(truncation error) is notNone
, a sixth entry is added tob
withb[5]
equal toterr
times the r.m.s. average of the other coefficients. Ignored ifterr=None
(default).
- qcdevol.GAMMA_MSB(nf, terr=None, alpha_qed=None, Q=None)¶
Compute
g[n]
fornf
-flavor MSbar mass:mu**2 dlog(m(mu))/dmu**2= -g[1]*al - g[2]*al**2 - g[3]*al**3 - g[4]*al**4 - g[5]*al**5
where
al = alpha_msb(mu)
.If
terr
(truncation error) is notNone
, a sixth entry is added tog
withg[5]
equal toterr
times the r.m.s. average of the other coefficients. Ignored ifterr=None
(default).
- qcdevol.ZETA2_G_MSB(nf, m, mu, terr=None)¶
Compute
z[i]
fornf
flavor MSbaralpha
where:alpha(mu,nf-1) = z * alpha(mu,nf) z = 1 + z[0]*alpha(mu,nf) + z[1]*alpha(mu,nf)**2 + z[2]*alpha(mu,nf)**3 + z[3]*alpha(mu,nf)**4
Here
m=mh(mu,nf)
is the MSbar for the quark flavor being removed.If
terr
is not none, a fifth entry is added toz
withz[4]
equal toterr
times the r.m.s. average of the other coefficients.
- qcdevol.ZETA_M_MSB(nf, m, mu, terr=None)¶
Compute
z[n]
fornf
flavor MSbar massm
where:m(mu,nf-1) = z * m(mu,nf) z = 1 + z[0]*alpha(mu,nf)**2 + z[1]*alpha(mu,nf)**3 + z[2]*alpha(mu,nf)**4
Here
m=mh(mu,nf)
is the MSbar mass for the quark flavor being removed.If
terr
is not none, a fourth entry is added toz
withz[3]
equal toterr
times the r.m.s. average of the other coefficients.
The dictionary qcd.SCHEMES
defines the coupling schemes
available in qcdevol
. qcd.SCHEMES[s](nf)
returns
a tuple (c,r)
where the coupling als
in
scheme s
is related to the MSbar coupling almsb
by:
als(mu) = almsb(r*mu) * (1 + c[0] * almsb(r*mu) + c[1] * almsb(r*mu)**2 + ...)
The schemes currently available are: s='msb'
for the MSbar
scheme (default scheme for most functions), and s='v'
or s='V'
the static-quark potential coupling, defined by:
V(Q) = -4 pi C_F alv(Q) / Q**2
through third order in the coupling.