diff --git a/Helicity/LorentzRSSpinor.h b/Helicity/LorentzRSSpinor.h --- a/Helicity/LorentzRSSpinor.h +++ b/Helicity/LorentzRSSpinor.h @@ -1,571 +1,571 @@ // -*- C++ -*- // // LorentzRSSpinor.h is a part of ThePEG - Toolkit for HEP Event Generation // Copyright (C) 2003-2017 Peter Richardson, Leif Lonnblad // // ThePEG is licenced under version 3 of the GPL, see COPYING for details. // Please respect the MCnet academic guidelines, see GUIDELINES for details. // #ifndef ThePEG_LorentzRSSpinor_H #define ThePEG_LorentzRSSpinor_H #include "ThePEG/Config/ThePEG.h" #include "ThePEG/Vectors/ThreeVector.h" #include "HelicityDefinitions.h" #include "LorentzRSSpinor.fh" #include "LorentzRSSpinorBar.h" #include "LorentzSpinorBar.h" #include "LorentzSpinor.h" #include "LorentzPolarizationVector.h" namespace ThePEG{ namespace Helicity{ /** * The LorentzRSSpinor class is designed to store a Rarita-Schwinger * spinor for a spin-3/2 particle. In addition to storing the * components of the spinor information is stored on the type of * spinor, for example u or v type. * * At the moment only one choice of the Dirac matrix representation * is supported. For high-energy calculations the choice made by the * HELAS collaboration is more efficient for numerical * calculations. In this representation * * \f[ * \gamma_{i=1,2,3}=\left(\begin{array}{cc} * 0 & \sigma_i \\ * -\sigma_i & 0 * \end{array}\right) * \quad * \gamma_0=\left(\begin{array}{cc} * 0 & 1 \\ * 1 & 0 * \end{array}\right) * \quad * \gamma_5=\left(\begin{array}{cc} * -1 & 0 \\ * 0 & 1 * \end{array}\right) * \f] * * The type of the spinor is also stored using the SpinorType * enumeration. There are three types supported SpinorType::u, * SpinorType::v, SpinorType::unknown. This information is intended * mainly for use in the case of Majorana particles where matrix * elements can be calculated with either u or v type spinors and * knowledge of which was used will be needed in order to give the * correct correlations. The SpinorType::unknowne is intended for * cases where either the spinor for an off-shell line in a matrix * element calculation or the information is genuinely unknown. * * The LorentzRSSpinorBar class is also provided to store the barred * spinor. * * * @see HelicityDefinitions * @see LorentzRSSpinorBar * * \author Peter Richardson * */ template class LorentzRSSpinor { public: /** @name Standard constructors. */ //@{ /** * Default zero constructor, optionally specifying \a t, the type. */ LorentzRSSpinor(SpinorType t = SpinorType::unknown) : _type(t), _spin() {} /** * Constructor with complex numbers specifying the components, * optionally specifying \a t, the type. */ LorentzRSSpinor(complex a1, complex b1, complex c1, complex d1, complex a2, complex b2, complex c2, complex d2, complex a3, complex b3, complex c3, complex d3, complex a4, complex b4, complex c4, complex d4, SpinorType t=SpinorType::unknown) : _type(t), _spin{{ {{a1,b1,c1,d1}}, {{a2,b2,c2,d2}}, {{a3,b3,c3,d3}}, {{a4,b4,c4,d4}} }} {} template LorentzRSSpinor(const LorentzRSSpinor & other) : _type(other._type), _spin(other._spin) {} //@} /** @name Access the components. */ //@{ /** * Subscript operator to return spinor components */ complex operator()(int i, int j) const { assert( i >= 0 && i <= 3 && j>=0 && j<=3); return _spin[i][j]; } /** * Set components by index */ complex & operator () (int i, int j) { assert( i >= 0 && i <= 3 && j>=0 && j<=3); return _spin[i][j]; } /** * Get first spinor component for the x vector */ complex xs1() const {return _spin[0][0];} /** * Get second spinor component for the x vector */ complex xs2() const {return _spin[0][1];} /** * Get third spinor component for the x vector */ complex xs3() const {return _spin[0][2];} /** * Get fourth spinor component for the x vector */ complex xs4() const {return _spin[0][3];} /** * Get first spinor component for the y vector */ complex ys1() const {return _spin[1][0];} /** * Get second spinor component for the y vector */ complex ys2() const {return _spin[1][1];} /** * Get third spinor component for the y vector */ complex ys3() const {return _spin[1][2];} /** * Get fourth spinor component for the y vector */ complex ys4() const {return _spin[1][3];} /** * Get first spinor component for the z vector */ complex zs1() const {return _spin[2][0];} /** * Get second spinor component for the z vector */ complex zs2() const {return _spin[2][1];} /** * Get third spinor component for the z vector */ complex zs3() const {return _spin[2][2];} /** * Get fourth spinor component for the z vector */ complex zs4() const {return _spin[2][3];} /** * Get first spinor component for the t vector */ complex ts1() const {return _spin[3][0];} /** * Get second spinor component for the t vector */ complex ts2() const {return _spin[3][1];} /** * Get third spinor component for the t vector */ complex ts3() const {return _spin[3][2];} /** * Get fourth spinor component for the t vector */ complex ts4() const {return _spin[3][3];} /** * Set first spinor component for the x vector */ void setXS1(complex in) {_spin[0][0]=in;} /** * Set second spinor component for the x vector */ void setXS2(complex in) {_spin[0][1]=in;} /** * Set third spinor component for the x vector */ void setXS3(complex in) {_spin[0][2]=in;} /** * Set fourth spinor component for the x vector */ void setXS4(complex in) {_spin[0][3]=in;} /** * Set first spinor component for the y vector */ void setYS1(complex in) {_spin[1][0]=in;} /** * Set second spinor component for the y vector */ void setYS2(complex in) {_spin[1][1]=in;} /** * Set third spinor component for the y vector */ void setYS3(complex in) {_spin[1][2]=in;} /** * Set fourth spinor component for the y vector */ void setYS4(complex in) {_spin[1][3]=in;} /** * Set first spinor component for the z vector */ void setZS1(complex in) {_spin[2][0]=in;} /** * Set second spinor component for the z vector */ void setZS2(complex in) {_spin[2][1]=in;} /** * Set third spinor component for the z vector */ void setZS3(complex in) {_spin[2][2]=in;} /** * Set fourth spinor component for the z vector */ void setZS4(complex in) {_spin[2][3]=in;} /** * Set first spinor component for the t vector */ void setTS1(complex in) {_spin[3][0]=in;} /** * Set second spinor component for the t vector */ void setTS2(complex in) {_spin[3][1]=in;} /** * Set third spinor component for the t vector */ void setTS3(complex in) {_spin[3][2]=in;} /** * Set fourth spinor component for the t vector */ void setTS4(complex in) {_spin[3][3]=in;} //@} /// @name Mathematical assignment operators. //@{ template LorentzRSSpinor & operator+=(const LorentzRSSpinor & a) { for(unsigned int ix=0;ix<4;++ix) for(unsigned int iy=0;iy<4;++iy) _spin[ix][iy] += a._spin[ix][iy]; return *this; } template LorentzRSSpinor & operator-=(const LorentzRSSpinor & a) { for(unsigned int ix=0;ix<4;++ix) for(unsigned int iy=0;iy<4;++iy) _spin[ix][iy] -= a._spin[ix][iy]; return *this; } LorentzRSSpinor & operator*=(double a) { for(unsigned int ix=0;ix<4;++ix) for(unsigned int iy=0;iy<4;++iy) _spin[ix][iy] *=a; return *this; } LorentzRSSpinor & operator/=(double a) { for(unsigned int ix=0;ix<4;++ix) for(unsigned int iy=0;iy<4;++iy) _spin[ix][iy] /=a; return *this; } //@} /** @name Arithmetic operators. */ //@{ /** * dot product with a polarization vector */ LorentzSpinor dot(const LorentzPolarizationVector & vec) const { LorentzSpinor output(_type); complex temp; unsigned int ix; for(ix=0;ix<4;++ix) { temp = _spin[3][ix]*vec.t(); temp -= _spin[0][ix]*vec.x(); temp -= _spin[1][ix]*vec.y(); temp -= _spin[2][ix]*vec.z(); output[ix]=temp; } return output; } /** * dot product with a 4-vector */ LorentzSpinor dot(const LorentzMomentum & invec) const { LorentzSpinor output(_type); complex temp; LorentzVector vec = UnitRemoval::InvE * invec; unsigned int ix; for(ix=0;ix<4;++ix) { temp = - ( _spin[0][ix]*vec.x() + _spin[1][ix]*vec.y()+ _spin[2][ix]*vec.z() ) + _spin[3][ix]*vec.t(); output[ix]=temp; } return output; } //@} /** @name Transformations. */ //@{ /** * return the barred spinor */ LorentzRSSpinorBar bar() const; /** * Standard Lorentz boost specifying the components of the beta vector. */ LorentzRSSpinor & boost(double,double,double); /** * Standard Lorentz boost specifying the beta vector. */ LorentzRSSpinor & boost(const Boost &); /** * General transform */ LorentzRSSpinor & transform(const LorentzRotation &); //@} /** @name Functions related to type. */ //@{ /** * Return the type of the spinor. */ SpinorType Type() const {return _type;} //@} /** * Scalar product \f$\bar{f}^\alpha(c_LP_L+c_RP_R)f_\alpha\f$for general couplings * @param fbar The barred spinor * @param left The left-handed coupling, \f$c_L\f$. * @param right The right-handed coupling, \f$c_R\f$. */ template auto generalScalar(LorentzRSSpinorBar& fbar, Complex left, Complex right) - -> decltype(left*fbar(3,0)*ts1()) + -> decltype(left*fbar(3,0)*this->ts1()) { decltype(left*fbar(3,0)*ts1()) output; unsigned int iz; output = left*(fbar(3,0)*_spin[3][0]+fbar(3,1)*_spin[3][1]) +right*(fbar(3,2)*_spin[3][2]+fbar(3,3)*_spin[3][3]); for(iz=0;iz<3;++iz) { output -= left*(fbar(iz,0)*_spin[iz][0]+fbar(iz,1)*_spin[iz][1]) +right*(fbar(iz,2)*_spin[iz][2]+fbar(iz,3)*_spin[iz][3]); } return output; } /** * Current \f$\bar{f}(c_LP_L+c_RP_R)f^\alpha\f$ for general couplings. * @param fbar The barred spinor * @param left The left-handed coupling, \f$c_L\f$. * @param right The right-handed coupling, \f$c_R\f$. */ template auto generalCurrent(LorentzSpinorBar& fbar, Complex left, Complex right) - -> LorentzVector + -> LorentzVectorts1())> { typedef decltype(left*fbar.s1()*ts1()) ResultT; ResultT output[4]; for(size_t iz=0;iz<4;++iz) output[iz]= left*(fbar.s1()*_spin[iz][0]+fbar.s2()*_spin[iz][1]) +right*(fbar.s3()*_spin[iz][2]+fbar.s4()*_spin[iz][3]); return LorentzVector(output[0],output[1],output[2],output[3]); } private: /** * Type of spinor */ SpinorType _type; /** * Storage of the components. */ std::array,4>,4> _spin; }; /// @name Basic mathematical operations //@{ template inline LorentzRSSpinor operator/(const LorentzRSSpinor & v, Value a) { return LorentzRSSpinor(v.xs1()/a, v.xs2()/a, v.xs3()/a, v.xs4()/a, v.ys1()/a, v.ys2()/a, v.ys3()/a, v.ys4()/a, v.zs1()/a, v.zs2()/a, v.zs3()/a, v.zs4()/a, v.ts1()/a, v.ts2()/a, v.ts3()/a, v.ts4()/a, v.Type()); } inline LorentzRSSpinor operator/(const LorentzRSSpinor & v, Complex a) { return LorentzRSSpinor(v.xs1()/a, v.xs2()/a, v.xs3()/a, v.xs4()/a, v.ys1()/a, v.ys2()/a, v.ys3()/a, v.ys4()/a, v.zs1()/a, v.zs2()/a, v.zs3()/a, v.zs4()/a, v.ts1()/a, v.ts2()/a, v.ts3()/a, v.ts4()/a, v.Type()); } template inline LorentzRSSpinor operator-(const LorentzRSSpinor & v) { return LorentzRSSpinor(-v.xs1(),-v.xs2(),-v.xs3(),-v.xs4(), -v.ys1(),-v.ys2(),-v.ys3(),-v.ys4(), -v.zs1(),-v.zs2(),-v.zs3(),-v.zs4(), -v.ts1(),-v.ts2(),-v.ts3(),-v.ts4(), v.Type()); } template inline LorentzRSSpinor operator+(LorentzRSSpinor a, const LorentzRSSpinor & b) { return a += b; } template inline LorentzRSSpinor operator-(LorentzRSSpinor a, const LorentzRSSpinor & b) { return a -= b; } template inline LorentzRSSpinor operator*(const LorentzRSSpinor & a, double b) { return LorentzRSSpinor(a.xs1()*b, a.xs2()*b, a.xs3()*b, a.xs4()*b, a.ys1()*b, a.ys2()*b, a.ys3()*b, a.ys4()*b, a.zs1()*b, a.zs2()*b, a.zs3()*b, a.zs4()*b, a.ts1()*b, a.ts2()*b, a.ts3()*b, a.ts4()*b,a.Type()); } template inline LorentzRSSpinor operator*(double b, LorentzRSSpinor a) { return a *= b; } template inline LorentzRSSpinor operator*(const LorentzRSSpinor & a, Complex b) { return LorentzRSSpinor(a.xs1()*b, a.xs2()*b, a.xs3()*b, a.xs4()*b, a.ys1()*b, a.ys2()*b, a.ys3()*b, a.ys4()*b, a.zs1()*b, a.zs2()*b, a.zs3()*b, a.zs4()*b, a.ts1()*b, a.ts2()*b, a.ts3()*b, a.ts4()*b,a.Type()); } template inline auto operator*(complex a, const LorentzRSSpinor & v) -> LorentzRSSpinor { return {a*v.xs1(), a*v.xs2(), a*v.xs3(), a*v.xs4(), a*v.ys1(), a*v.ys2(), a*v.ys3(), a*v.ys4(), a*v.zs1(), a*v.zs2(), a*v.zs3(), a*v.zs4(), a*v.ts1(), a*v.ts2(), a*v.ts3(), a*v.ts4(),v.Type()}; } template inline auto operator*(const LorentzRSSpinor & v, complex b) -> LorentzRSSpinor { return b*v; } template inline auto operator/(const LorentzRSSpinor & v, complex b) -> LorentzRSSpinor { return {v.xs1()/b, v.xs2()/b, v.xs3()/b, v.xs4()/b, v.ys1()/b, v.ys2()/b, v.ys3()/b, v.ys4()/b, v.zs1()/b, v.zs2()/b, v.zs3()/b, v.zs4()/b, v.ts1()/b, v.ts2()/b, v.ts3()/b, v.ts4()/b,v.Type()}; } template inline auto operator*(ValueB a, const LorentzRSSpinor & v) -> LorentzRSSpinor { return {a*v.xs1(), a*v.xs2(), a*v.xs3(), a*v.xs4(), a*v.ys1(), a*v.ys2(), a*v.ys3(), a*v.ys4(), a*v.zs1(), a*v.zs2(), a*v.zs3(), a*v.zs4(), a*v.ts1(), a*v.ts2(), a*v.ts3(), a*v.ts4(),v.Type()}; } template inline auto operator*(const LorentzRSSpinor & v, ValueB b) -> LorentzRSSpinor { return b*v; } template inline auto operator/(const LorentzRSSpinor & v, ValueB b) -> LorentzRSSpinor { return {v.xs1()/b, v.xs2()/b, v.xs3()/b, v.xs4()/b, v.ys1()/b, v.ys2()/b, v.ys3()/b, v.ys4()/b, v.zs1()/b, v.zs2()/b, v.zs3()/b, v.zs4()/b, v.ts1()/b, v.ts2()/b, v.ts3()/b, v.ts4()/b,v.Type()}; } //@} } } #ifndef ThePEG_TEMPLATES_IN_CC_FILE #include "LorentzRSSpinor.tcc" #endif #endif diff --git a/Helicity/LorentzRSSpinorBar.h b/Helicity/LorentzRSSpinorBar.h --- a/Helicity/LorentzRSSpinorBar.h +++ b/Helicity/LorentzRSSpinorBar.h @@ -1,506 +1,506 @@ // -*- C++ -*- // // LorentzRSSpinorBar.h is a part of ThePEG - Toolkit for HEP Event Generation // Copyright (C) 2003-2017 Peter Richardson, Leif Lonnblad // // ThePEG is licenced under version 3 of the GPL, see COPYING for details. // Please respect the MCnet academic guidelines, see GUIDELINES for details. // #ifndef ThePEG_LorentzRSSpinorBar_H #define ThePEG_LorentzRSSpinorBar_H // This is the declaration of the LorentzRSSpinorBar class. #include "ThePEG/Config/ThePEG.h" #include "ThePEG/Vectors/ThreeVector.h" #include "HelicityDefinitions.h" #include "LorentzRSSpinor.fh" #include "LorentzRSSpinorBar.fh" #include "LorentzSpinorBar.h" #include "LorentzSpinor.h" #include "LorentzPolarizationVector.h" namespace ThePEG { namespace Helicity { /** * The LorentzRSSpinorBar class implements the storage of a * barred Lorentz Rarita-Schwinger Spinor for a spin-3/2 particle. * The design is based on that of the * LorentzRSSpinor class and the details of the implemented * are discussed in more detail in the header file for that class. * * @see LorentzRSSpinor * * \author Peter Richardson * */ template class LorentzRSSpinorBar { public: /** @name Standard constructors. */ //@{ /** * Default zero constructor, optionally specifying \a t, the type. */ LorentzRSSpinorBar(SpinorType t = SpinorType::unknown) : _type(t), _spin() {} /** * Constructor with complex numbers specifying the components, * optionally specifying \a t, the type. */ LorentzRSSpinorBar(complex a1,complex b1, complex c1,complex d1, complex a2,complex b2, complex c2,complex d2, complex a3,complex b3, complex c3,complex d3, complex a4,complex b4, complex c4,complex d4, SpinorType t=SpinorType::unknown) : _type(t), _spin{{ {{a1,b1,c1,d1}}, {{a2,b2,c2,d2}}, {{a3,b3,c3,d3}}, {{a4,b4,c4,d4}} }} {} template LorentzRSSpinorBar(const LorentzRSSpinorBar & other) : _type(other._type), _spin(other._spin) {} //@} /** @name Access the components. */ //@{ /** * Subscript operator to return spinor components */ complex operator()(int i, int j) const { assert( i >= 0 && i <= 3 && j>=0 && j<=3 ); return _spin[i][j]; } /** * Set components by index */ complex & operator () (int i, int j) { assert( i >= 0 && i <= 3 && j>=0 && j<=3 ); return _spin[i][j]; } /** * Get first spinor component for the x vector */ complex xs1() const {return _spin[0][0];} /** * Get second spinor component for the x vector */ complex xs2() const {return _spin[0][1];} /** * Get third spinor component for the x vector */ complex xs3() const {return _spin[0][2];} /** * Get fourth spinor component for the x vector */ complex xs4() const {return _spin[0][3];} /** * Get first spinor component for the y vector */ complex ys1() const {return _spin[1][0];} /** * Get second spinor component for the y vector */ complex ys2() const {return _spin[1][1];} /** * Get third spinor component for the y vector */ complex ys3() const {return _spin[1][2];} /** * Get fourth spinor component for the y vector */ complex ys4() const {return _spin[1][3];} /** * Get first spinor component for the z vector */ complex zs1() const {return _spin[2][0];} /** * Get second spinor component for the z vector */ complex zs2() const {return _spin[2][1];} /** * Get third spinor component for the z vector */ complex zs3() const {return _spin[2][2];} /** * Get fourth spinor component for the z vector */ complex zs4() const {return _spin[2][3];} /** * Get first spinor component for the t vector */ complex ts1() const {return _spin[3][0];} /** * Get second spinor component for the t vector */ complex ts2() const {return _spin[3][1];} /** * Get third spinor component for the t vector */ complex ts3() const {return _spin[3][2];} /** * Get fourth spinor component for the t vector */ complex ts4() const {return _spin[3][3];} /** * Set first spinor component for the x vector */ void setXS1(complex in) {_spin[0][0]=in;} /** * Set second spinor component for the x vector */ void setXS2(complex in) {_spin[0][1]=in;} /** * Set third spinor component for the x vector */ void setXS3(complex in) {_spin[0][2]=in;} /** * Set fourth spinor component for the x vector */ void setXS4(complex in) {_spin[0][3]=in;} /** * Set first spinor component for the y vector */ void setYS1(complex in) {_spin[1][0]=in;} /** * Set second spinor component for the y vector */ void setYS2(complex in) {_spin[1][1]=in;} /** * Set third spinor component for the y vector */ void setYS3(complex in) {_spin[1][2]=in;} /** * Set fourth spinor component for the y vector */ void setYS4(complex in) {_spin[1][3]=in;} /** * Set first spinor component for the z vector */ void setZS1(complex in) {_spin[2][0]=in;} /** * Set second spinor component for the z vector */ void setZS2(complex in) {_spin[2][1]=in;} /** * Set third spinor component for the z vector */ void setZS3(complex in) {_spin[2][2]=in;} /** * Set fourth spinor component for the z vector */ void setZS4(complex in) {_spin[2][3]=in;} /** * Set first spinor component for the t vector */ void setTS1(complex in) {_spin[3][0]=in;} /** * Set second spinor component for the t vector */ void setTS2(complex in) {_spin[3][1]=in;} /** * Set third spinor component for the t vector */ void setTS3(complex in) {_spin[3][2]=in;} /** * Set fourth spinor component for the t vector */ void setTS4(complex in ) {_spin[3][3]=in;} //@} /// @name Mathematical assignment operators. //@{ template LorentzRSSpinorBar & operator+=(const LorentzRSSpinorBar & a) { for(unsigned int ix=0;ix<4;++ix) for(unsigned int iy=0;iy<4;++iy) _spin[ix][iy] += a._spin[ix][iy]; return *this; } template LorentzRSSpinorBar & operator-=(const LorentzRSSpinorBar & a) { for(unsigned int ix=0;ix<4;++ix) for(unsigned int iy=0;iy<4;++iy) _spin[ix][iy] -= a._spin[ix][iy]; return *this; } LorentzRSSpinorBar & operator*=(double a) { for(unsigned int ix=0;ix<4;++ix) for(unsigned int iy=0;iy<4;++iy) _spin[ix][iy] *=a; return *this; } LorentzRSSpinorBar & operator/=(double a) { for(unsigned int ix=0;ix<4;++ix) for(unsigned int iy=0;iy<4;++iy) _spin[ix][iy] /=a; return *this; } //@} /** @name Arithmetic operators. */ //@{ /** * dot product with a polarization vector */ LorentzSpinorBar dot(const LorentzPolarizationVector & vec) const { LorentzSpinorBar output(_type); for(unsigned int ix=0;ix<4;++ix) { output[ix]=_spin[3][ix]*vec.t()-_spin[0][ix]*vec.x() -_spin[1][ix]*vec.y()-_spin[2][ix]*vec.z(); } return output; } /** * dot product with a 4-momentum */ LorentzSpinorBar dot(const LorentzMomentum & invec) const { LorentzSpinorBar output(_type); LorentzVector vec = UnitRemoval::InvE * invec; unsigned int ix; for(ix=0;ix<4;++ix) { output[ix]=_spin[3][ix]*vec.t()-_spin[0][ix]*vec.x() -_spin[1][ix]*vec.y()-_spin[2][ix]*vec.z(); } return output; } //@} /** @name Transformations. */ //@{ /** * return the barred spinor */ LorentzRSSpinor bar() const; /** * Standard Lorentz boost specifying the components of the beta vector. */ LorentzRSSpinorBar & boost(double,double,double); /** * Standard Lorentz boost specifying the beta vector. */ LorentzRSSpinorBar & boost(const Boost &); /** * General transform */ LorentzRSSpinorBar & transform(const LorentzRotation &); //@} /** @name Functions related to type. */ //@{ /** * Return the type of the spinor. */ SpinorType Type() const {return _type;} //@} /** * Current \f$\bar{f}^\alpha(c_LP_L+c_RP_R)f\f$ for general couplings. * @param f The unbarred spinor * @param left The left-handed coupling, \f$c_L\f$. * @param right The right-handed coupling, \f$c_R\f$. */ template auto generalCurrent(LorentzSpinor& f, Complex left, Complex right) - -> LorentzVector + -> LorentzVectorts1()*f.s1())> { typedef decltype(left*ts1()*f.s1()) ResultT; ResultT output[4]; unsigned int iz; for(iz=0;iz<4;++iz){ output[iz]= left*(_spin[iz][0]*f.s1()+_spin[iz][1]*f.s2()) +right*(_spin[iz][2]*f.s3()+_spin[iz][3]*f.s4()); } return LorentzVector(output[0],output[1], output[2],output[3]); } private: /** * Type of spinor. */ SpinorType _type; /** * Storage of the components. */ std::array,4>,4> _spin; }; /// @name Basic mathematical operations //@{ template inline LorentzRSSpinorBar operator/(const LorentzRSSpinorBar & v, Value a) { return LorentzRSSpinorBar(v.xs1()/a, v.xs2()/a, v.xs3()/a, v.xs4()/a, v.ys1()/a, v.ys2()/a, v.ys3()/a, v.ys4()/a, v.zs1()/a, v.zs2()/a, v.zs3()/a, v.zs4()/a, v.ts1()/a, v.ts2()/a, v.ts3()/a, v.ts4()/a, v.Type()); } inline LorentzRSSpinorBar operator/(const LorentzRSSpinorBar & v, Complex a) { return LorentzRSSpinorBar(v.xs1()/a, v.xs2()/a, v.xs3()/a, v.xs4()/a, v.ys1()/a, v.ys2()/a, v.ys3()/a, v.ys4()/a, v.zs1()/a, v.zs2()/a, v.zs3()/a, v.zs4()/a, v.ts1()/a, v.ts2()/a, v.ts3()/a, v.ts4()/a, v.Type()); } template inline LorentzRSSpinorBar operator-(const LorentzRSSpinorBar & v) { return LorentzRSSpinorBar(-v.xs1(),-v.xs2(),-v.xs3(),-v.xs4(), -v.ys1(),-v.ys2(),-v.ys3(),-v.ys4(), -v.zs1(),-v.zs2(),-v.zs3(),-v.zs4(), -v.ts1(),-v.ts2(),-v.ts3(),-v.ts4(), v.Type()); } template inline LorentzRSSpinorBar operator+(LorentzRSSpinorBar a, const LorentzRSSpinorBar & b) { return a += b; } template inline LorentzRSSpinorBar operator-(LorentzRSSpinorBar a, const LorentzRSSpinorBar & b) { return a -= b; } template inline LorentzRSSpinorBar operator*(const LorentzRSSpinorBar & a, double b) { return LorentzRSSpinorBar(a.xs1()*b, a.xs2()*b, a.xs3()*b, a.xs4()*b, a.ys1()*b, a.ys2()*b, a.ys3()*b, a.ys4()*b, a.zs1()*b, a.zs2()*b, a.zs3()*b, a.zs4()*b, a.ts1()*b, a.ts2()*b, a.ts3()*b, a.ts4()*b,a.Type()); } template inline LorentzRSSpinorBar operator*(double b, LorentzRSSpinorBar a) { return a *= b; } template inline LorentzRSSpinorBar operator*(const LorentzRSSpinorBar & a, Complex b) { return LorentzRSSpinorBar(a.xs1()*b, a.xs2()*b, a.xs3()*b, a.xs4()*b, a.ys1()*b, a.ys2()*b, a.ys3()*b, a.ys4()*b, a.zs1()*b, a.zs2()*b, a.zs3()*b, a.zs4()*b, a.ts1()*b, a.ts2()*b, a.ts3()*b, a.ts4()*b,a.Type()); } template inline auto operator*(complex a, const LorentzRSSpinorBar & v) -> LorentzRSSpinorBar { return {a*v.xs1(), a*v.xs2(), a*v.xs3(), a*v.xs4(), a*v.ys1(), a*v.ys2(), a*v.ys3(), a*v.ys4(), a*v.zs1(), a*v.zs2(), a*v.zs3(), a*v.zs4(), a*v.ts1(), a*v.ts2(), a*v.ts3(), a*v.ts4(),v.Type()}; } template inline auto operator*(const LorentzRSSpinorBar & v, complex b) -> LorentzRSSpinorBar { return b*v; } template inline auto operator/(const LorentzRSSpinorBar & v, complex b) -> LorentzRSSpinorBar { return {v.xs1()/b, v.xs2()/b, v.xs3()/b, v.xs4()/b, v.ys1()/b, v.ys2()/b, v.ys3()/b, v.ys4()/b, v.zs1()/b, v.zs2()/b, v.zs3()/b, v.zs4()/b, v.ts1()/b, v.ts2()/b, v.ts3()/b, v.ts4()/b,v.Type()}; } template inline auto operator*(ValueB a, const LorentzRSSpinorBar & v) -> LorentzRSSpinorBar { return {a*v.xs1(), a*v.xs2(), a*v.xs3(), a*v.xs4(), a*v.ys1(), a*v.ys2(), a*v.ys3(), a*v.ys4(), a*v.zs1(), a*v.zs2(), a*v.zs3(), a*v.zs4(), a*v.ts1(), a*v.ts2(), a*v.ts3(), a*v.ts4(),v.Type()}; } template inline auto operator*(const LorentzRSSpinorBar & v, ValueB b) -> LorentzRSSpinorBar { return b*v; } template inline auto operator/(const LorentzRSSpinorBar & v, ValueB b) -> LorentzRSSpinorBar { return {v.xs1()/b, v.xs2()/b, v.xs3()/b, v.xs4()/b, v.ys1()/b, v.ys2()/b, v.ys3()/b, v.ys4()/b, v.zs1()/b, v.zs2()/b, v.zs3()/b, v.zs4()/b, v.ts1()/b, v.ts2()/b, v.ts3()/b, v.ts4()/b,v.Type()}; } //@} } } #ifndef ThePEG_TEMPLATES_IN_CC_FILE #include "LorentzRSSpinorBar.tcc" #endif #endif diff --git a/Helicity/LorentzSpinor.h b/Helicity/LorentzSpinor.h --- a/Helicity/LorentzSpinor.h +++ b/Helicity/LorentzSpinor.h @@ -1,626 +1,626 @@ // -*- C++ -*- // // LorentzSpinor.h is a part of ThePEG - Toolkit for HEP Event Generation // Copyright (C) 2003-2017 Peter Richardson, Leif Lonnblad // // ThePEG is licenced under version 3 of the GPL, see COPYING for details. // Please respect the MCnet academic guidelines, see GUIDELINES for details. // #ifndef ThePEG_LorentzSpinor_H #define ThePEG_LorentzSpinor_H // This is the declaration of the LorentzSpinor class. #include "ThePEG/Config/ThePEG.h" #include "ThePEG/Vectors/LorentzRotation.h" #include "ThePEG/Vectors/ThreeVector.h" #include "HelicityDefinitions.h" #include "LorentzSpinor.fh" #include "LorentzSpinorBar.h" #include "LorentzPolarizationVector.h" #include "LorentzTensor.h" #include namespace ThePEG{ namespace Helicity{ /** * The LorentzSpinor class is designed to store a spinor. In addition * to storing the components of the spinor, information is stored on * the representation of the type of spinor, for example u or v type. * * At the moment only one choice of the Dirac matrix representation * is supported. For high-energy calculations the choice made by the * HELAS collaboration is more efficient for numerical * calculations. In this representation * * \f[ * \gamma_{i=1,2,3}=\left(\begin{array}{cc} * 0 & \sigma_i \\ * -\sigma_i & 0 * \end{array}\right) * \quad * \gamma_0=\left(\begin{array}{cc} * 0 & 1 \\ * 1 & 0 * \end{array}\right) * \quad * \gamma_5=\left(\begin{array}{cc} * -1 & 0 \\ * 0 & 1 * \end{array}\right) * \f] * * The type of the spinor is also stored using the SpinorType * enumeration. There are three types supported SpinorType::u, * SpinorType::v, SpinorType::unknown. This information is intended * mainly for use in the case of Majorana particles where matrix * elements can be calculated with either u or v type spinors and * knowledge of which was used will be needed in order to give the * correct correlations. The SpinorType::unknowne is intended for * cases where either the spinor for an off-shell line in a matrix * element calculation or the information is genuinely unknown. * * The LorentzSpinorBar class is also provided to store the barred * spinor. * * @see LorentzSpinorBar * * @author Peter Richardson * */ template class LorentzSpinor { public: /** @name Standard constructors. */ //@{ /** * Default zero constructor, optionally specifying \a t, the type. */ LorentzSpinor(SpinorType t = SpinorType::unknown) : _type(t), _spin() {} /** * Constructor with complex numbers specifying the components, * optionally specifying \a s, the type. */ LorentzSpinor(complex a,complex b, complex c,complex d, SpinorType s = SpinorType::unknown) : _type(s), _spin{{a,b,c,d}} {} //@} template LorentzSpinor(const LorentzSpinor & other) : _type(other._type), _spin(other._spin) {} /** @name Access the components. */ //@{ /** * Subscript operator to return spinor components */ complex operator[](int i) const { assert( i >= 0 && i <= 3 ); return _spin[i]; } /** * Subscript operator to return spinor components */ complex operator()(int i) const { assert( i >= 0 && i <= 3 ); return _spin[i]; } /** * Set components by index. */ complex & operator()(int i) { assert( i >= 0 && i <= 3 ); return _spin[i]; } /** * Set components by index. */ complex & operator[](int i) { assert( i >= 0 && i <= 3 ); return _spin[i]; } /** * Get first component. */ complex s1() const {return _spin[0];} /** * Get second component. */ complex s2() const {return _spin[1];} /** * Get third component. */ complex s3() const {return _spin[2];} /** * Get fourth component. */ complex s4() const {return _spin[3];} /** * Set first component. */ void setS1(complex in) {_spin[0]=in;} /** * Set second component. */ void setS2(complex in) {_spin[1]=in;} /** * Set third component. */ void setS3(complex in) {_spin[2]=in;} /** * Set fourth component. */ void setS4(complex in) {_spin[3]=in;} //@} /// @name Mathematical assignment operators. //@{ template LorentzSpinor & operator+=(const LorentzSpinor & a) { for(unsigned int ix=0;ix<4;++ix) _spin[ix] += a._spin[ix]; return *this; } template LorentzSpinor & operator-=(const LorentzSpinor & a) { for(unsigned int ix=0;ix<4;++ix) _spin[ix] -= a._spin[ix]; return *this; } LorentzSpinor & operator*=(double a) { for(unsigned int ix=0;ix<4;++ix) _spin[ix] *=a; return *this; } LorentzSpinor & operator/=(double a) { for(unsigned int ix=0;ix<4;++ix) _spin[ix] /=a; return *this; } //@} /** @name Transformations. */ //@{ /** * Return the barred spinor */ LorentzSpinorBar bar() const; /** * Return the conjugated spinor \f$u_c=C\bar{u}^T\f$. This operation * transforms u-spinors to v-spinors and vice-versa and is required when * dealing with majorana particles. */ LorentzSpinor conjugate() const; /** * Standard Lorentz boost specifying the components of the beta vector. */ LorentzSpinor & boost(double,double,double); /** * Standard Lorentz boost specifying the beta vector. */ LorentzSpinor & boost(const Boost &); /** * General Lorentz transformation */ LorentzSpinor & transform(const SpinHalfLorentzRotation & ); /** * General Lorentz transformation */ LorentzSpinor & transform(const LorentzRotation & r) { transform(r.half()); return *this; } //@} /** @name Functions related to type. */ //@{ /** * Return the type of the spinor. */ SpinorType Type() const {return _type;} //@} /** * @name Functions to apply the projection operator */ //@{ /** * Apply \f$p\!\!\!\!\!\not\,\,\,+m\f$ */ template auto projectionOperator(const LorentzVector & p, const ValueB & m) const -> LorentzSpinor { LorentzSpinor spin; static const Complex ii(0.,1.); complex p0pp3=p.t()+p.z(); complex p0mp3=p.t()-p.z(); complex p1pp2=p.x()+ii*p.y(); complex p1mp2=p.x()-ii*p.y(); spin.setS1(m*s1()+p0mp3*s3()-p1mp2*s4()); spin.setS2(m*s2()+p0pp3*s4()-p1pp2*s3()); spin.setS3(m*s3()+p0pp3*s1()+p1mp2*s2()); spin.setS4(m*s4()+p0mp3*s2()+p1pp2*s1()); return spin; } /** * Apply \f$g^LP_L+g^RP_R\f$ */ LorentzSpinor helicityProjectionOperator(const Complex & gL, const Complex & gR) const { LorentzSpinor spin; spin.setS1(gL*s1()); spin.setS2(gL*s2()); spin.setS3(gR*s3()); spin.setS4(gR*s4()); return spin; } //@} /** @name Functions to calculate certain currents. */ //@{ /** * Apply \f$p\!\!\!\!\!\not\f$ */ template auto slash(const LorentzVector & p) const -> LorentzSpinor { LorentzSpinor spin; static const Complex ii(0.,1.); complex p0pp3=p.t()+p.z(); complex p0mp3=p.t()-p.z(); complex p1pp2=p.x()+ii*p.y(); complex p1mp2=p.x()-ii*p.y(); spin.setS1(p0mp3*s3()-p1mp2*s4()); spin.setS2(p0pp3*s4()-p1pp2*s3()); spin.setS3(p0pp3*s1()+p1mp2*s2()); spin.setS4(p0mp3*s2()+p1pp2*s1()); return spin; } /** * Apply \f$p\!\!\!\!\!\not\f$ */ template auto slash(const LorentzVector > & p) const -> LorentzSpinor { LorentzSpinor spin; static const Complex ii(0.,1.); complex p0pp3=p.t()+p.z(); complex p0mp3=p.t()-p.z(); complex p1pp2=p.x()+ii*p.y(); complex p1mp2=p.x()-ii*p.y(); spin.setS1(p0mp3*s3()-p1mp2*s4()); spin.setS2(p0pp3*s4()-p1pp2*s3()); spin.setS3(p0pp3*s1()+p1mp2*s2()); spin.setS4(p0mp3*s2()+p1pp2*s1()); return spin; } /** * Calculate the left-handed current \f$\bar{f}\gamma^\mu P_Lf\f$. * @param fb The barred spinor. */ template auto leftCurrent(const LorentzSpinorBar& fb) const - -> LorentzVector + -> LorentzVectors2())> { typedef decltype(fb.s3()*s2()) ResultT; LorentzVector vec; Complex ii(0.,1.); ResultT p1(fb.s3()*s2()),p2(fb.s4()*s1()); vec.setX( -(p1+p2) ); vec.setY( ii*(p1-p2) ); p1 = fb.s3()*s1();p2 = fb.s4()*s2(); vec.setZ( -(p1-p2) ); vec.setT( (p1+p2) ); return vec; } /** * Calculate the right-handed current \f$\bar{f}\gamma^\mu P_Rf\f$. * @param fb The barred spinor. */ template auto rightCurrent(const LorentzSpinorBar& fb) const - -> LorentzVector + -> LorentzVectors4())> { typedef decltype(fb.s1()*s4()) ResultT; LorentzVector vec; Complex ii(0.,1.); ResultT p1(fb.s1()*s4()),p2(fb.s2()*s3()); vec.setX( (p1+p2)); vec.setY( -ii*(p1-p2)); p1 = fb.s1()*s3();p2 = fb.s2()*s4(); vec.setZ( (p1-p2)); vec.setT( (p1+p2)); return vec; } /** * Calculate the vector current \f$\bar{f}\gamma^\mu f\f$ * @param fb The barred spinor. */ template auto vectorCurrent(const LorentzSpinorBar& fb) const - -> LorentzVector + -> LorentzVectors4())> { - typedef decltype(fb.s1()*s4()) ResultT; + typedef decltype(fb.s1()*this->s4()) ResultT; LorentzVector vec; Complex ii(0.,1.); ResultT s1s4(fb.s1()*s4()),s2s3(fb.s2()*s3()), s3s2(fb.s3()*s2()),s4s1(fb.s4()*s1()), s1s3(fb.s1()*s3()),s2s4(fb.s2()*s4()), s3s1(fb.s3()*s1()),s4s2(fb.s4()*s2()); vec.setX( s1s4+s2s3-s3s2-s4s1 ); vec.setY( -ii*(s1s4-s2s3-s3s2+s4s1)); vec.setZ( s1s3-s2s4-s3s1+s4s2 ); vec.setT( s1s3+s2s4+s3s1+s4s2); return vec; } /** * Calculate a general current with arbitary left and right couplings, * i.e. \f$\bar{f}\gamma^\mu(c_lP_L+c_RP_R)f\f$ * @param fb The barred spinor. * @param left The left coupling, \f$c_L\f$. * @param right The right coupling, \f$c_R\f$. */ template auto generalCurrent(const LorentzSpinorBar& fb, Complex left, Complex right) const - -> LorentzVector + -> LorentzVectors2())> { - typedef decltype(fb.s3()*s2()) ResultT; + typedef decltype(fb.s3()*this->s2()) ResultT; LorentzVector vec; Complex ii(0.,1.); ResultT p1(fb.s3()*s2()),p2(fb.s4()*s1()); vec.setX( -left*(p1+p2)); vec.setY( ii*left*(p1-p2)); p1 = fb.s3()*s1();p2 = fb.s4()*s2(); vec.setZ( -left*(p1-p2)); vec.setT( left*(p1+p2)); p1=fb.s1()*s4();p2=fb.s2()*s3(); vec.setX(vec.x()+right*(p1+p2)); vec.setY(vec.y()-ii*right*(p1-p2)); p1 = fb.s1()*s3();p2 = fb.s2()*s4(); vec.setZ(vec.z()+right*(p1-p2)); vec.setT(vec.t()+right*(p1+p2)); return vec; } //@} /** @name Functions to calculate certain scalars. */ //@{ /** * Calculate the left-handed scalar \f$\bar{f}P_Lf\f$. * @param fb The barred spinor. */ template auto leftScalar(const LorentzSpinorBar& fb) const - -> decltype(fb.s1()*s1()) + -> decltype(fb.s1()*this->s1()) { return fb.s1()*s1()+fb.s2()*s2(); } /** * Calculate the right-handed scalar \f$\bar{f}P_Rf\f$. * @param fb The barred spinor. */ template auto rightScalar(const LorentzSpinorBar& fb) const - -> decltype(fb.s3()*s3()) + -> decltype(fb.s3()*this->s3()) { return fb.s3()*s3()+fb.s4()*s4(); } /** * Calculate the scalar \f$\bar{f}f\f$. * @param fb The barred spinor. */ template auto scalar(const LorentzSpinorBar& fb) const - -> decltype(fb.s1()*s1()) + -> decltype(fb.s1()*this->s1()) { return fb.s1()*s1()+fb.s2()*s2() +fb.s3()*s3()+fb.s4()*s4(); } /** * Calculate the pseudoscalar \f$\bar{f}\gamma_5f\f$. * @param fb The barred spinor. */ template auto pseudoScalar(const LorentzSpinorBar& fb) const - -> decltype(fb.s1()*s1()) + -> decltype(fb.s1()*this->s1()) { return -fb.s1()*s1()-fb.s2()*s2() +fb.s3()*s3()+fb.s4()*s4(); } /** * Calculate a general scalar product with arbitary left and right couplings, * i.e. \f$\bar{f}c_lP_L+c_RP_Rf\f$ * @param fb The barred spinor. * @param left The left coupling, \f$c_L\f$. * @param right The right coupling, \f$c_R\f$. */ template auto generalScalar(const LorentzSpinorBar& fb, Complex left, Complex right) const - -> decltype(left*fb.s1()*s1()) + -> decltype(left*fb.s1()*this->s1()) { return left*(fb.s1()*s1()+fb.s2()*s2()) + right*(fb.s3()*s3()+fb.s4()*s4()); } //@} /** * Calculate the product with \f$\sigma^{\mu\nu}\f$, i.e. * \f$\bar{f}\sigma^{\mu\nu}f\f$ */ template auto sigma(const LorentzSpinorBar& fb) const - -> LorentzTensor + -> LorentzTensors1())> { - typedef decltype(fb.s1()*s1()) ResultT; + typedef decltype(fb.s1()*this->s1()) ResultT; LorentzTensor output; ResultT s11(fb.s1()*s1()),s22(fb.s2()*s2()), s33(fb.s3()*s3()),s44(fb.s4()*s4()), s12(fb.s1()*s2()),s21(fb.s2()*s1()), s34(fb.s3()*s4()),s43(fb.s4()*s3()); Complex ii(0.,1.); ResultT zero; zero = ZERO; output.setTT( zero ); output.setTX(-ii*( s12+s21-s34-s43)); output.setTY( -s12+s21+s34-s43 ); output.setTZ(-ii*( s11-s22-s33+s44)); output.setXT( -output.tx() ); output.setXX( zero ); output.setXY( s11-s22+s33-s44 ); output.setXZ(-ii*(-s12+s21-s34+s43)); output.setYT( -output.ty() ); output.setYX( -output.xy() ); output.setYY( zero ); output.setYZ( s12+s21+s34+s43 ); output.setZT( -output.tz() ); output.setZX( -output.xz() ); output.setZY( -output.yz() ); output.setZZ( zero ); return output; } private: /** * Type of spinor */ SpinorType _type; /** * Storage of the components. */ std::array,4> _spin; }; /// @name Basic mathematical operations //@{ template inline LorentzSpinor operator/(const LorentzSpinor & v, Value a) { return LorentzSpinor(v.s1()/a, v.s2()/a, v.s3()/a, v.s4()/a,v.Type()); } inline LorentzSpinor operator/(const LorentzSpinor & v, Complex a) { return LorentzSpinor(v.s1()/a, v.s2()/a, v.s3()/a, v.s4()/a,v.Type()); } template inline LorentzSpinor operator-(const LorentzSpinor & v) { return LorentzSpinor(-v.s1(),-v.s2(),-v.s3(),-v.s4(),v.Type()); } template inline LorentzSpinor operator+(LorentzSpinor a, const LorentzSpinor & b) { return a += b; } template inline LorentzSpinor operator-(LorentzSpinor a, const LorentzSpinor & b) { return a -= b; } template inline LorentzSpinor operator*(const LorentzSpinor & a, double b) { return LorentzSpinor(a.s1()*b, a.s2()*b, a.s3()*b, a.s4()*b,a.Type()); } template inline LorentzSpinor operator*(double b, LorentzSpinor a) { return a *= b; } template inline LorentzSpinor operator*(const LorentzSpinor & a, Complex b) { return LorentzSpinor(a.s1()*b, a.s2()*b, a.s3()*b, a.s4()*b,a.Type()); } template inline auto operator*(complex a, const LorentzSpinor & v) -> LorentzSpinor { return {a*v.s1(), a*v.s2(), a*v.s3(), a*v.s4(),v.Type()}; } template inline auto operator*(const LorentzSpinor & v, complex b) -> LorentzSpinor { return b*v; } template inline auto operator/(const LorentzSpinor & v, complex b) -> LorentzSpinor { return {v.s1()/b, v.s2()/b, v.s3()/b, v.s4()/b,v.Type()}; } template inline auto operator*(ValueB a, const LorentzSpinor & v) -> LorentzSpinor { return {a*v.s1(), a*v.s2(), a*v.s3(), a*v.s4(),v.Type()}; } template inline auto operator*(const LorentzSpinor & v, ValueB b) -> LorentzSpinor { return b*v; } template inline auto operator/(const LorentzSpinor & v, ValueB b) -> LorentzSpinor { return {v.s1()/b, v.s2()/b, v.s3()/b, v.s4()/b,v.Type()}; } } } #ifndef ThePEG_TEMPLATES_IN_CC_FILE #include "LorentzSpinor.tcc" #endif #endif diff --git a/Helicity/LorentzTensor.h b/Helicity/LorentzTensor.h --- a/Helicity/LorentzTensor.h +++ b/Helicity/LorentzTensor.h @@ -1,533 +1,533 @@ // -*- C++ -*- // // LorentzTensor.h is a part of ThePEG - Toolkit for HEP Event Generation // Copyright (C) 2003-2017 Peter Richardson, Leif Lonnblad // // ThePEG is licenced under version 3 of the GPL, see COPYING for details. // Please respect the MCnet academic guidelines, see GUIDELINES for details. // #ifndef ThePEG_LorentzTensor_H #define ThePEG_LorentzTensor_H // This is the declaration of the LorentzTensor class. #include "ThePEG/Config/PhysicalQtyComplex.h" #include "ThePEG/Config/ThePEG.h" #include "LorentzPolarizationVector.h" namespace ThePEG { namespace Helicity { // compiler magic needs these pre-declarations to make friend templates work template class LorentzTensor; /** * The LorentzTensor class is designed to implement the storage of a * complex tensor to be used to representation the wavefunction of a * spin-2 particle. * * At the moment it only implements the storage of the tensor * components but it is envisaged that it will be extended to include * boost methods etc. * * @author Peter Richardson * */ template class LorentzTensor { public: /** @name Standard constructors and destructors. */ //@{ /** * Default zero constructor. */ LorentzTensor() = default; /** * Constructor specifyign all components. */ LorentzTensor(complex xx, complex xy, complex xz, complex xt, complex yx, complex yy, complex yz, complex yt, complex zx, complex zy, complex zz, complex zt, complex tx, complex ty, complex tz, complex tt) : _tensor{{ {{xx,xy,xz,xt}}, {{yx,yy,yz,yt}}, {{zx,zy,zz,zt}}, {{tx,ty,tz,tt}} }} {} /** * Constructor in terms of two polarization vectors. */ LorentzTensor(const LorentzPolarizationVector & p, const LorentzPolarizationVector & q) { setXX(p.x() * q.x()); setYX(p.y() * q.x()); setZX(p.z() * q.x()); setTX(p.t() * q.x()); setXY(p.x() * q.y()); setYY(p.y() * q.y()); setZY(p.z() * q.y()); setTY(p.t() * q.y()); setXZ(p.x() * q.z()); setYZ(p.y() * q.z()); setZZ(p.z() * q.z()); setTZ(p.t() * q.z()); setXT(p.x() * q.t()); setYT(p.y() * q.t()); setZT(p.z() * q.t()); setTT(p.t() * q.t()); } //@} /** @name Access individual components. */ //@{ /** * Get x,x component. */ complex xx() const {return _tensor[0][0];} /** * Get y,x component. */ complex yx() const {return _tensor[1][0];} /** * Get z,x component. */ complex zx() const {return _tensor[2][0];} /** * Get t,x component. */ complex tx() const {return _tensor[3][0];} /** * Get x,y component. */ complex xy() const {return _tensor[0][1];} /** * Get y,y component. */ complex yy() const {return _tensor[1][1];} /** * Get z,y component. */ complex zy() const {return _tensor[2][1];} /** * Get t,y component. */ complex ty() const {return _tensor[3][1];} /** * Get x,z component. */ complex xz() const {return _tensor[0][2];} /** * Get y,z component. */ complex yz() const {return _tensor[1][2];} /** * Get z,z component. */ complex zz() const {return _tensor[2][2];} /** * Get t,z component. */ complex tz() const {return _tensor[3][2];} /** * Get x,t component. */ complex xt() const {return _tensor[0][3];} /** * Get y,t component. */ complex yt() const {return _tensor[1][3];} /** * Get z,t component. */ complex zt() const {return _tensor[2][3];} /** * Get t,t component. */ complex tt() const {return _tensor[3][3];} /** * Set x,x component. */ void setXX(complex a) {_tensor[0][0]=a;} /** * Set y,x component. */ void setYX(complex a) {_tensor[1][0]=a;} /** * Set z,x component. */ void setZX(complex a) {_tensor[2][0]=a;} /** * Set t,x component. */ void setTX(complex a) {_tensor[3][0]=a;} /** * Set x,y component. */ void setXY(complex a) {_tensor[0][1]=a;} /** * Set y,y component. */ void setYY(complex a) {_tensor[1][1]=a;} /** * Set z,y component. */ void setZY(complex a) {_tensor[2][1]=a;} /** * Set t,y component. */ void setTY(complex a) {_tensor[3][1]=a;} /** * Set x,z component. */ void setXZ(complex a) {_tensor[0][2]=a;} /** * Set y,z component. */ void setYZ(complex a) {_tensor[1][2]=a;} /** * Set z,z component. */ void setZZ(complex a) {_tensor[2][2]=a;} /** * Set t,z component. */ void setTZ(complex a) {_tensor[3][2]=a;} /** * Set x,t component. */ void setXT(complex a) {_tensor[0][3]=a;} /** * Set y,t component. */ void setYT(complex a) {_tensor[1][3]=a;} /** * Set z,t component. */ void setZT(complex a) {_tensor[2][3]=a;} /** * Set t,t component. */ void setTT(complex a) {_tensor[3][3]=a;} /** * Get components by indices. */ complex operator () (int i, int j) const { assert( i>=0 && i<=3 && j>=0 && j<=3); return _tensor[i][j]; } /** * Set components by indices. */ complex & operator () (int i, int j) { assert( i>=0 && i<=3 && j>=0 && j<=3); return _tensor[i][j]; } //@} /** @name Transformations. */ //@{ /** * Standard Lorentz boost specifying the components of the beta vector. */ LorentzTensor & boost(double,double,double); /** * Standard Lorentz boost specifying the beta vector. */ LorentzTensor & boost(const Boost & b) { return boost(b.x(), b.y(), b.z()); } /** * General Lorentz transformation */ LorentzTensor & transform(const SpinOneLorentzRotation & r){ unsigned int ix,iy,ixa,iya; LorentzTensor output; complex temp; for(ix=0;ix<4;++ix) { for(iy=0;iy<4;++iy) { temp=complex(); for(ixa=0;ixa<4;++ixa) { for(iya=0;iya<4;++iya) temp+=r(ix,ixa)*r(iy,iya)*(*this)(ixa,iya); } output(ix,iy)=temp; } } *this=output; return *this; } /** * Return the complex conjugate. */ LorentzTensor conjugate() { return LorentzTensor(conj(xx()), conj(xy()), conj(xz()), conj(xt()), conj(yx()), conj(yy()), conj(yz()), conj(yt()), conj(zx()), conj(zy()), conj(zz()), conj(zt()), conj(tx()), conj(ty()), conj(tz()), conj(tt())); } //@} /** @name Arithmetic operators. */ //@{ /** * Scaling with a complex number */ LorentzTensor operator*=(Complex a) { for(int ix=0;ix<4;++ix) for(int iy=0;iy<4;++iy) _tensor[ix][iy]*=a; return *this; } /** * Scalar product with other tensor */ template friend auto operator*(const LorentzTensor & t, const LorentzTensor & u) -> decltype(t.xx()*u.xx()) ; /** * Addition. */ LorentzTensor operator+(const LorentzTensor & in) const { return LorentzTensor(xx()+in.xx(),xy()+in.xy(),xz()+in.xz(),xt()+in.xt(), yx()+in.yx(),yy()+in.yy(),yz()+in.yz(),yt()+in.yt(), zx()+in.zx(),zy()+in.zy(),zz()+in.zz(),zt()+in.zt(), tx()+in.tx(),ty()+in.ty(),tz()+in.tz(),tt()+in.tt()); } /** * Subtraction. */ LorentzTensor operator-(const LorentzTensor & in) const { return LorentzTensor(xx()-in.xx(),xy()-in.xy(),xz()-in.xz(),xt()-in.xt(), yx()-in.yx(),yy()-in.yy(),yz()-in.yz(),yt()-in.yt(), zx()-in.zx(),zy()-in.zy(),zz()-in.zz(),zt()-in.zt(), tx()-in.tx(),ty()-in.ty(),tz()-in.tz(),tt()-in.tt()); } /** * Trace */ complex trace() const { return _tensor[3][3]-_tensor[0][0]-_tensor[1][1]-_tensor[2][2]; } //@} /** * Various dot products */ //@{ /** * First index dot product with polarization vector */ LorentzVector > preDot (const LorentzPolarizationVector & vec) const { LorentzVector > output; output.setX(vec.t()*_tensor[3][0]-vec.x()*_tensor[0][0]- vec.y()*_tensor[1][0]-vec.z()*_tensor[2][0]); output.setY(vec.t()*_tensor[3][1]-vec.x()*_tensor[0][1]- vec.y()*_tensor[1][1]-vec.z()*_tensor[2][1]); output.setZ(vec.t()*_tensor[3][2]-vec.x()*_tensor[0][2]- vec.y()*_tensor[1][2]-vec.z()*_tensor[2][2]); output.setT(vec.t()*_tensor[3][3]-vec.x()*_tensor[0][3]- vec.y()*_tensor[1][3]-vec.z()*_tensor[2][3]); return output; } /** * Second index dot product with polarization vector */ LorentzVector > postDot(const LorentzPolarizationVector & vec) const { LorentzVector > output; output.setX(vec.t()*_tensor[0][3]-vec.x()*_tensor[0][0]- vec.y()*_tensor[0][1]-vec.z()*_tensor[0][2]); output.setY(vec.t()*_tensor[1][3]-vec.x()*_tensor[1][0]- vec.y()*_tensor[1][1]-vec.z()*_tensor[1][2]); output.setZ(vec.t()*_tensor[2][3]-vec.x()*_tensor[2][0]- vec.y()*_tensor[2][1]-vec.z()*_tensor[2][2]); output.setT(vec.t()*_tensor[3][3]-vec.x()*_tensor[3][0]- vec.y()*_tensor[3][1]-vec.z()*_tensor[3][2]); return output; } /** * First index dot product with momentum */ auto preDot (const Lorentz5Momentum & vec) const - -> LorentzVector + -> LorentzVectorxx())> { - LorentzVector output; + LorentzVectorxx())> output; output.setX(vec.t()*_tensor[3][0]-vec.x()*_tensor[0][0]- vec.y()*_tensor[1][0]-vec.z()*_tensor[2][0]); output.setY(vec.t()*_tensor[3][1]-vec.x()*_tensor[0][1]- vec.y()*_tensor[1][1]-vec.z()*_tensor[2][1]); output.setZ(vec.t()*_tensor[3][2]-vec.x()*_tensor[0][2]- vec.y()*_tensor[1][2]-vec.z()*_tensor[2][2]); output.setT(vec.t()*_tensor[3][3]-vec.x()*_tensor[0][3]- vec.y()*_tensor[1][3]-vec.z()*_tensor[2][3]); return output; } /** * Second index dot product with momentum */ auto postDot(const Lorentz5Momentum & vec) const - -> LorentzVector + -> LorentzVectorxx())> { - LorentzVector output; + LorentzVectorxx())> output; output.setX(vec.t()*_tensor[0][3]-vec.x()*_tensor[0][0]- vec.y()*_tensor[0][1]-vec.z()*_tensor[0][2]); output.setY(vec.t()*_tensor[1][3]-vec.x()*_tensor[1][0]- vec.y()*_tensor[1][1]-vec.z()*_tensor[1][2]); output.setZ(vec.t()*_tensor[2][3]-vec.x()*_tensor[2][0]- vec.y()*_tensor[2][1]-vec.z()*_tensor[2][2]); output.setT(vec.t()*_tensor[3][3]-vec.x()*_tensor[3][0]- vec.y()*_tensor[3][1]-vec.z()*_tensor[3][2]); return output; } //@} private: /** * The components. */ std::array,4>,4> _tensor; }; /** * Multiplication by a complex number. */ template inline auto operator*(complex a, const LorentzTensor & t) -> LorentzTensor { return {a*t.xx(), a*t.xy(), a*t.xz(), a*t.xt(), a*t.yx(), a*t.yy(), a*t.yz(), a*t.yt(), a*t.zx(), a*t.zy(), a*t.zz(), a*t.zt(), a*t.tx(), a*t.ty(), a*t.tz(), a*t.tt()}; } /** * Multiplication by a complex number. */ template inline auto operator*(const LorentzTensor & t,complex a) -> LorentzTensor { return {a*t.xx(), a*t.xy(), a*t.xz(), a*t.xt(), a*t.yx(), a*t.yy(), a*t.yz(), a*t.yt(), a*t.zx(), a*t.zy(), a*t.zz(), a*t.zt(), a*t.tx(), a*t.ty(), a*t.tz(), a*t.tt()}; } /** * Multiply a LorentzVector by a LorentzTensor. */ template inline auto operator*(const LorentzVector & v, const LorentzTensor & t) -> LorentzVector { LorentzVector outvec; outvec.setX( v.t()*t(3,0)-v.x()*t(0,0) -v.y()*t(1,0)-v.z()*t(2,0)); outvec.setY( v.t()*t(3,1)-v.x()*t(0,1) -v.y()*t(1,1)-v.z()*t(2,1)); outvec.setZ( v.t()*t(3,2)-v.x()*t(0,2) -v.y()*t(1,2)-v.z()*t(2,2)); outvec.setT( v.t()*t(3,3)-v.x()*t(0,3) -v.y()*t(1,3)-v.z()*t(2,3)); return outvec; } /** * Multiply a LorentzTensor by a LorentzVector. */ template inline auto operator*(const LorentzTensor & t, const LorentzVector & v) -> LorentzVector { LorentzVector outvec; outvec.setX( v.t()*t(0,3)-v.x()*t(0,0) -v.y()*t(0,1)-v.z()*t(0,2)); outvec.setY( v.t()*t(1,3)-v.x()*t(1,0) -v.y()*t(1,1)-v.z()*t(1,2)); outvec.setZ( v.t()*t(2,3)-v.x()*t(2,0) -v.y()*t(2,1)-v.z()*t(2,2)); outvec.setT( v.t()*t(3,3)-v.x()*t(3,0) -v.y()*t(3,1)-v.z()*t(3,2)); return outvec; } /** * Multiply a LorentzTensor by a LorentzTensor */ template inline auto operator*(const LorentzTensor & t, const LorentzTensor & u) -> decltype(t.xx()*u.xx()) { using RetT = decltype(t.xx()*u.xx()); RetT output=RetT(),temp; for(unsigned int ix=0;ix<4;++ix) { temp = t._tensor[ix][3]*u._tensor[ix][3]; for(unsigned int iy=0;iy<3;++iy) { temp+= t._tensor[ix][iy]*u._tensor[ix][iy]; } if(ix<3) output-=temp; else output+=temp; } return output; } } } #ifndef ThePEG_TEMPLATES_IN_CC_FILE #include "LorentzTensor.tcc" #endif #endif