herwig is hosted by Hepforge, IPPP Durham
Herwig  7.2.1
SpinorHelicity.h
1 // -*- C++ -*-
2 //
3 // SpinorHelicity.h is a part of Herwig - A multi-purpose Monte Carlo event generator
4 // Copyright (C) 2002-2019 The Herwig Collaboration
5 //
6 // Herwig is licenced under version 3 of the GPL, see COPYING for details.
7 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
8 //
9 #ifndef HERWIG_SpinorHelicity_H
10 #define HERWIG_SpinorHelicity_H
11 
12 #include "ThePEG/Config/Complex.h"
14 
15 #include <boost/operators.hpp>
16 
17 namespace Herwig {
18 
19 using namespace ThePEG;
20 
21 namespace SpinorHelicity {
22 
30  struct PlusSpinorTag {};
31 
39  struct MinusSpinorTag {};
40 
49 
58 
66  template<class Value>
68 
69  typedef decltype(sqr(std::declval<Value>())) ResultType;
70  typedef complex<ResultType> ComplexResultType;
72 
73  };
74 
82  template<class Type>
84 
85  // specialize for |p>
86  template<>
88 
89  template<class Value, class MValue>
90  static pair<complex<Value>,complex<Value> >
91  components(const LorentzVector<MValue>& p) {
92  if ( p.t() < ZERO ) {
93  pair<complex<Value>,complex<Value> > res =
94  components<Value,MValue>(-p);
95  // do not revert to *=, breaks with XCode 5.1
96  res.first = res.first * Complex(0.,1.);
97  res.second = res.second * Complex(0.,1.);
98  return res;
99  }
100  Energy pPlus = p.t() + p.x();
101  if ( abs(pPlus) < 1.e-10 * GeV ) {
102  return make_pair(complex<Value>(ZERO),
103  complex<Value>(sqrt(2.*p.t())));
104  }
105  return make_pair(complex<Value>(sqrt(pPlus)),
106  complex<Value>(p.z()/sqrt(pPlus),p.y()/sqrt(pPlus)));
107  }
108 
109  };
110 
111  // specialize for |p]
112  template<>
114 
115  template<class Value, class MValue>
116  static pair<complex<Value>,complex<Value> >
117  components(const LorentzVector<MValue>& p) {
118  if ( p.t() < ZERO ) {
119  pair<complex<Value>,complex<Value> > res =
120  components<Value,MValue>(-p);
121  // do not revert to *=, breaks with XCode 5.1
122  res.first = res.first * Complex(0.,1.);
123  res.second = res.second * Complex(0.,1.);
124  return res;
125  }
126  Energy pPlus = p.t() + p.x();
127  if ( abs(pPlus) < 1.e-10 * GeV ) {
128  return make_pair(complex<Value>(sqrt(2.*p.t())),
129  complex<Value>(ZERO));
130  }
131  return make_pair(complex<Value>(p.z()/sqrt(pPlus),-p.y()/sqrt(pPlus)),
132  -complex<Value>(sqrt(pPlus)));
133  }
134 
135  };
136 
137  // specialize for <p|
138  template<>
140 
143 
144  template<class Value, class MValue>
145  static pair<complex<Value>,complex<Value> >
146  components(const LorentzVector<MValue>& p) {
147  pair<complex<Value>,complex<Value> > res =
148  WeylSpinorTraits<PlusSpinorTag>::template components<Value>(p);
149  res.first = -res.first;
150  swap(res.first,res.second);
151  return res;
152  }
153 
154  };
155 
156  // specialize for [p|
157  template<>
159 
161  typedef PlusSpinorTag BarSpinorTag;
162 
163  template<class Value, class MValue>
164  static pair<complex<Value>,complex<Value> >
165  components(const LorentzVector<MValue>& p) {
166  pair<complex<Value>,complex<Value> > res =
167  WeylSpinorTraits<MinusSpinorTag>::template components<Value>(p);
168  res.second = -res.second;
169  swap(res.first,res.second);
170  return res;
171  }
172 
173  };
174 
182  template<class Type, class Value>
183  class WeylSpinor {
184 
185  public:
186 
187  typedef complex<Value> ComplexType;
188  typedef pair<ComplexType,ComplexType> ComponentsType;
189  typedef Type Tag;
191  typedef Value ValueType;
192 
193  private:
194 
198  ComponentsType theComponents;
199 
200  public:
201 
205  explicit WeylSpinor(const ComponentsType& c = ComponentsType())
206  : theComponents(c) {}
207 
211  template<class MValue>
213  : theComponents(Traits::template components<Value>(p)) {}
214 
218  const ComponentsType& components() const { return theComponents; }
219 
223  const ComplexType& s1() const { return theComponents.first; }
224 
228  const ComplexType& s2() const { return theComponents.second; }
229 
230  };
231 
234 
237 
240 
243 
251  template<class Type, class Value>
253  : public boost::addable<SpinorProduct<Type,Value> >,
254  public boost::subtractable<SpinorProduct<Type,Value> >,
255  public boost::multipliable<SpinorProduct<Type,Value>, double>,
256  public boost::multipliable<SpinorProduct<Type,Value>, complex<double> > {
257 
258  public:
259 
260  typedef typename SpinorMultiplicationTraits<Value>::ComplexResultType ResultType;
262  typedef typename WeylSpinorTraits<Type>::ConjugateSpinorTag RightSpinorTag;
264 
265  private:
266 
270  ResultType theResult;
271 
272  public:
273 
279  explicit SpinorProduct(const LeftSpinorType& left,
280  const RightSpinorType& right)
281  : theResult(left.s1()*right.s1()+left.s2()*right.s2()) {}
282 
286  operator ResultType() const { return theResult; }
287 
291  ResultType eval() const { return theResult; }
292 
293  public:
294 
295  SpinorProduct& operator+= (const SpinorProduct& other) {
296  theResult += other.theResult;
297  return *this;
298  }
299 
300  SpinorProduct& operator-= (const SpinorProduct& other) {
301  theResult -= other.theResult;
302  return *this;
303  }
304 
305  SpinorProduct& operator*= (double x) {
306  theResult *= x;
307  return *this;
308  }
309 
310  SpinorProduct& operator*= (complex<double> x) {
311  theResult *= x;
312  return *this;
313  }
314 
315  };
316 
319 
322 
330  template<class Type, class Value>
332  : public boost::addable<SpinorCurrent<Type,Value> >,
333  public boost::subtractable<SpinorCurrent<Type,Value> >,
334  public boost::multipliable<SpinorCurrent<Type,Value>, double>,
335  public boost::multipliable<SpinorCurrent<Type,Value>, complex<double> > {
336 
337  public:
338 
341  typedef typename WeylSpinorTraits<Type>::BarSpinorTag RightSpinorTag;
343 
344  private:
345 
346  ResultType theResult;
347 
352  const WeylSpinor<PlusSpinorTag,Value>& right) {
353  return
354  ResultType(right.s1()*left.s1()-right.s2()*left.s2(),
355  complex<double>(0.,1.)*(right.s1()*left.s2()-right.s2()*left.s1()),
356  right.s1()*left.s2()+right.s2()*left.s1(),
357  right.s1()*left.s1()+right.s2()*left.s2());
358  }
359 
364  const WeylSpinor<MinusSpinorTag,Value>& right) {
365  return
366  ResultType(-right.s1()*left.s1()+right.s2()*left.s2(),
367  -complex<double>(0.,1.)*(right.s1()*left.s2()-right.s2()*left.s1()),
368  -right.s1()*left.s2()-right.s2()*left.s1(),
369  right.s1()*left.s1()+right.s2()*left.s2());
370  }
371 
372  public:
373 
378  explicit SpinorCurrent(const LeftSpinorType& left,
379  const RightSpinorType& right)
380  : theResult(evaluate(left,right)) {}
381 
385  operator ResultType() const { return theResult; }
386 
390  ResultType eval() const { return theResult; }
391 
392  public:
393 
394  SpinorCurrent& operator+= (const SpinorCurrent& other) {
395  theResult += other.theResult;
396  return *this;
397  }
398 
399  SpinorCurrent& operator-= (const SpinorCurrent& other) {
400  theResult -= other.theResult;
401  return *this;
402  }
403 
404  SpinorCurrent& operator*= (double x) {
405  theResult *= x;
406  return *this;
407  }
408 
409  SpinorCurrent& operator*= (complex<double> x) {
410  theResult *= x;
411  return *this;
412  }
413 
414  };
415 
418 
421 
425  template<class T>
426  auto abs2(const complex<T>& x) -> decltype((x*conj(x)).real())
427  {
428  return (x*conj(x)).real();
429  }
430 
431 }
432 
433 }
434 
435 #endif // HERWIG_SpinorHelicity_H
WeylSpinor(const ComponentsType &c=ComponentsType())
Construct from components.
double sqrt(int x)
std::complex< double > Complex
WeylSpinor(const LorentzVector< MValue > &p)
Construct from momentum.
SpinorProduct(const LeftSpinorType &left, const RightSpinorType &right)
Construct from two spinors; note that the spinor metric is included, when constructing spinors...
constexpr auto sqr(const T &x) -> decltype(x *x)
ResultType evaluate(const WeylSpinor< MinusConjugateSpinorTag, Value > &left, const WeylSpinor< PlusSpinorTag, Value > &right)
Calculate [p|^|q>
ResultType eval() const
Return result.
Base class for Weyl spinors.
ComponentsType theComponents
The components.
SpinorCurrent(const LeftSpinorType &left, const RightSpinorType &right)
Construct from two spinors.
ResultType eval() const
Return result.
const ComponentsType & components() const
Return the components.
Helpers for commonly encountered types.
const ComplexType & s1() const
Return the first component.
-*- C++ -*-
const ComplexType & s2() const
Return the second component.
ResultType evaluate(const WeylSpinor< PlusConjugateSpinorTag, Value > &left, const WeylSpinor< MinusSpinorTag, Value > &right)
Calculate <p|^|q].
constexpr ZeroUnit ZERO