The Spacetime Class Library

Joe Boudreau 2017

Table of Contents


1 Introduction

The Spacetime Class Library is a small collection of classes intended to facilitate numerical computations in nonrelativistic and relativistic quantum mechanics. The library contains rotations and Lorentz transformation classes, instances of which can provide matrix representations of the corresponding group element. Objects such as vectors, spinors, and Dirac spinors (among others) which belong to some representation of the rotation group (or the Lorentz group) are transformed accordingly. The design allows other types of covariant objects to be introduced into the library at a later time. Presently, the library contains vectors (v⃗⃗), four-vectors (pμ), spinors (|↑⟩, |↓⟩), Dirac spinors (u,v,u, v) and Weyl spinors, and special operators such as Pauli matrices (σ⃗ = σx, σyσz), gamma matrices, γ0, γ1, γ2, γ3, and the matrices Sμν = (i)/(4)[γμ, γν]. The usual operations are defined on these classes. The library depends upon the Eigen library for matrix manipulation [A]  [A] http://eigen.tuxfamily.org. This documentation is the reference manual for the Spacetime library.

2 Rotations and Lorentz Tranformations

2.1 class Rotation

#include “Spacetime/Rotation.h”
Description
A rotation is specfied by an axis and an angle, interpreted as an active right-handed rotation. For example, a rotation about the +z axis carries the x axis into the y-axis and the y axis into the  − x-axis. A rule of composition applies to rotations and is expressed with the multiplication symbol (operator *). The rotation also provides a matrix in any desired D = 2j + 1 dimensional irreducible representation of SU(2), expressed in the basis{e⃗0 = |j, j, e⃗1 = |j, j − 1⟩.. e⃗2j = |j,  − j⟩}. Here, j is the angular momentum quantum number, integer or half-integer. Rotations can also act on other objects like vectors and spinors, but this action is defined in free subroutines, which are discussed in conjunction with vectors, spinors, and other mathematical objects upon which the action of the rotation is defined. These are described later.
Methods
Constructor. Constructs an identity rotation:
Constructor. Constructs a rotation from a vector. The length of the vector represents the angle of rotation, in radians. the direction represents the axis of rotation.
Composition
Retrieve the rotation matrix in the D-dimensional representation of SU(2). Here D = 2j + 1 where j is the angular momentum quantum number, integer or half-integer.
Get the rotation vector, whose length is equal to the angle of rotation in radians and whose direction represents the axis of rotation (right-handed, active).
Return the inverse rotation:


2.2 class LorentzTransformation

#include “Spacetime/LorentzTransformation.h
Description
The LorentzTransformation class represents a proper, orthochronous Lorentz transformation. It is specified by two vectors, the first representing the boost and the second representing the rotation vector. Both are to be considered as active, i.e. acting on the body rather than the coordinate system. The boost vector points along the direction of the boost and has a magnitude of η = tanh − 1β (Achtung!) where β = v ⁄ c, v is the velocity and c is the speed of light. The rotation vector is described in section 2.1↑. Irreducible representations of the proper orthochronous Lorentz group are equivalent to SU(2) × SU(2),  and accordingly two matrices are given, one for each of the two SU(2) subgroups. Composition of LorentzTransformations can be carried out with the * operator.
Methods
Constructor. Constructs an identity element.
Constructor. Takes the rapidity vector and the rotation vector
Composition.
Get the D-dimensional representation, matrix 1 and matrix 2
Get the rotation vector and the rapidity vector
Get the inverse.


3 Vectors, Four-Vectors, and Tensors

This section describes vectors in Cartesian 3-space, 4-vectors in Minkokski space, and the assymetric four-tensor, which is an tensor having simple properties under Lorentz tranformations. These objects may be transformed by rotations and/or Lorentz transformations. The work is done in free subroutines, which are also described in this section.

3.1 template <class T> BasicThreeVector

#include “Spacetime/ThreeVector.h”
Description
A BasicThreeVector is the classic vector in Cartesian 3-space. It is affected by rotations. It is implemented as a template class because in addition to real vectors, we sometimes need complex vectors, to describe certain polarization states, for example. Two parameterized classes are also defined in the header, ThreeVector and ComplexThreeVector.
Methods
Constructor. Constructs a zero vector
Constructor. Constructs from components
Access to elements
Access to elements:
Cross and dot:
Compound operations
Unary minus
Conjugate
Norm
Squared norm:
Normalized version of four-vector.
Parameterized types
Other operations
Rotation
Formatted I/O
Vector addition and subtraction
Scaling:


3.2 template class<T> BasicFourVector

#include “Spacetime/FourVector.h”
Description
A BasicFourVector is a four-component vector in Minkowski space. It is affected by Lorentz transformations. It is implemented as a template class because in addition to real four-vectors, we sometimes need complex four-vectors, to describe circular polarization states. Two parameterized classes are also defined in the header, FourVector and ComplexFourVector.
Methods
Construct a null four vector
Construct from four components
Access to elements
Access to elements
Return the space portion
Returns the invariant interval
Compound operations
Unary minus
Conjugate
Norm
Squared norm:
Parameterized types
Other operations
Lorentz Transformations
Formatted I/O
Four-vector addition and subtraction
Scaling:


3.3 class AntisymmetricFourTensor

#include “Spacetime/AntisymmetricFourTensor.h”
Description
An antisymmetric four-tensor is a tensor in Minkowski space with simple transformation properties under Lorentz transformations. This class has the property that modifying any one of its elements Tij automatically changes Tji to  − Tij; this works thanks to the helper class ddouble which you can be used exactly as a double. This works because of a cast operator that permits automatic type conversion. The best known example of an antisymmetric four-tensor is the Maxwell field tensor.
Methods
Construct a null tensor:
Access to elements. Class ddouble is used just like a double.
Compound operations
Unary minus:
Other operations
Formatted output:
Arithmetic operations:
Lorentz Transformation:


4 Spinors

This section describes three types of spinors (in a D-dimensional Hilbert space), Dirac spinors, and Weyl spinors.

4.1 template <unsigned int D=2> Spinor

#include “Spacetime/Spinor.h”
Description
This class represents a spinor in a D = 2j + 1 dimensional space. By default D = 2 and the spinor describes spin j = 1 ⁄ 2 particles.
Methods
Construct a null spinor:
Construct a spinor from an initializer list of std::complex<double>:
Construct a spinor from other elgible datatypes from the Eigen library:
Assign Eigen expressions to the spinor:
Other operations
Action of Rotation upon Spinor:
Recover the spin from the Spinor:


4.2 Weyl spinor classes

#include “Spacetime/WeylSpinor.h”
Description
Weyl spinors are two-component spinors representing massless particles, transforming under rotations and Lorentz transformations. The following datatypes are defined within the namespace WeylSpinor:
Operations
Action of Rotation upon the two types of spinor:
Return the Four-momentum of this Weyl Spinor
Return the Spin (magnitude and direction)

4.2.1 template <int T> class BasicSpinor

The basic 2-component Weyl spinor
Inherits
Eigen::Vector2cd
Methods
Constructs a null spinor:
Constructor from momentum vector:
Construct a spinor from other elgible datatypes from the Eigen library:
Assign a spinor from other elgible datatypes from the Eigen library:


4.3 Dirac spinor classes

#include “Spacetime/DiracSpinor.h”
Description
Dirac spinors describe four-component complex column vectors with specific transformation properties under rotations and Lorentz transformations. In this library the internal representation is the Weyl, or chiral, representation. The following datatypes are defined in the namespace Dirac:
Operations
Take the spinor adjoint (spinor-bar):
Action of Lorentz transformation upon the spinor
Return the four-momentum of this Dirac spinor
Return the spin (magnitude and direction)

4.3.1 template<int T> BasicSpinor

The basic 4-component column spinor.
Inherits
Eigen::Vector4cd
Methods:
Construct a null Dirac spinor:
Construct a Dirac spinor of a particle at rest with spin described by the Spinor s and mass m :
Construct a Dirac spinor of a particle with spin described by the Spinor s and with four-momentum p:
Construct a Dirac spinor for a massless particle from the corresponding Weyl spinor:
Construct a Dirac spinor from other elgible datatypes from the Eigen library:
Assign a Dirac spinor from other elgible datatypes from the Eigen library:

4.3.2 class Bar

The basic 4-component row spinor.
Inherits
Eigen::RowVector4cd
Methods:
Construct a null Dirac spinor-bar:
Construct a Dirac spinor-bar from four numbers:
Construct a Dirac spinor-bar from other elgible datatypes from the Eigen library:
Assign a Dirac spinor-bar from other elgible datatypes from the Eigen library:

4.3.3 class Any

Represents a generic Dirac Spinor without a specific flavor (u-type or v-type):
Inherits:
Eigen::Vector4cd
Methods:
Construct a generic Dirac spinor from other elgible datatypes from the Eigen library:
Assign a generic Dirac spinor from other elgible datatypes from the Eigen library:


5 Operators

Operators are defined in two header files, Operator4.h and SpecialOperators.h. Operator4.h defines some generic interfaces to four-vector operators like (γ0, γ1, γ2, γ3), and four-tensor operators like σμν ≡ (i)/(2)[γu, γν]. SpecialOperators.h defines objects like the Pauli matrices, the γ matrices, the Feynman slash function, etc., and two classes (Gamma and Sigma4x4) providing implementations of the aforementioned “generic interfaces”.

5.1 class AbsOperator4

#include “Spacetime/Operator4.h”
Description
An AbsOperator4 in this library is a four-vector operator. Examples of this include
AbsOperator4 is an abstract base class for these operators. It specifies two operations, access to elements, and virtual copy constructor (clone).
Methods
Readonly access to elements:
Clone:
Other operations:


5.2 class Operator4

#include “Spacetime/Operator4.h”
Description
Operator4 is an implementation of AbsOperator4, which stores each of the four four-dimensional complex operators in an array of complex matrices. The Weyl representation is used.
Methods
Readonly access to elements:
Read/write access to the elements:
Clone:


5.3 class AbsOperator4x4

#include “Spacetime/Operator4.h”
Description
AbsOperator4x4 is an abstract base class for tensor operators, such as σμν ≡ (i)/(2)[γu, γν], which currently is implemented by class Sigma4x4 which lives in SpecialOperators.h, and is the only subclass of AbsOperator4x4 implemented so far.
Methods
Readonly access to elements:
Clone:
Take a dot product with a four vector to return an Operator 4, eg sigma.dot(q)


5.4 class SU2Generator

#include “Spacetime/SU2Generator.h”
Description
SU2Generator is a class with no instances, only static member functions. It provides the generators Jx, Jy, and Jz in D = 2j + 1 dimensions. It also includes facilities to obtain the rotation matrix, or Drehung-matrix, or D-matrix, for a rotation about any specified axis, and to go backwards from a D-matrix to its logarithm.
static member functions
Access to the generators, held in cache:
Exponentiates linear combinations of SU2 generators in an arbitrary dimensional space to form a rotation (Drehung) matrix, which is a rotation around the axis nHat by the angle phi.
Takes the logarithm of a rotation (Drehung) matrix


5.5 Special Operators

#include “Spacetime/SpecialOperators.h”
Description
The header file SpecialOperators.h defines several objects like Pauli spin matrices, gamma matrices, and also a few objects that give a friendlier interface to those; plus functions for the Feynman slash notation. These are scoped within the namespace SpecialOperators.
const objects defined with the namespace SpecialOperators
Pauli spin matrices:
Gamma matrices (Weyl representation is used):
Array of first four gamma matrices. See also class Gamma:
Projection operators to left-and right-handed chirality states:
4x4 array of sigma matrices and spin matricies (S = σ ⁄ 2). See also class Sigma4x4.
Other operations
Feynman slash notation:


5.5.1 class Gamma

#include “SpecialOperators.h”
description
Access to gamma matrices γ0, γ1, γ2 and γ3. The clone method allows the construction of four-currents.
inherits
AbsOperator4
methods
Clone:
Return one of the gamma matrices


5.5.2 class Sigma4x4

#include “SpecialOperators.h”
description
Provides access to the matrices σμν ≡ (1)/(2)[γμ, γν]. The clone method allows the construction of four-currents.
inherits
Operator4x4
methods
Clone:
Return one of the gamma matrices:
Contract with a four-vector to obtain a four-vector operator:


5.6 class Current4

#include “Spacetime/Current4.h”
description
The Current4 class facilitates the construction of expressions like u(k)γμu(k), u(k)σμνγνu(k), etc. In many cases it is necessary to explicitly create any instances of Current4 because temporary Current4 objects will be created automatically when expresssions such as the above examples are written naturally. The header file must nonetheless be included. This class has no methods, and two friends:
operations:


6 A typical expression

With few classes, you can do a lot of fundamental physics. In this section, we give a brief example. In Compton scattering, one has to evaluate the following matrix elements:
 − is = u ϵν(ieγν)(i( +  + me))/((p + k)2 − m2e)(ieγμ)ϵμu
(s-channel matrix element), and
 − it = u ϵμ(ieγμ)(i ( −  + me))/((p − k)2 − m2e)(ieγν)ϵνu, 
(t-channel matrix element); these can be expressed in the Spacetime library as

DiracSpinor uIn, uOut;
FourVector p, k, kPrime;
ComplexFourVector epsIn, epsOut;
.
.
.
std::complex<double> MS =                                        
  e*e*
  ((bar(uOut)*                                       
  (slash(epsOut.conjugate())*                        
  (slash(p+k)+M)/((p+k).squaredNorm()-m*m)*            
  slash(epsIn))*                                      
  uIn)                                                 
  (0)); 
​
std::complex<double>  MT =                                          
  e*e*
  ((bar(uOut)*                                      
  (slash(epsIn)*                                       
  (slash(p-kPrime)+M)/((p-kPrime).squaredNorm()-m*m)*  
  slash(epsOut.conjugate()))*                        
  uIn)                                                
(0));

after proper initialization of all DiracSpinors, FourVectors, and ComplexFourVectors entering into the expression. Notice that Current4 class temporaries are created so the code will not compile unless the Current4 header file is included.