QEDbase
Documentation for QEDbase. This library provides (more or less) simple implementations of Lorentz vectors, Dirac spinors and Dirac matrices. These are usually used in order to do calculations in particle physics, e.g. they are part of the set of Feynman rules for a given quantum field theory (see [1]).
QEDbase.AbstractCoordinateSystem
— TypeAbstractCoordinateSystem
TBW
QEDbase.AbstractDefinitePolarization
— TypeAbstract base type for definite polarization of particles with is_boson
.
Concrete types are PolarizationX
and PolarizationY
.
QEDbase.AbstractDefiniteSpin
— TypeAbstract base type for definite spins of particles with is_fermion
.
QEDbase.AbstractDiracMatrix
— Typeabstract type AbstractDiracMatrix{T} <: StaticArraysCore.FieldMatrix{4, 4, T}
Abstract type for Dirac matrices, i.e. matrix representations for linear mappings from a spinor space into another.
QEDbase.AbstractDiracVector
— Typeabstract type AbstractDiracVector{T} <: StaticArraysCore.FieldVector{4, T}
Abstract type for Dirac vectors, e.g. four dimensional vectors from a spinor space.
QEDbase.AbstractFourMomentum
— TypeAbstractFourMomentum
Abstract base type for four-momentas, representing one energy and three spacial components.
Also see: QEDcore.SFourMomentum
, QEDcore.MFourMomentum
QEDbase.AbstractFrameOfReference
— TypeAbstractFrameOfReference
TBW
QEDbase.AbstractInPhaseSpacePoint
— TypeAbstractInPhaseSpacePoint
A partial type specialization on AbstractPhaseSpacePoint
which can be used for dispatch in functions requiring only the in channel of the phase space to exist, for example implementations of _incident_flux
. No restrictions are imposed on the out-channel, which may or may not exist.
See also: AbstractOutPhaseSpacePoint
QEDbase.AbstractIndefinitePolarization
— TypeAbstract base type for indefinite polarization of particles with is_boson
.
One concrete type is AllPolarization
.
QEDbase.AbstractIndefiniteSpin
— TypeAbstract base type for indefinite spins of particles with is_fermion
.
One concrete type is AllSpin
.
QEDbase.AbstractInvalidInputException
— TypeAbstract base type for exceptions indicating invalid input. See InvalidInputError
for a simple concrete implementation. Concrete implementations should at least implement
Base.showerror(io::IO, err::CustomInvalidError) where {CustomInvalidError<:AbstractInvalidInputException}
QEDbase.AbstractLorentzVector
— Typeabstract type AbstractLorentzVector{T} <: StaticArraysCore.FieldVector{4, T}
Abstract type to model generic Lorentz vectors, i.e. elements of a minkowski-like space, where the component space is arbitray.
QEDbase.AbstractModelDefinition
— TypeAbstract base type for all compute model definitions in the context of scattering processes. Every subtype of AbstractModelDefinition
is associated with a fundamental interaction. Therefore, one needs to implement the following soft interface function
fundamental_interaction_type(::AbstractModelDefinition)
QEDbase.AbstractOutPhaseSpacePoint
— TypeAbstractOutPhaseSpacePoint
A partial type specialization on AbstractPhaseSpacePoint
which can be used for dispatch in functions requiring only the out channel of the phase space to exist. No restrictions are imposed on the in-channel, which may or may not exist.
See also: AbstractInPhaseSpacePoint
QEDbase.AbstractParticle
— TypeAbstract base type for every type which might be considered a particle in the context of QED.jl
. For every (concrete) subtype of AbstractParticle
, there are two kinds of interface functions implemented: static functions and property functions. The static functions provide information on what kind of particle it is (defaults are written in square brackets)
is_fermion(::AbstractParticle)::Bool [= false]
is_boson(::AbstractParticle)::Bool [= false]
is_particle(::AbstractParticle)::Bool [= true]
is_anti_particle(::AbstractParticle)::Bool [= false]
If the output of those functions differ from the defaults for a subtype of AbstractParticle
, these functions need to be overwritten. The second type of functions define a hard interface for AbstractParticle
:
mass(::AbstractParticle)::Real
charge(::AbstractParticle)::Real
These functions must be implemented in order to have the subtype of AbstractParticle
work with the functionalities of QEDprocesses.jl
.
QEDbase.AbstractParticleStateful
— TypeAbstractParticleStateful <: QEDbase.AbstractParticle
Abstract base type for the representation of a particle with a state. It requires the following interface functions to be provided:
particle_direction
: Returning the particle's direction.particle_species
: Returning the particle's species.momentum
: Returning the particle's momentum.
Implementations for is_fermion
, is_boson
, is_particle
, is_anti_particle
, is_incoming
, is_outgoing
, mass
, and charge
are automatically provided using the interface functions above to fulfill the QEDbase.AbstractParticle
interface.
QEDbase.AbstractParticleType
— TypeAbstractParticleType <: AbstractParticle
This is the abstract base type for every species of particles. All functionalities defined on subtypes of AbstractParticleType
should be static, i.e. known at compile time. For adding runtime information, e.g. four-momenta or particle states, to a particle, consider implementing a concrete subtype of AbstractParticle
instead, which may have a type parameter P<:AbstractParticleType
.
Concrete built-in subtypes of AbstractParticleType
are available in QEDcore.jl
and should always be singletons..
QEDbase.AbstractPhaseSpacePoint
— TypeAbstractPhaseSpacePoint{PROC, MODEL, PSDEF, IN_PARTICLES, OUT_PARTICLES}
Representation of a point in the phase space of a process. It has several template arguments:
PROC <:
AbstractProcessDefinition
MODEL <:
AbstractModelDefinition
PSDEF <:
AbstractPhasespaceDefinition
IN_PARTICLES <:
Tuple{Vararg{AbstractParticleStateful}}: The tuple type of all the incoming [
AbstractParticleStateful`](@ref)s.OUT_PARTICLES <:
Tuple{Vararg{AbstractParticleStateful}}: The tuple type of all the outgoing [
AbstractParticleStateful`](@ref)s.
The following interface functions must be provided:
Base.getindex(psp::AbstractPhaseSpacePoint, dir::ParticleDirection, n::Int)
: Return the nthAbstractParticleStateful
of the given direction. ThrowBoundsError
for invalid indices.particles(psp::AbstractPhaseSpacePoint, dir::ParticleDirection)
: Return the particle tuple (typeIN_PARTICLES
orOUT_PARTICLES
depending ondir
)process
: Return the process.model
: Return the model.phase_space_definition
: Return the phase space definition.
From this, the following functions are automatically derived:
momentum(psp::AbstractPhaseSpacePoint, dir::ParticleDirection, n::Int)
: Return the momentum of the nthAbstractParticleStateful
of the given direction.momenta(psp::PhaseSpacePoint, ::ParticleDirection)
: Return aTuple
of all the momenta of the given direction.
Furthermore, an implementation of an AbstractPhaseSpacePoint
has to verify on construction that it is valid, i.e., the following conditions are fulfilled:
IN_PARTICLES
must matchincoming_particles(::PROC)
in length, order, and type or be an emptyTuple
.OUT_PARTICLES
must match theoutgoing_particles(::PROC)
in length, order, and type, or be an emptyTuple
.IN_PARTICLES
andOUT_PARTICLES
may not both be empty.
If IN_PARTICLES
is non-empty, AbstractPhaseSpacePoint <: AbstractInPhaseSpacePoint
is true. Likewise, if OUT_PARTICLES
is non-empty, AbstractPhaseSpacePoint <: AbstractOutPhaseSpacePoint
is true. Consequently, if both IN_PARTICLES
and OUT_PARTICLES
are non-empty, both <:
statements are true.
QEDbase.AbstractPhasespaceDefinition
— TypeAbstractPhasespaceDefinition
TBW
QEDbase.AbstractPolarization
— TypeAbstract base type for the polarization of particles with is_boson
.
QEDbase.AbstractProcessDefinition
— TypeAbstract base type for definitions of scattering processes. It is the root type for the process interface, which assumes that every subtype of AbstractProcessDefinition
implements at least
incoming_particles(proc_def::AbstractProcessDefinition)
outgoing_particles(proc_def::AbstractProcessDefinition)
which return a tuple of the incoming and outgoing particles, respectively.
Furthermore, to calculate scattering probabilities and differential cross sections, the following interface functions need to be implemented for every combination of CustomProcess<:AbstractProcessDefinition
, CustomModel<:AbstractModelDefinition
, and CustomPhasespaceDefinition<:AbstractPhasespaceDefinition
.
_incident_flux(psp::InPhaseSpacePoint{CustomProcess,CustomModel})
_matrix_element(psp::PhaseSpacePoint{CustomProcess,CustomModel})
_averaging_norm(proc::CustomProcess)
_is_in_phasespace(psp::PhaseSpacePoint{CustomProcess,CustomModel})
_phase_space_factor(psp::PhaseSpacePoint{CustomProcess,CustomModel,CustomPhasespaceDefinition})
Optional is the implementation of
_total_probability(psp::PhaseSpacePoint{CustomProcess,CustomModel,CustomPhasespaceDefinition})
to enable the calculation of total probabilities and cross sections.
QEDbase.AbstractSpin
— TypeAbstract base type for the spin of particles with is_fermion
.
QEDbase.AbstractSpinOrPolarization
— TypeAbstract base type for the spin or polarization of particles with is_fermion
or is_boson
, respectively.
QEDbase.AllPolarization
— TypeConcrete type indicating that a particle with is_boson
has an indefinite polarization and the differential cross section calculation should average or sum over all polarizations, depending on the direction (Incoming
or Outgoing
) of the particle in question.
julia> using QEDbase
julia> AllPol()
all polarizations
There is a built-in alias for AllPolarization
:
julia> using QEDbase
julia> AllPolarization === AllPol
true
QEDbase.AllSpin
— TypeConcrete type indicating that a particle with is_fermion
has an indefinite spin and the differential cross section calculation should average or sum over all spins, depending on the direction (Incoming
or Outgoing
) of the particle in question.
julia> using QEDbase
julia> AllSpin()
all spins
QEDbase.Incoming
— TypeIncoming <: ParticleDirection
Concrete implementation of a ParticleDirection
to indicate that a particle is incoming in the context of a given process. Mostly used for dispatch.
julia> using QEDbase
julia> Incoming()
incoming
Besides being a subtype of ParticleDirection
, Incoming
has
is_incoming(::Incoming) = true
is_outgoing(::Incoming) = false
QEDbase.InvalidInputError
— TypeInvalidInputError(msg::String)
Exception which is thrown if a function input is invalid.
QEDbase.Outgoing
— TypeOutgoing <: ParticleDirection
Concrete implementation of a ParticleDirection
to indicate that a particle is outgoing in the context of a given process. Mostly used for dispatch.
julia> using QEDbase
julia> Outgoing()
outgoing
Besides being a subtype of ParticleDirection
, Outgoing
has
is_incoming(::Outgoing) = false
is_outgoing(::Outgoing) = true
QEDbase.ParticleDirection
— TypeAbstract base type for the directions of particles in the context of processes, i.e. either they are incoming or outgoing. Subtypes of this are mostly used for dispatch.
QEDbase.PolarizationX
— TypeConcrete type which indicates, that a particle with is_boson
has polarization in $x$-direction.
julia> using QEDbase
julia> PolX()
x-polarized
The notion of axes, e.g. $x$- and $y$-direction is just to distinguish two orthogonal polarization directions. However, if the three-momentum of the particle with is_boson
is aligned to the $z$-axis of a coordinate system, the polarization axes define the $x$- or $y$-axis, respectively.
There is a built-in alias for PolarizationX
:
julia> using QEDbase
julia> PolarizationX === PolX
true
QEDbase.PolarizationY
— TypeConcrete type which indicates, that a particle with is_boson
has polarization in $y$-direction.
julia> using QEDbase
julia> PolY()
y-polarized
The notion of axes, e.g. $x$- and $y$-direction is just to distinguish two orthogonal polarization directions. However, if the three-momentum of the particle with is_boson
is aligned to the $z$-axis of a coordinate system, the polarization axes define the $x$- or $y$-axis, respectively.
There is a built-in alias for PolarizationY
:
julia> using QEDbase
julia> PolarizationY === PolY
true
QEDbase.RegistryError
— Typestruct RegistryError <: Exception
Exception raised, if a certain type target_type
can not be registed for a certain interface, since there needs the function func
to be impleemnted.
Fields
func::Function
target_type::Any
QEDbase.SpinDown
— TypeConcrete type indicating that a particle with is_fermion
has spin-down.
julia> using QEDbase
julia> SpinDown()
spin down
QEDbase.SpinUp
— TypeConcrete type indicating that a particle with is_fermion
has spin-up.
julia> using QEDbase
julia> SpinUp()
spin up
QEDbase.UnknownDirection
— TypeUnknownDirection <: ParticleDirection
Concrete implementation of a ParticleDirection
to indicate that a particle has an unknown direction. This can mean that a specific direction does not make sense in the given context, that the direction is unavailable, or that it is unnecessary.
julia> using QEDbase
julia> UnknownDirection()
unknown direction
Besides being a subtype of ParticleDirection
, UnknownDirection
has
is_incoming(::UnknownDirection) = false
is_outgoing(::UnknownDirection) = false
QEDbase._as_svec
— Function_as_svec(x)
Accepts a single object, an SVector
of objects or a tuple of objects, and returns them in a single "layer" of SVector.
Useful with base_state
.
QEDbase._averaging_norm
— Function_averaging_norm(proc::AbstractProcessDefinition)
Interface function, which returns a normalization for the averaging of the squared matrix elements over spins and polarizations.
QEDbase._generate_incoming_momenta
— Function_generate_incoming_momenta(
proc::AbstractProcessDefinition,
model::AbstractModelDefinition,
phase_space_def::AbstractPhasespaceDefinition,
in_phase_space::NTuple{N,T},
) where {N,T<:Real}
Interface function to generate the four-momenta of the incoming particles from coordinates for a given phase-space definition.
QEDbase._generate_momenta
— Method_generate_momenta(
proc::AbstractProcessDefinition,
model::AbstractModelDefinition,
phase_space_def::AbstractPhasespaceDefinition,
in_phase_space::NTuple{N,T},
out_phase_space::NTuple{M,T},
) where {N,M,T<:Real}
Return four-momenta for incoming and outgoing particles for given coordinate based phase-space points.
QEDbase._generate_outgoing_momenta
— Function_generate_outgoing_momenta(
proc::AbstractProcessDefinition,
model::AbstractModelDefinition,
phase_space_def::AbstractPhasespaceDefinition,
in_phase_space::NTuple{N,T},
out_phase_space::NTuple{M,T},
) where {N,M,T<:Real}
Interface function to generate the four-momenta of the outgoing particles from coordinates for a given phase-space definition.
QEDbase._hasmethod_registry
— Method_hasmethod_registry(fun, _)
Wrapper around Base.hasmethod
with a more meaningful error message in the context of function registration.
QEDbase._incident_flux
— Function_incident_flux(in_psp::InPhaseSpacePoint{PROC,MODEL}) where {
PROC <: AbstractProcessDefinition,
MODEL <: AbstractModelDefinition,
}
Interface function which returns the incident flux of the given scattering process for a given InPhaseSpacePoint
.
QEDbase._is_in_phasespace
— Function_is_in_phasespace(PhaseSpacePoint{PROC,MODEL}) where {
PROC <: AbstractProcessDefinition,
MODEL <: AbstractModelDefinition,
}
Interface function which returns true
if the combination of the given incoming and outgoing phase space is physical, i.e. all momenta are on-shell and some sort of energy-momentum conservation holds.
QEDbase._matrix_element
— Function_matrix_element(PhaseSpacePoint{PROC,MODEL}) where {
PROC <: AbstractProcessDefinition,
MODEL <: AbstractModelDefinition,
}
Interface function which returns a tuple of scattering matrix elements for each spin and polarization combination of proc
.
QEDbase._phase_space_factor
— Function_phase_space_factor(PhaseSpacePoint{PROC,MODEL,PSDEF}) where {
PROC <: AbstractProcessDefinition,
MODEL <: AbstractModelDefinition
PSDEF <: AbstractPhasespaceDefinition,
}
Interface function, which returns the pre-differential factor of the invariant phase space intergral measure.
It is assumed, that this function returns the value of
\[\mathrm{d}\Pi_n:= \prod_{i=1}^N \frac{\mathrm{d}^3p_i}{(2\pi)^3 2 p_i^0} H(P_t, p_1, \dots, p_N),\]
where $H(\dots)$ is a characteristic function (or distribution) which constrains the phase space, e.g. $\delta^{(4)}(P_t - \sum_i p_i)$.
QEDbase._spin_index
— Method_spin_index(_::AbstractDefiniteSpin) -> Int64
Utility function, which converts subtypes of AbstractDefiniteSpin
into an integer representation, e.g. an index:
SpinUp
$\rightarrow 1$SpinDown
$\rightarrow 2$
This is useful in the implementation of some of base_state
's overloads.
QEDbase._total_probability
— Function_total_probability(in_psp::InPhaseSpacePoint{PROC,MODEL}) where {
PROC <: AbstractProcessDefinition,
MODEL <: AbstractModelDefinition,
}
Interface function for the combination of a scattering process and a physical model. Return the total of a given process and model for a passed InPhaseSpacePoint
.
Given an implementation of this method and _incident_flux
, the respective function for the total cross section total_cross_section
is also available.
QEDbase.assert_onshell
— Methodassert_onshell(mom, mass)
Assertion if a FourMomentum mom
is on-shell w.r.t a given mass mass
.
The precision of this functions is explained in isonshell
.
QEDbase.base_state
— Function base_state(
particle::AbstractParticleType,
direction::ParticleDirection,
momentum::QEDbase.AbstractFourMomentum,
[spin_or_pol::AbstractSpinOrPolarization]
)
Return the base state of a directed on-shell particle with a given four-momentum. For internal usage only.
Input
particle
– the type of the particle, i.e., an instance of anAbstractParticleType
, e.g.QEDcore.Electron
,QEDcore.Positron
, orQEDcore.Photon
.direction
– the direction of the particle, i.e.Incoming
orOutgoing
.momentum
– the four-momentum of the particle[spin_or_pol]
– if given, the spin or polarization of the particle, e.g.SpinUp
/SpinDown
orPolarizationX
/PolarizationY
.
Output
The output type of base_state
depends on wether the spin or polarization of the particle passed in is specified or not.
If spin_or_pol
is passed, the output of base_state
is
base_state(::Fermion, ::Incoming, mom, spin_or_pol) # -> BiSpinor
base_state(::AntiFermion, ::Incoming, mom, spin_or_pol) # -> AdjointBiSpinor
base_state(::Fermion, ::Outgoing, mom, spin_or_pol) # -> AdjointBiSpinor
base_state(::AntiFermion, ::Outgoing, mom, spin_or_pol) # -> BiSpinor
base_state(::Photon, ::Incoming, mom, spin_or_pol) # -> SLorentzVector{ComplexF64}
base_state(::Photon, ::Outgoing, mom, spin_or_pol) # -> SLorentzVector{ComplexF64}
If spin_or_pol
is of type AllPolarization
or AllSpin
, the output is an SVector
with both spin/polarization alignments:
base_state(::Fermion, ::Incoming, mom) # -> SVector{2,BiSpinor}
base_state(::AntiFermion, ::Incoming, mom) # -> SVector{2,AdjointBiSpinor}
base_state(::Fermion, ::Outgoing, mom) # -> SVector{2,AdjointBiSpinor}
base_state(::AntiFermion, ::Outgoing, mom) # -> SVector{2,BiSpinor}
base_state(::Photon, ::Incoming, mom) # -> SVector{2,SLorentzVector{ComplexF64}}
base_state(::Photon, ::Outgoing, mom) # -> SVector{2,SLorentzVector{ComplexF64}}
Example
using QEDbase
mass = 1.0 # set electron mass to 1.0
px,py,pz = rand(3) # generate random spatial components
E = sqrt(px^2 + py^2 + pz^2 + mass^2) # compute energy, i.e. the electron is on-shell
mom = SFourMomentum(E, px, py, pz) # initialize the four-momentum of the electron
# compute the state of an incoming electron with spin = SpinUp
# note: base_state is not exported!
electron_state = base_state(QEDcore.Electron(), Incoming(), mom, SpinUp())
julia> using QEDbase; using QEDcore
julia> mass = 1.0; px,py,pz = (0.1, 0.2, 0.3); E = sqrt(px^2 + py^2 + pz^2 + mass^2); mom = SFourMomentum(E, px, py, pz)
4-element SFourMomentum with indices SOneTo(4):
1.0677078252031311
0.1
0.2
0.3
julia> electron_state = base_state(Electron(), Incoming(), mom, SpinUp())
4-element BiSpinor with indices SOneTo(4):
1.4379526505428235 + 0.0im
0.0 + 0.0im
-0.20862995724285552 + 0.0im
-0.06954331908095185 - 0.1390866381619037im
julia> electron_states = base_state(Electron(), Incoming(), mom, AllSpin())
2-element StaticArraysCore.SVector{2, BiSpinor} with indices SOneTo(2):
[1.4379526505428235 + 0.0im, 0.0 + 0.0im, -0.20862995724285552 + 0.0im, -0.06954331908095185 - 0.1390866381619037im]
[0.0 + 0.0im, 1.4379526505428235 + 0.0im, -0.06954331908095185 + 0.1390866381619037im, 0.20862995724285552 + 0.0im]
The returned objects of base_state
can be consistently wrapped in an SVector
for iteration using _as_svec
.
This way, a loop like the following becomes possible when spin
may be definite or indefinite.
for state in QEDbase._as_svec(base_state(Electron(), Incoming(), momentum, spin))
# ...
end
For an incoming fermion with momentum $p$, we use the explicit formula:
\[u_\sigma(p) = \frac{\gamma^\mu p_\mu + m}{\sqrt{\vert p_0\vert + m}} \eta_\sigma,\]
where the elementary base spinors are given as
\[\eta_1 = (1, 0, 0, 0)^T\\ \eta_2 = (0, 1, 0, 0)^T\]
For an outgoing anti-fermion with momentum $p$, we use the explicit formula:
\[v_\sigma(p) = \frac{-\gamma^\mu p_\mu + m}{\sqrt{\vert p_0\vert + m}} \chi_\sigma,\]
where the elementary base spinors are given as
\[\chi_1 = (0, 0, 1, 0)^T\\ \chi_2 = (0, 0, 0, 1)^T\]
For outgoing fermions and incoming anti-fermions with momentum $p$, the base state is given as the Dirac-adjoint of the respective incoming fermion or outgoing anti-fermion state:
\[\overline{u}_\sigma(p) = u_\sigma^\dagger \gamma^0\\ \overline{v}_\sigma(p) = v_\sigma^\dagger \gamma^0\]
where $v_\sigma$ is the base state of the respective outgoing anti-fermion.
For a photon with four-momentum $k^\mu = \omega (1, \cos\phi \sin\theta, \sin\phi \sin\theta, \cos\theta)$, the two polarization vectors are given as
\[\begin{align*} \epsilon^\mu_1 &= (0, \cos\theta \cos\phi, \cos\theta \sin\phi, -\sin\theta),\\ \epsilon^\mu_2 &= (0, -\sin\phi, \cos\phi, 0). \end{align*}\]
In the current implementation there are no checks built-in, which verify the passed arguments whether they describe on-shell particles, i.e. p*p≈mass^2
. Using base_state
with off-shell particles will cause unpredictable behavior.
QEDbase.charge
— Functioncharge(::AbstractParticle)::Real
Interface function for particles. Return the electric charge of a particle (in units of the elementary electric charge).
This needs to be implemented for each concrete subtype of AbstractParticle
.
QEDbase.differential_cross_section
— Methoddifferential_cross_section(phase_space_point::PhaseSpacePoint)
If the given phase spaces are physical, return differential cross section evaluated on a phase space point. Zero otherwise.
QEDbase.differential_probability
— Methoddifferential_probability(phase_space_point::AbstractPhaseSpacePoint)
If the given phase spaces are physical, return differential probability evaluated on a phase space point. Zero otherwise.
QEDbase.fundamental_interaction_type
— Functionfundamental_interaction_type(models_def::AbstractModelDefinition)
Return the fundamental interaction associated with the passed model definition.
QEDbase.getBeta
— MethodgetBeta(lv)
Return magnitude of the beta vector for a given LorentzVectorLike
, i.e. the magnitude of the LorentzVectorLike
divided by its 0-component.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to sqrt(px^2 + py^2 + pz^2)/E
.
QEDbase.getCosPhi
— MethodgetCosPhi(lv)
Return the cosine of the phi angle for a given LorentzVectorLike
.
This is an equivalent but faster version of cos(getPhi(lv))
; see getPhi
.
QEDbase.getCosTheta
— MethodgetCosTheta(lv)
Return the cosine of the theta angle for a given LorentzVectorLike
.
This is an equivalent but faster version of cos(getTheta(lv))
; see getTheta
.
QEDbase.getE
— MethodgetE(lv)
Return the energy component of a given LorentzVectorLike
, i.e. its 0-component.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to E
.
QEDbase.getEnergy
— FunctionFunction alias for getE
.
QEDbase.getGamma
— MethodgetGamma(lv)
Return the relativistic gamma factor for a given LorentzVectorLike
, i.e. the inverse square root of one minus the beta vector squared.
If (E,px,py,pz)
is a LorentzVectorLike
with beta vector β
, this is equivalent to 1/sqrt(1- β^2)
.
QEDbase.getInvariantMass
— MethodgetInvariantMass(lv)
Return the the invariant mass of a given LorentzVectorLike
, i.e. the square root of the minkowski dot with itself.
If (t,x,y,z)
is a LorentzVectorLike
, this is equivalent to sqrt(t^2 - (x^2 + y^2 + z^2))
.
If the squared invariant mass m2
is negative, -sqrt(-m2)
is returned.
QEDbase.getInvariantMass2
— MethodgetInvariantMass2(lv)
Return the squared invariant mass of a given LorentzVectorLike
, i.e. the minkowski dot with itself.
If (t,x,y,z)
is a LorentzVectorLike
, this is equivalent to t^2 - (x^2 + y^2 + z^2)
.
QEDbase.getMag
— FunctionFunctiom alias for getMagnitude
.
QEDbase.getMag2
— FunctionFunctiom alias for getMagnitude2
.
QEDbase.getMagnitude
— MethodgetMagnitude(lv)
Return the magnitude of a given LorentzVectorLike
, i.e. the euklidian norm spatial components.
If (t,x,y,z)
is a LorentzVectorLike
, this is equivalent to sqrt(x^2 + y^2 + z^2)
.
This function differs from a similar function for the TLorentzVector
used in the famous ROOT
library.
QEDbase.getMagnitude2
— MethodgetMagnitude2(lv)
Return the square of the magnitude of a given LorentzVectorLike
, i.e. the sum of the squared spatial components.
If (t,x,y,z)
is a LorentzVectorLike
, this is equivalent to x^2+ y^2 + z^2
.
This function differs from a similar function for the TLorentzVector
used in the famous ROOT
library.
QEDbase.getMass
— FunctionFunction alias for getInvariantMass
.
QEDbase.getMass2
— FunctionFunction alias for getInvariantMass2
QEDbase.getMinus
— MethodgetMinus(lv)
Return the minus component for a given LorentzVectorLike
in light-cone coordinates.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to (E-pz)/2
.
The light-cone coordinates are defined w.r.t. to the 3-axis.
The definition `p^- := (E - p_z)/2
differs from the usual definition of light-cone coordinates in general relativity.
QEDbase.getMt
— FunctionFunction alias for getTransverseMass
QEDbase.getMt2
— FunctionFunction alias for getTransverseMass2
QEDbase.getPerp
— FunctionFunction alias for getTransverseMomentum
.
QEDbase.getPerp2
— FunctionFunction alias for getTransverseMomentum2
.
QEDbase.getPhi
— MethodgetPhi(lv)
Return the phi angle for a given LorentzVectorLike
, i.e. the azimuthal angle of its spatial components in spherical coordinates.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to atan(py,px)
.
The spherical coordinates are defined w.r.t. to the 3-axis.
QEDbase.getPlus
— MethodgetPlus(lv)
Return the plus component for a given LorentzVectorLike
in light-cone coordinates.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to (E+pz)/2
.
The light-cone coordinates are defined w.r.t. to the 3-axis.
The definition `p^+ := (E + p_z)/2
differs from the usual definition of light-cone coordinates in general relativity.
QEDbase.getPt
— FunctionFunction alias for getTransverseMomentum
.
QEDbase.getPt2
— FunctionFunction alias for getTransverseMomentum2
.
QEDbase.getPx
— MethodgetPx(lv)
Return the $p_x$ component of a given LorentzVectorLike
, i.e. its 1-component.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to px
.
QEDbase.getPy
— MethodgetPy(lv)
Return the $p_y$ component of a given LorentzVectorLike
, i.e. its 2-component.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to py
.
QEDbase.getPz
— MethodgetPz(lv)
Return the $p_z$ component of a given LorentzVectorLike
, i.e. its 3-component.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to pz
.
QEDbase.getRapidity
— MethodgetRapidity(lv)
Return the rapidity for a given LorentzVectorLike
.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to 0.5*log((E+pz)/(E-pz))
.
The transverse components are defined w.r.t. to the 3-axis.
QEDbase.getRho
— FunctionFunction alias for getMagnitude
QEDbase.getRho2
— FunctionFunction alias for getMagnitude2
QEDbase.getSinPhi
— MethodgetSinPhi(lv)
Return the sine of the phi angle for a given LorentzVectorLike
.
This is an equivalent but faster version of sin(getPhi(lv))
; see getPhi
.
QEDbase.getT
— FunctiongetT(lv)
Return the 0-component of a given LorentzVectorLike
.
If (t,x,y,z)
is a LorentzVectorLike
, this is equivalent to t
.
QEDbase.getTheta
— MethodgetTheta(lv)
Return the theta angle for a given LorentzVectorLike
, i.e. the polar angle of its spatial components in spherical coordinates.
If (E,px,py,pz)
is a LorentzVectorLike
with magnitude rho
, this is equivalent to arccos(pz/rho)
, which is also equivalent to arctan(sqrt(px^2+py^2)/pz)
.
The spherical coordinates are defined w.r.t. to the 3-axis.
QEDbase.getTransverseMass
— MethodgetTransverseMass(lv)
Return the transverse momentum for a given LorentzVectorLike
, i.e. the square root of its squared transverse mass.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to sqrt(E^2 - pz^2)
.
The transverse components are defined w.r.t. to the 3-axis.
If the squared transverse mass mT2
is negative, -sqrt(-mT2)
is returned.
QEDbase.getTransverseMass2
— MethodgetTransverseMass2(lv)
Return the squared transverse mass for a given LorentzVectorLike
, i.e. the difference of its squared 0- and 3-component.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to E^2 - pz^2
.
The transverse components are defined w.r.t. to the 3-axis.
QEDbase.getTransverseMomentum
— MethodgetTransverseMomentum(lv)
Return the transverse momentum for a given LorentzVectorLike
, i.e. the magnitude of its transverse components.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to sqrt(px^2 + py^2)
.
The transverse components are defined w.r.t. to the 3-axis.
QEDbase.getTransverseMomentum2
— MethodgetTransverseMomentum2(lv)
Return the squared transverse momentum for a given LorentzVectorLike
, i.e. the sum of its squared 1- and 2-component.
If (E,px,py,pz)
is a LorentzVectorLike
, this is equivalent to px^2 + py^2
.
The transverse components are defined w.r.t. to the 3-axis.
QEDbase.getX
— FunctiongetX(lv)
Return the 1-component of a given LorentzVectorLike
.
If (t,x,y,z)
is a LorentzVectorLike
, this is equivalent to x
.
QEDbase.getY
— FunctiongetY(lv)
Return the 2-component of a given LorentzVectorLike
.
If (t,x,y,z)
is a LorentzVectorLike
, this is equivalent to y
.
QEDbase.getZ
— FunctiongetZ(lv)
Return the 3-component of a given LorentzVectorLike
.
If (t,x,y,z)
is a LorentzVectorLike
, this is equivalent to z
.
QEDbase.in_phase_space_dimension
— Functionin_phase_space_dimension(
proc::AbstractProcessDefinition,
model::AbstractModelDefinition,
)
TBW
QEDbase.incoming_particles
— Functionincoming_particles(proc_def::AbstractProcessDefinition)
Interface function for scattering processes. Return a tuple of the incoming particles for the given process definition. This function needs to be given to implement the scattering process interface.
QEDbase.is_anti_particle
— Methodis_anti_particle(_::AbstractParticle) -> Bool
Interface function for particles. Return true if the passed subtype of AbstractParticle
can be considered an anti particle as distinct from their particle counterpart, and false
otherwise. The default implementation of is_anti_particle
for every subtype of AbstractParticle
will always return false
.
QEDbase.is_boson
— Methodis_boson(_::AbstractParticle)
Interface function for particles. Return true
if the passed subtype of AbstractParticle
can be considered a boson in the sense of particle statistics, and false
otherwise. The default implementation of is_boson
for every subtype of AbstractParticle
will always return false
.
QEDbase.is_fermion
— Methodis_fermion(_::AbstractParticle) -> Bool
Interface function for particles. Return true
if the passed subtype of AbstractParticle
can be considered a fermion in the sense of particle statistics, and false
otherwise. The default implementation of is_fermion
for every subtype of AbstractParticle
will always return false
.
QEDbase.is_incoming
— Functionis_incoming(dir::ParticleDirection)
is_incoming(particle::AbstractParticleStateful)
Convenience function that returns true for Incoming
and incoming AbstractParticleStateful
and false otherwise.
QEDbase.is_outgoing
— Functionis_outgoing(dir::ParticleDirection)
is_outgoing(particle::AbstractParticleStateful)
Convenience function that returns true for Outgoing
and outgoing AbstractParticleStateful
and false otherwise.
QEDbase.is_particle
— Methodis_particle(_::AbstractParticle) -> Bool
Interface function for particles. Return true
if the passed subtype of AbstractParticle
can be considered a particle as distinct from anti-particles, and false
otherwise. The default implementation of is_particle
for every subtype of AbstractParticle
will always return true
.
QEDbase.isonshell
— Methodisonshell(mom, mass)
On-shell check of a given four-momentum mom
w.r.t. a given mass mass
.
For AbstactFourMomentum
, the element type is fixed to Float64
, limiting the precision of comparisons between elements. The current implementation has been tested within the boundaries for energy scales E
with 1e-9 <= E <= 1e5
. In those bounds, the mass error, which is correctly detected as off-shell, is 1e-4
times the mean value of the components, but has at most the value 0.01
, e.g. at the high energy end. The energy scales correspond to 0.5 meV
for the lower bound and 50 GeV
for the upper bound.
- if
AbstractFourMomentum
is updated to elementtypesT<:Real
, theAbstractLorentzVector
should be updated with theAbstractFourMomentum
.
QEDbase.mass
— Functionmass(particle::AbstractParticle)::Real
Interface function for particles. Return the rest mass of a particle (in units of the electron mass).
This needs to be implemented for each concrete subtype of AbstractParticle
.
QEDbase.mdot
— FunctionFunction alias for minkowski_dot
.
QEDbase.minkowski_dot
— Methodminkowski_dot(v1,v2)
Return the Minkowski dot product of two LorentzVectorLike
.
If (t1,x1,y1,z1)
and (t2,x2,y2,z2)
are two LorentzVectorLike
, this is equivalent to
t1*t2 - (x1*x2 + y1*y2 + z1*z2)
We use the mostly minus metric.
QEDbase.model
— Functionmodel(psp::AbstractPhaseSpacePoint)
Return the phase space point's set model.
QEDbase.momenta
— Methodmomenta(psp::AbstractPhaseSpacePoint, ::ParticleDirection)
Return a Tuple
of all the particles' momenta for the given ParticleDirection
.
QEDbase.momentum
— Functionmomentum(part::AbstractParticleStateful)
Interface function that must return the particle's AbstractFourMomentum
.
QEDbase.momentum
— Methodmomentum(psp::AbstractPhaseSpacePoint, dir::ParticleDirection, species::AbstractParticleType, n::Int)
Returns the momentum of the n
th particle in the given AbstractPhaseSpacePoint
which has direction dir
and species species
. If n
is outside the valid range for this phase space point, a BoundsError
is thrown.
This function accepts n as an Int
value. If n
is a compile-time constant (for example, a literal 1
or 2
), you can use Val(n)
instead to call a zero overhead version of this function.
QEDbase.momentum
— Methodmomentum(psp::AbstractPhaseSpacePoint, dir::ParticleDirection, species::AbstractParticleType)
Returns the momentum of the particle in the given AbstractPhaseSpacePoint
with dir
and species
, if there is only one such particle. If there are multiple or none, an InvalidInputError
is thrown.
QEDbase.momentum
— Methodmomentum(psp::AbstractPhaseSpacePoint, dir::ParticleDirection, n::Int)
Returns the momentum of the n
th particle in the given AbstractPhaseSpacePoint
which has direction dir
. If n
is outside the valid range for this phase space point, a BoundsError
is thrown.
QEDbase.momentum
— Methodmomentum(psp::AbstractPhaseSpacePoint, dir::ParticleDirection, species::AbstractParticleType, n::Val{N})
Returns the momentum of the n
th particle in the given AbstractPhaseSpacePoint
which has direction dir
and species species
. If n
is outside the valid range for this phase space point, a BoundsError
is thrown.
This function accepts n as a Val{N}
type, i.e., a compile-time constant value (for example a literal 1
or 2
). This allows this function to add zero overhead, but only if N
is actually known at compile time. If it is not, use the overload of this function that uses n::Int
instead. That function is faster than calling this one with Val(n)
.
QEDbase.mul
— Methodmul(
DM::Union{AbstractDiracMatrix, AbstractDiracVector},
L::AbstractLorentzVector
) -> Any
Product of generic Lorentz vector with a Dirac tensor from the left. Basically, the multiplication is piped to the components from the Lorentz vector.
This also overloads the *
operator for this types.
QEDbase.mul
— Methodmul(
L::AbstractLorentzVector,
DM::Union{AbstractDiracMatrix, AbstractDiracVector}
) -> Any
Product of generic Lorentz vector with a Dirac tensor from the right. Basically, the multiplication is piped to the components from the Lorentz vector.
This also overloads the *
operator for this types.
QEDbase.multiplicity
— Functionmultiplicity(spin_or_pol)
Return the number of spins or polarizations respresented by spin_or_pol
, e.g. multiplicity(SpinUp()) == 1
, but multiplicity(AllSpin()) = 2
.
QEDbase.number_incoming_particles
— Methodnumber_incoming_particles(proc_def::AbstractProcessDefinition)
Return the number of incoming particles of a given process.
QEDbase.number_outgoing_particles
— Methodnumber_outgoing_particles(proc_def::AbstractProcessDefinition)
Return the number of outgoing particles of a given process.
QEDbase.number_particles
— Methodnumber_particles(proc_def::AbstractProcessDefinition, particle::AbstractParticleStateful)
Return the number of particles of the given particle's direction and species in the given process definition.
QEDbase.number_particles
— Methodnumber_particles(proc_def::AbstractProcessDefinition, dir::ParticleDirection)
Convenience function dispatching to number_incoming_particles
or number_outgoing_particles
depending on the given direction, returning the number of incoming or outgoing particles, respectively.
QEDbase.number_particles
— Methodnumber_particles(proc_def::AbstractProcessDefinition, dir::ParticleDirection, species::AbstractParticleType)
Return the number of particles of the given direction and species in the given process definition.
QEDbase.out_phase_space_dimension
— Functionout_phase_space_dimension(
proc::AbstractProcessDefinition,
model::AbstractModelDefinition,
)
TBW
QEDbase.outgoing_particles
— Functionoutgoing_particles(proc_def::AbstractProcessDefinition)
Interface function for scattering processes. Return the tuple of outgoing particles for the given process definition. This function needs to be given to implement the scattering process interface.
QEDbase.particle_direction
— Functionparticle_direction(part::AbstractParticleStateful)
Interface function that must return the particle's ParticleDirection
.
QEDbase.particle_species
— Functionparticle_species(part::AbstractParticleStateful)
Interface function that must return the particle's AbstractParticleType
, e.g. QEDcore.Electron()
.
QEDbase.particles
— Methodparticles(proc_def::AbstractProcessDefinition, ::ParticleDirection)
Convenience function dispatching to incoming_particles
or outgoing_particles
depending on the given direction.
QEDbase.phase_space_definition
— Functionphase_space_definition(psp::AbstractPhaseSpacePoint)
Return the phase space point's set phase space definition.
QEDbase.process
— Functionprocess(psp::AbstractPhaseSpacePoint)
Return the phase space point's set process.
QEDbase.propagator
— Functionpropagator(particle::AbstractParticleType, mom::QEDbase.AbstractFourMomentum, [mass::Real])
Return the propagator of a particle for a given four-momentum. If mass
is passed, the respective propagator for massive particles is used, if not, it is assumed the particle passed in is massless.
There are two types of implementations for propagators given in QEDProcesses
: For a BosonLike
particle with four-momentum $k$ and mass $m$, the propagator is given as
\[D(k) = \frac{1}{k^2 - m^2}.\]
For a FermionLike
particle with four-momentum $p$ and mass $m$, the propagator is given as
\[S(p) = \frac{\gamma^\mu p_\mu + mass}{p^2 - m^2}.\]
This function does not throw when the given particle is off-shell. If an off-shell particle is passed, the function propagator
returns Inf
.
QEDbase.register_LorentzVectorLike
— Methodregister_LorentzVectorLike(T)
Function to register a custom type as a LorentzVectorLike.
Ensure the passed custom type has implemented at least the function getT, getX, getY, getZ
and enables getter functions of the lorentz vector library for the given type. If additionally the functions setT!, setX!, setY!, setZ!
are implemened for the passed custom type, also the setter functions of the Lorentz vector interface are enabled.
QEDbase.setCosTheta!
— MethodsetCosTheta!(lv,value)
Sets the cosine of the theta angle of a LorentzVectorLike
to a given value
.
The value
set with setCosTheta!
is then returned by getCosTheta
. Since the cosine of the theta angle is computed on the call of getCosTheta
, the setter setCosTheta!
changes several components of the given LorentzVectorLike
.
QEDbase.setE!
— MethodsetE!(lv,value)
Sets the energy component of a given LorentzVectorLike
to a given value
.
The value
set with setE!
is then returned by getE
.
QEDbase.setEnergy!
— FunctionFunction alias for setE!
.
QEDbase.setMinus!
— MethodsetMinus!(lv,value)
Sets the minus component of a LorentzVectorLike
to a given value
.
The value
set with setMinus!
is then returned by getMinus
. Since the minus component is computed on the call of getMinus
, the setter setMinus!
changes several components of the given LorentzVectorLike
.
QEDbase.setMt!
— FunctionFunction alias for setTransverseMass!
.
QEDbase.setPerp!
— FunctionFunction alias for setTransverseMomentum!
.
QEDbase.setPhi!
— MethodsetPhi!(lv,value)
Sets the phi angle of a LorentzVectorLike
to a given value
.
The value
set with setPhi!
is then returned by getPhi
. Since the phi angle is computed on the call of getPhi
, the setter setPhi!
changes several components of the given LorentzVectorLike
.
QEDbase.setPlus!
— MethodsetPlus!(lv,value)
Sets the plus component of a LorentzVectorLike
to a given value
.
The value
set with setPlus!
is then returned by getPlus
. Since the plus component is computed on the call of getPlus
, the setter setPlus!
changes several components of the given LorentzVectorLike
.
QEDbase.setPt!
— FunctionFunction alias for setTransverseMomentum!
.
QEDbase.setPx!
— MethodsetPx!(lv,value)
Sets the 1-component of a given LorentzVectorLike
to a given value
.
The value
set with setPx!
is then returned by getPx
.
QEDbase.setPy!
— MethodsetPy!(lv,value)
Sets the 2-component of a given LorentzVectorLike
to a given value
.
The value
set with setPy!
is then returned by getPy
.
QEDbase.setPz!
— MethodsetPz!(lv,value)
Sets the 3-component of a given LorentzVectorLike
to a given value
.
The value
set with setPz!
is then returned by getPz
.
QEDbase.setRapidity!
— MethodsetRapidity!(lv,value)
Sets the rapidity of a LorentzVectorLike
to a given value
.
The value
set with setRapidity!
is then returned by getRapidity
. Since the rapidity is computed on the call of setRapidity
, the setter setRapidity!
changes several components of the given LorentzVectorLike
.
QEDbase.setRho!
— MethodsetRho!(lv,value)
Sets the magnitude of a LorentzVectorLike
to a given value
.
The value
set with setRho!
is then returned by getRho
. Since the magnitude is computed on the call of getRho
, the setter setRho!
changes several components of the given LorentzVectorLike
.
QEDbase.setT!
— FunctionsetT!(lv,value)
Sets the 0-component of a given LorentzVectorLike
to the given value
.
QEDbase.setTheta!
— MethodsetTheta!(lv,value)
Sets the theta angle of a LorentzVectorLike
to a given value
.
The value
set with setTheta!
is then returned by getTheta
. Since the theta angle is computed on the call of getTheta
, the setter setTheta!
changes several components of the given LorentzVectorLike
.
QEDbase.setTransverseMass!
— MethodsetTransverseMass!(lv,value)
Sets the transverse mass of a LorentzVectorLike
to a given value
.
The value
set with setTransverseMass!
is then returned by getTransverseMass
. Since the transverse mass is computed on the call of getTransverseMass
, the setter setTransverseMass!
changes several components of the given LorentzVectorLike
.
QEDbase.setTransverseMomentum!
— MethodsetTransverseMomentum!(lv,value)
Sets the transverse momentum of a LorentzVectorLike
to a given value
.
The value
set with setTransverseMomentum!
is then returned by getTransverseMomentum
. Since the transverse momentum is computed on the call of getTransverseMomentum
, the setter setTransverseMomentum!
changes several components of the given LorentzVectorLike
.
QEDbase.setX!
— FunctionsetX!(lv,value)
Sets the 1-component of a given LorentzVectorLike
to the given value
.
QEDbase.setY!
— FunctionsetY!(lv,value)
Sets the 2-component of a given LorentzVectorLike
to the given value
.
QEDbase.setZ!
— FunctionsetZ!(lv,value)
Sets the 3-component of a given LorentzVectorLike
to the given value
.
QEDbase.total_cross_section
— Methodtotal_cross_section(in_psp::AbstractInPhaseSpacePoint)
Return the total cross section for a given AbstractInPhaseSpacePoint
.
QEDbase.total_probability
— Methodtotal_probability(in_psp::AbstractInPhaseSpacePoint)
Return the total probability of a given AbstractInPhaseSpacePoint
.
QEDbase.unsafe_differential_cross_section
— Methodunsafe_differential_cross_section(phase_space_point::AbstractPhaseSpacePoint)
Return the differential cross section evaluated on a phase space point without checking if the given phase space is physical.
QEDbase.unsafe_differential_probability
— Methodunsafe_differential_probability(phase_space_point::AbstractPhaseSpacePoint)
Return differential probability evaluated on a phase space point without checking if the given phase space(s) are physical.