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.AbstractDiracMatrixType
abstract 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.

source
QEDbase.AbstractDiracVectorType
abstract type AbstractDiracVector{T} <: StaticArraysCore.FieldVector{4, T}

Abstract type for Dirac vectors, e.g. four dimensional vectors from a spinor space.

source
QEDbase.AbstractFourMomentumType
AbstractFourMomentum

Abstract base type for four-momentas, representing one energy and three spacial components.

Also see: QEDcore.SFourMomentum, QEDcore.MFourMomentum

source
QEDbase.AbstractInvalidInputExceptionType

Abstract 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}
source
QEDbase.AbstractLorentzVectorType
abstract 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.

source
QEDbase.AbstractModelDefinitionType

Abstract 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)
source
QEDbase.AbstractParticleType

Abstract 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.

source
QEDbase.AbstractParticleStatefulType
AbstractParticleStateful <: QEDbase.AbstractParticle

Abstract base type for the representation of a particle with a state. It requires the following interface functions to be provided:

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.

source
QEDbase.AbstractParticleTypeType
AbstractParticleType <: 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..

source
QEDbase.AbstractPhaseSpacePointType
AbstractPhaseSpacePoint{PROC, MODEL, PSDEF, IN_PARTICLES, OUT_PARTICLES}

Representation of a point in the phase space of a process. It has several template arguments:

The following interface functions must be provided:

  • Base.getindex(psp::AbstractPhaseSpacePoint, dir::ParticleDirection, n::Int): Return the nth AbstractParticleStateful of the given direction. Throw BoundsError for invalid indices.
  • particles(psp::AbstractPhaseSpacePoint, dir::ParticleDirection): Return the particle tuple (type IN_PARTICLES or OUT_PARTICLES depending on dir)
  • 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 nth AbstractParticleStateful of the given direction.
  • momenta(psp::PhaseSpacePoint, ::ParticleDirection): Return a Tuple 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 match incoming_particles(::PROC) in length, order, and type or be an empty Tuple.
  • OUT_PARTICLES must match the outgoing_particles(::PROC) in length, order, and type, or be an empty Tuple.
  • IN_PARTICLES and OUT_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.

source
QEDbase.AbstractProcessDefinitionType

Abstract 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.

source
QEDbase.AllPolarizationType

Concrete 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
Alias

There is a built-in alias for AllPolarization:

julia> using QEDbase

julia> AllPolarization === AllPol
true
source
QEDbase.AllSpinType

Concrete 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
source
QEDbase.IncomingType
Incoming <: 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
ParticleDirection Interface

Besides being a subtype of ParticleDirection, Incoming has

is_incoming(::Incoming) = true
is_outgoing(::Incoming) = false
source
QEDbase.OutgoingType
Outgoing <: 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
ParticleDirection Interface

Besides being a subtype of ParticleDirection, Outgoing has

is_incoming(::Outgoing) = false
is_outgoing(::Outgoing) = true
source
QEDbase.ParticleDirectionType

Abstract 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.

source
QEDbase.PolarizationXType

Concrete type which indicates, that a particle with is_boson has polarization in $x$-direction.

julia> using QEDbase

julia> PolX()
x-polarized
Coordinate axes

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.

Alias

There is a built-in alias for PolarizationX:

julia> using QEDbase

julia> PolarizationX === PolX
true
source
QEDbase.PolarizationYType

Concrete type which indicates, that a particle with is_boson has polarization in $y$-direction.

julia> using QEDbase

julia> PolY()
y-polarized
Coordinate axes

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.

Alias

There is a built-in alias for PolarizationY:

julia> using QEDbase

julia> PolarizationY === PolY
true
source
QEDbase.RegistryErrorType
struct 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

source
QEDbase.UnknownDirectionType
UnknownDirection <: 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
ParticleDirection Interface

Besides being a subtype of ParticleDirection, UnknownDirection has

is_incoming(::UnknownDirection) = false
is_outgoing(::UnknownDirection) = false
source
QEDbase._as_svecFunction
_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.

source
QEDbase._averaging_normFunction
_averaging_norm(proc::AbstractProcessDefinition)

Interface function, which returns a normalization for the averaging of the squared matrix elements over spins and polarizations.

source
QEDbase._generate_incoming_momentaFunction
_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.

source
QEDbase._generate_momentaMethod
_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.

source
QEDbase._generate_outgoing_momentaFunction
_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.

source
QEDbase._hasmethod_registryMethod
_hasmethod_registry(fun, _)

Wrapper around Base.hasmethod with a more meaningful error message in the context of function registration.

source
QEDbase._incident_fluxFunction
_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.

source
QEDbase._is_in_phasespaceFunction
_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.

source
QEDbase._matrix_elementFunction
_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.

source
QEDbase._phase_space_factorFunction
_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.

Convention

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)$.

source
QEDbase._spin_indexMethod
_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.

source
QEDbase._total_probabilityFunction
_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.

total cross section

Given an implementation of this method and _incident_flux, the respective function for the total cross section total_cross_section is also available.

source
QEDbase.assert_onshellMethod
assert_onshell(mom, mass)

Assertion if a FourMomentum mom is on-shell w.r.t a given mass mass.

See also

The precision of this functions is explained in isonshell.

source
QEDbase.base_stateFunction
    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 an AbstractParticleType, e.g. QEDcore.Electron, QEDcore.Positron, or QEDcore.Photon.
  • direction – the direction of the particle, i.e. Incoming or Outgoing.
  • momentum – the four-momentum of the particle
  • [spin_or_pol] – if given, the spin or polarization of the particle, e.g. SpinUp/SpinDown or PolarizationX/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]
Iterator convenience

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
Conventions

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*}\]

Warning

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.

source
QEDbase.chargeFunction
charge(::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.

source
QEDbase.differential_cross_sectionMethod
differential_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.

source
QEDbase.differential_probabilityMethod
differential_probability(phase_space_point::AbstractPhaseSpacePoint)

If the given phase spaces are physical, return differential probability evaluated on a phase space point. Zero otherwise.

source
QEDbase.getBetaMethod
getBeta(lv)

Return magnitude of the beta vector for a given LorentzVectorLike, i.e. the magnitude of the LorentzVectorLike divided by its 0-component.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to sqrt(px^2 + py^2 + pz^2)/E.

source
QEDbase.getCosPhiMethod
getCosPhi(lv)

Return the cosine of the phi angle for a given LorentzVectorLike.

Note

This is an equivalent but faster version of cos(getPhi(lv)); see getPhi.

source
QEDbase.getCosThetaMethod
getCosTheta(lv)

Return the cosine of the theta angle for a given LorentzVectorLike.

Note

This is an equivalent but faster version of cos(getTheta(lv)); see getTheta.

source
QEDbase.getEMethod
getE(lv)

Return the energy component of a given LorentzVectorLike, i.e. its 0-component.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to E.

source
QEDbase.getGammaMethod
getGamma(lv)

Return the relativistic gamma factor for a given LorentzVectorLike, i.e. the inverse square root of one minus the beta vector squared.

Example

If (E,px,py,pz) is a LorentzVectorLike with beta vector β, this is equivalent to 1/sqrt(1- β^2).

source
QEDbase.getInvariantMassMethod
getInvariantMass(lv)

Return the the invariant mass of a given LorentzVectorLike, i.e. the square root of the minkowski dot with itself.

Example

If (t,x,y,z) is a LorentzVectorLike, this is equivalent to sqrt(t^2 - (x^2 + y^2 + z^2)).

Note

If the squared invariant mass m2 is negative, -sqrt(-m2) is returned.

source
QEDbase.getInvariantMass2Method
getInvariantMass2(lv)

Return the squared invariant mass of a given LorentzVectorLike, i.e. the minkowski dot with itself.

Example

If (t,x,y,z) is a LorentzVectorLike, this is equivalent to t^2 - (x^2 + y^2 + z^2).

source
QEDbase.getMagnitudeMethod
getMagnitude(lv)

Return the magnitude of a given LorentzVectorLike, i.e. the euklidian norm spatial components.

Example

If (t,x,y,z) is a LorentzVectorLike, this is equivalent to sqrt(x^2 + y^2 + z^2).

Warning

This function differs from a similar function for the TLorentzVector used in the famous ROOT library.

source
QEDbase.getMagnitude2Method
getMagnitude2(lv)

Return the square of the magnitude of a given LorentzVectorLike, i.e. the sum of the squared spatial components.

Example

If (t,x,y,z) is a LorentzVectorLike, this is equivalent to x^2+ y^2 + z^2.

Warning

This function differs from a similar function for the TLorentzVector used in the famous ROOT library.

source
QEDbase.getPhiMethod
getPhi(lv)

Return the phi angle for a given LorentzVectorLike, i.e. the azimuthal angle of its spatial components in spherical coordinates.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to atan(py,px).

Note

The spherical coordinates are defined w.r.t. to the 3-axis.

source
QEDbase.getPxMethod
getPx(lv)

Return the $p_x$ component of a given LorentzVectorLike, i.e. its 1-component.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to px.

source
QEDbase.getPyMethod
getPy(lv)

Return the $p_y$ component of a given LorentzVectorLike, i.e. its 2-component.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to py.

source
QEDbase.getPzMethod
getPz(lv)

Return the $p_z$ component of a given LorentzVectorLike, i.e. its 3-component.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to pz.

source
QEDbase.getRapidityMethod
getRapidity(lv)

Return the rapidity for a given LorentzVectorLike.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to 0.5*log((E+pz)/(E-pz)).

Note

The transverse components are defined w.r.t. to the 3-axis.

source
QEDbase.getSinPhiMethod
getSinPhi(lv)

Return the sine of the phi angle for a given LorentzVectorLike.

Note

This is an equivalent but faster version of sin(getPhi(lv)); see getPhi.

source
QEDbase.getTFunction
getT(lv)

Return the 0-component of a given LorentzVectorLike.

Example

If (t,x,y,z) is a LorentzVectorLike, this is equivalent to t.

source
QEDbase.getThetaMethod
getTheta(lv)

Return the theta angle for a given LorentzVectorLike, i.e. the polar angle of its spatial components in spherical coordinates.

Example

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).

Note

The spherical coordinates are defined w.r.t. to the 3-axis.

source
QEDbase.getTransverseMassMethod
getTransverseMass(lv)

Return the transverse momentum for a given LorentzVectorLike, i.e. the square root of its squared transverse mass.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to sqrt(E^2 - pz^2).

Note

The transverse components are defined w.r.t. to the 3-axis.

Note

If the squared transverse mass mT2 is negative, -sqrt(-mT2) is returned.

source
QEDbase.getTransverseMass2Method
getTransverseMass2(lv)

Return the squared transverse mass for a given LorentzVectorLike, i.e. the difference of its squared 0- and 3-component.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to E^2 - pz^2.

Note

The transverse components are defined w.r.t. to the 3-axis.

source
QEDbase.getTransverseMomentumMethod
getTransverseMomentum(lv)

Return the transverse momentum for a given LorentzVectorLike, i.e. the magnitude of its transverse components.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to sqrt(px^2 + py^2).

Note

The transverse components are defined w.r.t. to the 3-axis.

source
QEDbase.getTransverseMomentum2Method
getTransverseMomentum2(lv)

Return the squared transverse momentum for a given LorentzVectorLike, i.e. the sum of its squared 1- and 2-component.

Example

If (E,px,py,pz) is a LorentzVectorLike, this is equivalent to px^2 + py^2.

Note

The transverse components are defined w.r.t. to the 3-axis.

source
QEDbase.getXFunction
getX(lv)

Return the 1-component of a given LorentzVectorLike.

Example

If (t,x,y,z) is a LorentzVectorLike, this is equivalent to x.

source
QEDbase.getYFunction
getY(lv)

Return the 2-component of a given LorentzVectorLike.

Example

If (t,x,y,z) is a LorentzVectorLike, this is equivalent to y.

source
QEDbase.getZFunction
getZ(lv)

Return the 3-component of a given LorentzVectorLike.

Example

If (t,x,y,z) is a LorentzVectorLike, this is equivalent to z.

source
QEDbase.incoming_particlesFunction
incoming_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.

source
QEDbase.is_anti_particleMethod
is_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.

source
QEDbase.is_bosonMethod
is_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.

source
QEDbase.is_fermionMethod
is_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.

source
QEDbase.is_particleMethod
is_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.

source
QEDbase.isonshellMethod
isonshell(mom, mass)

On-shell check of a given four-momentum mom w.r.t. a given mass mass.

Precision

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.

FourMomenta with real entries
  • if AbstractFourMomentum is updated to elementtypes T<:Real, the AbstractLorentzVector should be updated with the AbstractFourMomentum.
source
QEDbase.massFunction
mass(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.

source
QEDbase.minkowski_dotMethod
minkowski_dot(v1,v2)

Return the Minkowski dot product of two LorentzVectorLike.

Example

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)
Note

We use the mostly minus metric.

source
QEDbase.modelFunction
model(psp::AbstractPhaseSpacePoint)

Return the phase space point's set model.

source
QEDbase.momentaMethod
momenta(psp::AbstractPhaseSpacePoint, ::ParticleDirection)

Return a Tuple of all the particles' momenta for the given ParticleDirection.

source
QEDbase.momentumMethod
momentum(psp::AbstractPhaseSpacePoint, dir::ParticleDirection, species::AbstractParticleType, n::Int)

Returns the momentum of the nth 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.

Note

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.

source
QEDbase.momentumMethod
momentum(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.

source
QEDbase.momentumMethod
momentum(psp::AbstractPhaseSpacePoint, dir::ParticleDirection, n::Int)

Returns the momentum of the nth 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.

source
QEDbase.momentumMethod
momentum(psp::AbstractPhaseSpacePoint, dir::ParticleDirection, species::AbstractParticleType, n::Val{N})

Returns the momentum of the nth 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.

Note

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).

source
QEDbase.mulMethod
mul(
    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.

Multiplication operator

This also overloads the * operator for this types.

source
QEDbase.mulMethod
mul(
    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.

Multiplication operator

This also overloads the * operator for this types.

source
QEDbase.multiplicityFunction
multiplicity(spin_or_pol)

Return the number of spins or polarizations respresented by spin_or_pol, e.g. multiplicity(SpinUp()) == 1, but multiplicity(AllSpin()) = 2.

source
QEDbase.number_particlesMethod
number_particles(proc_def::AbstractProcessDefinition, particle::AbstractParticleStateful)

Return the number of particles of the given particle's direction and species in the given process definition.

source
QEDbase.number_particlesMethod
number_particles(proc_def::AbstractProcessDefinition, dir::ParticleDirection, species::AbstractParticleType)

Return the number of particles of the given direction and species in the given process definition.

source
QEDbase.outgoing_particlesFunction
outgoing_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.

source
QEDbase.processFunction
process(psp::AbstractPhaseSpacePoint)

Return the phase space point's set process.

source
QEDbase.propagatorFunction
propagator(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.

Convention

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}.\]

Warning

This function does not throw when the given particle is off-shell. If an off-shell particle is passed, the function propagator returns Inf.

source
QEDbase.register_LorentzVectorLikeMethod
register_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.

source
QEDbase.setCosTheta!Method
setCosTheta!(lv,value)

Sets the cosine of the theta angle of a LorentzVectorLike to a given value.

Note

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.

source
QEDbase.setE!Method
setE!(lv,value)

Sets the energy component of a given LorentzVectorLike to a given value.

Note

The value set with setE! is then returned by getE.

source
QEDbase.setMinus!Method
setMinus!(lv,value)

Sets the minus component of a LorentzVectorLike to a given value.

Note

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.

source
QEDbase.setPhi!Method
setPhi!(lv,value)

Sets the phi angle of a LorentzVectorLike to a given value.

Note

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.

source
QEDbase.setPlus!Method
setPlus!(lv,value)

Sets the plus component of a LorentzVectorLike to a given value.

Note

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.

source
QEDbase.setPx!Method
setPx!(lv,value)

Sets the 1-component of a given LorentzVectorLike to a given value.

Note

The value set with setPx! is then returned by getPx.

source
QEDbase.setPy!Method
setPy!(lv,value)

Sets the 2-component of a given LorentzVectorLike to a given value.

Note

The value set with setPy! is then returned by getPy.

source
QEDbase.setPz!Method
setPz!(lv,value)

Sets the 3-component of a given LorentzVectorLike to a given value.

Note

The value set with setPz! is then returned by getPz.

source
QEDbase.setRapidity!Method
setRapidity!(lv,value)

Sets the rapidity of a LorentzVectorLike to a given value.

Note

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.

source
QEDbase.setRho!Method
setRho!(lv,value)

Sets the magnitude of a LorentzVectorLike to a given value.

Note

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.

source
QEDbase.setT!Function
setT!(lv,value)

Sets the 0-component of a given LorentzVectorLike to the given value.

source
QEDbase.setTheta!Method
setTheta!(lv,value)

Sets the theta angle of a LorentzVectorLike to a given value.

Note

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.

source
QEDbase.setTransverseMass!Method
setTransverseMass!(lv,value)

Sets the transverse mass of a LorentzVectorLike to a given value.

Note

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.

source
QEDbase.setTransverseMomentum!Method
setTransverseMomentum!(lv,value)

Sets the transverse momentum of a LorentzVectorLike to a given value.

Note

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.

source
QEDbase.setX!Function
setX!(lv,value)

Sets the 1-component of a given LorentzVectorLike to the given value.

source
QEDbase.setY!Function
setY!(lv,value)

Sets the 2-component of a given LorentzVectorLike to the given value.

source
QEDbase.setZ!Function
setZ!(lv,value)

Sets the 3-component of a given LorentzVectorLike to the given value.

source
QEDbase.unsafe_differential_cross_sectionMethod
unsafe_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.

source
QEDbase.unsafe_differential_probabilityMethod
unsafe_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.

source