herwig is hosted by Hepforge, IPPP Durham
Herwig  7.2.1
GeneralDecayMatrixElement.h
1 // -*- C++ -*-
2 //
3 // GeneralDecayMatrixElement.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_GeneralDecayMatrixElement_H
10 #define HERWIG_GeneralDecayMatrixElement_H
11 //
12 // This is the declaration of the GeneralDecayMatrixElement class.
13 
14 #include "DecayMatrixElement.h"
15 #include "GeneralDecayMatrixElement.fh"
16 
17 namespace Herwig {
18 
19 using namespace ThePEG;
20 
21 
33 
34 public:
35 
42 
50  PDT::Spin outspin1,PDT::Spin outspin2)
51  : DecayMatrixElement(2,spinin) {
52  outspin().push_back(outspin1);
53  outspin().push_back(outspin2);
54  setMESize();
55  }
56 
65  PDT::Spin outspin2,PDT::Spin outspin3)
66  : DecayMatrixElement(3,spinin) {
67  outspin().push_back(outspin1);
68  outspin().push_back(outspin2);
69  outspin().push_back(outspin3);
70  setMESize();
71  }
72 
82  PDT::Spin outspin3,PDT::Spin outspin4)
83  : DecayMatrixElement(4,spinin) {
84  outspin().push_back(outspin1);
85  outspin().push_back(outspin2);
86  outspin().push_back(outspin3);
87  outspin().push_back(outspin4);
88  setMESize();
89  }
90 
101  PDT::Spin outspin3,PDT::Spin outspin4,PDT::Spin outspin5)
102  : DecayMatrixElement(5,spinin) {
103  outspin().push_back(outspin1);
104  outspin().push_back(outspin2);
105  outspin().push_back(outspin3);
106  outspin().push_back(outspin4);
107  outspin().push_back(outspin5);
108  setMESize();
109  }
110 
122  PDT::Spin outspin3,PDT::Spin outspin4,PDT::Spin outspin5,
123  PDT::Spin outspin6)
124  : DecayMatrixElement(6,spinin) {
125  outspin().push_back(outspin1);
126  outspin().push_back(outspin2);
127  outspin().push_back(outspin3);
128  outspin().push_back(outspin4);
129  outspin().push_back(outspin5);
130  outspin().push_back(outspin6);
131  setMESize();
132  }
133 
139  GeneralDecayMatrixElement(PDT::Spin spinin,vector<PDT::Spin> spinout)
140  : DecayMatrixElement(spinout.size(),spinin) {
141  outspin() = spinout;
142  setMESize();
143  }
144 
149  GeneralDecayMatrixElement(vector<PDT::Spin> extspin)
150  : DecayMatrixElement(int(extspin.size())-1,extspin[0]) {
151  outspin() = vector<PDT::Spin>(extspin.begin()+1,extspin.end());
152  setMESize();
153  }
155 
156 public:
157 
166  RhoDMatrix calculateDMatrix(const vector<RhoDMatrix> & rhoout) const;
167 
174  RhoDMatrix calculateRhoMatrix(int ipart,const RhoDMatrix & rhoin,
175  const vector<RhoDMatrix> & rhoout) const;
176 
182  Complex contract(const RhoDMatrix & rhoin) const;
183 
190  Complex contract(const GeneralDecayMatrixElement & con,
191  const RhoDMatrix & rhoin);
193 
194 public:
195 
206  Complex operator () (unsigned int inhel,unsigned int outhel1,
207  unsigned int outhel2) const {
208  unsigned int iloc = inhel*constants_[1]+
209  outhel1*constants_[2]+outhel2*constants_[3];
210  assert(outspin().size()==2&&iloc<matrixElement_.size());
211  return matrixElement_[iloc];
212  }
213 
220  Complex & operator () (unsigned int inhel,unsigned int outhel1,
221  unsigned int outhel2) {
222  unsigned int iloc = inhel*constants_[1]+
223  outhel1*constants_[2]+outhel2*constants_[3];
224  assert(outspin().size()==2&&iloc<matrixElement_.size());
225  return matrixElement_[iloc];
226  }
227 
235  Complex operator () (unsigned int inhel,unsigned int outhel1,
236  unsigned int outhel2,unsigned int outhel3) const {
237  unsigned int iloc = inhel*constants_[1]+outhel1*constants_[2]+
238  outhel2*constants_[3]+outhel3*constants_[4];
239  assert(outspin().size()==3&&iloc<matrixElement_.size());
240  return matrixElement_[iloc];
241  }
242 
250  Complex & operator () (unsigned int inhel,unsigned int outhel1,
251  unsigned int outhel2,unsigned int outhel3) {
252  unsigned int iloc = inhel*constants_[1]+outhel1*constants_[2]+
253  outhel2*constants_[3]+outhel3*constants_[4];
254  assert(outspin().size()==3&&iloc<matrixElement_.size());
255  return matrixElement_[iloc];
256  }
257 
266  Complex operator () (unsigned int inhel,unsigned int outhel1,
267  unsigned int outhel2,unsigned int outhel3,
268  unsigned int outhel4) const {
269  vector<unsigned int> itemp(5);
270  itemp[0]=inhel ; itemp[1]=outhel1;
271  itemp[2]=outhel2; itemp[3]=outhel3;
272  itemp[4]=outhel4 ;
273  return (*this)(itemp);
274  }
275 
284  Complex & operator () (unsigned int inhel,unsigned int outhel1,
285  unsigned int outhel2,unsigned int outhel3,
286  unsigned int outhel4) {
287  vector<unsigned int> itemp(5);
288  itemp[0]=inhel ; itemp[1]=outhel1;
289  itemp[2]=outhel2; itemp[3]=outhel3;
290  itemp[4]=outhel4;
291  return (*this)(itemp);
292  }
293 
303  Complex operator () (unsigned int inhel,unsigned int outhel1,
304  unsigned int outhel2,unsigned int outhel3,
305  unsigned int outhel4,unsigned int outhel5) const {
306  vector<unsigned int> itemp(6);
307  itemp[0]=inhel ; itemp[1]=outhel1;
308  itemp[2]=outhel2; itemp[3]=outhel3;
309  itemp[4]=outhel4 ;itemp[5]=outhel5;
310  return (*this)(itemp);
311  }
312 
322  Complex & operator () (unsigned int inhel,unsigned int outhel1,
323  unsigned int outhel2,unsigned int outhel3,
324  unsigned int outhel4,unsigned int outhel5) {
325  vector<unsigned int> itemp(6);
326  itemp[0]=inhel ; itemp[1]=outhel1;
327  itemp[2]=outhel2; itemp[3]=outhel3;
328  itemp[4]=outhel4; itemp[5]=outhel5;
329  return (*this)(itemp);
330  }
331 
342  Complex operator () (unsigned int inhel,unsigned int outhel1,
343  unsigned int outhel2,unsigned int outhel3,
344  unsigned int outhel4,unsigned int outhel5,
345  unsigned int outhel6) const {
346  vector<unsigned int> itemp(7);
347  itemp[0]=inhel ; itemp[1]=outhel1;
348  itemp[2]=outhel2; itemp[3]=outhel3;
349  itemp[4]=outhel4 ;itemp[5]=outhel5;
350  itemp[6]=outhel6;
351  return (*this)(itemp);
352  }
353 
364  Complex & operator () (unsigned int inhel,unsigned int outhel1,
365  unsigned int outhel2,unsigned int outhel3,
366  unsigned int outhel4,unsigned int outhel5,
367  unsigned int outhel6) {
368  vector<unsigned int> itemp(7);
369  itemp[0]=inhel ; itemp[1]=outhel1;
370  itemp[2]=outhel2; itemp[3]=outhel3;
371  itemp[4]=outhel4; itemp[5]=outhel5;
372  itemp[6]=outhel6;
373  return (*this)(itemp);
374  }
375 
380  Complex operator () (const vector<unsigned int> & in) const {
381  assert(in.size()==outspin().size()+1);
382  // contribution for the incoming particle
383  unsigned int iloc(in[0]*constants_[1]);
384  // contributions for the outgoing particles
385  for(unsigned int ix=1;ix<in.size();++ix){iloc+=in[ix]*constants_[ix+1];}
386  assert(iloc<matrixElement_.size());
387  return matrixElement_[iloc];
388  }
389 
394  Complex & operator () (const vector<unsigned int> & in) {
395  assert(in.size()==outspin().size()+1);
396  // contribution for the incoming particle
397  unsigned int iloc(in[0]*constants_[1]);
398  // contributions for the outgoing particles
399  for(unsigned int ix=1;ix<in.size();++ix){iloc+=in[ix]*constants_[ix+1];}
400  assert(iloc<matrixElement_.size());
401  return matrixElement_[iloc];
402  }
404 
408  void zero() {
409  for(unsigned int ix=0;ix<matrixElement_.size();++ix)
410  matrixElement_[ix]=0.;
411  }
412 
416  void setMESize() {
417  int isize = inspin();
418  for(unsigned int ix=0;ix<outspin().size();++ix) isize*=outspin()[ix];
419  matrixElement_.resize(isize,0.);
420  // set up the constants for the mapping of helicity to vector index
421  constants_.resize(outspin().size()+2);
422  int temp=1;
423  for(unsigned int ix=outspin().size();ix>0;--ix) {
424  temp*=outspin()[ix-1];constants_[ix]=temp;
425  }
426  temp *= inspin();
427  constants_[0]=temp;
428  constants_[outspin().size()+1]=1;
429  }
430 
431 private:
432 
436  mutable vector<Complex> matrixElement_;
437 
441  mutable vector<unsigned int> constants_;
442 
443 };
444 
445 }
446 
447 #endif /* HERWIG_GeneralDecayMatrixElement_H */
GeneralDecayMatrixElement(PDT::Spin spinin, PDT::Spin outspin1, PDT::Spin outspin2, PDT::Spin outspin3)
Constructor for three body decay.
GeneralDecayMatrixElement()
Default constructor.
std::complex< double > Complex
GeneralDecayMatrixElement(PDT::Spin spinin, vector< PDT::Spin > spinout)
Constructor for arbitary body decay.
GeneralDecayMatrixElement(PDT::Spin spinin, PDT::Spin outspin1, PDT::Spin outspin2)
Constructor for two body decay.
void zero()
Member to zero all the elements for the matrix element.
GeneralDecayMatrixElement(PDT::Spin spinin, PDT::Spin outspin1, PDT::Spin outspin2, PDT::Spin outspin3, PDT::Spin outspin4, PDT::Spin outspin5, PDT::Spin outspin6)
Constructor for six body decay.
vector< unsigned int > constants_
Constants needed to map the index of the vector to a helicity structure.
GeneralDecayMatrixElement(vector< PDT::Spin > extspin)
Constructor for arbitary body decay.
Here is the documentation of the DecayMatrixElement class.
-*- C++ -*-
void setMESize()
Set the size of the vector containing the matrix element.
GeneralDecayMatrixElement(PDT::Spin spinin, PDT::Spin outspin1, PDT::Spin outspin2, PDT::Spin outspin3, PDT::Spin outspin4, PDT::Spin outspin5)
Constructor for five body decay.
vector< Complex > matrixElement_
Storage of the matrix element, a vector is better for memory usage.
GeneralDecayMatrixElement(PDT::Spin spinin, PDT::Spin outspin1, PDT::Spin outspin2, PDT::Spin outspin3, PDT::Spin outspin4)
Constructor for four body decay.