herwig is hosted by Hepforge, IPPP Durham
Herwig  7.2.1
MEQCD2to2Fast.h
1 // -*- C++ -*-
2 //
3 // MEQCD2to2Fast.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_MEQCD2to2Fast_H
10 #define HERWIG_MEQCD2to2Fast_H
11 //
12 // This is the declaration of the MEQCD2to2Fast class.
13 //
14 
15 #include "Herwig/MatrixElement/HwMEBase.h"
16 #include "ThePEG/Repository/UseRandom.h"
17 
18 namespace Herwig {
19 using namespace ThePEG;
20 using namespace ThePEG::Helicity;
21 
31 class MEQCD2to2Fast: public HwMEBase {
32 
33 public:
34 
38  MEQCD2to2Fast() :_maxflavour(5),_process(0),_strictFlavourScheme(false) {
39  massOption(vector<unsigned int>(2,0));
40  }
41 
48  virtual unsigned int orderInAlphaS() const;
49 
54  virtual unsigned int orderInAlphaEW() const;
55 
63  virtual double me2() const;
64 
68  virtual Energy2 scale() const;
69 
73  virtual void getDiagrams() const;
74 
83  virtual Selector<DiagramIndex> diagrams(const DiagramVector & dv) const;
84 
93  colourGeometries(tcDiagPtr diag) const;
95 
96 
97 public:
98 
105  void persistentOutput(PersistentOStream & os) const;
106 
112  void persistentInput(PersistentIStream & is, int version);
114 
121  static void Init();
122 
123 protected:
124 
132  double gg2ggME() const {
133  Energy2 u(uHat()),t(tHat()),s(sHat());
134  double output = 9./4.*(3.-t*u/s/s-s*u/t/t-s*t/u/u);
135  double flow[3]={(1.-u*t/s/s-s*t/u/u+t*t/s/u),
136  (1.-t*u/s/s-s*u/t/t+u*u/s/t),
137  (1.-t*s/u/u-u*s/t/t+s*s/u/t)};
138  _flow = 1+UseRandom::rnd3(flow[0],flow[1],flow[2]);
139  double diag[3]={(sqr(u)+sqr(t))/sqr(s),
140  (sqr(s)+sqr(u))/sqr(t),
141  (sqr(s)+sqr(t))/sqr(u)};
142  if(_flow==1) diag[1]=0;
143  else if(_flow==2) diag[2]=0;
144  else if(_flow==3) diag[0]=0;
145  _diagram=1+UseRandom::rnd3(diag[0],diag[1],diag[2]);
146  return output;
147  }
148 
152  double gg2qqbarME() const {
153  Energy2 u(uHat()),t(tHat()),s(sHat());
154  Energy4 u2(sqr(u)),t2(sqr(t)),s2(sqr(s));
155  double output =(1./6./u/t-3./8./s2)*(t2+u2);
156  double flow[2]={u2/(u2+t2),t2/(u2+t2)};
157  _flow = 1+UseRandom::rnd2(flow[0],flow[1]);
158  _diagram=3+_flow;
159  return output;
160  }
161 
165  double qqbar2ggME() const {
166  Energy2 u(uHat()),t(tHat()),s(sHat());
167  Energy4 s2(sqr(s)),u2(sqr(u)),t2(sqr(t));
168  double output = 0.5*(32./27./u/t-8./3./s2)*(t2+u2);
169  double flow[2] = {u2/(u2+t2),t2/(t2+u2)};
170  _flow=1+UseRandom::rnd2(flow[0],flow[1]);
171  _diagram=6+_flow;
172  return output;
173  }
174 
178  double qg2qgME() const {
179  Energy2 u(uHat()),t(tHat()),s(sHat());
180  Energy4 s2(sqr(s)),u2(sqr(u)),t2(sqr(t));
181  double output = (-4./9./s/u+1./t2)*(s2+u2);
182  double flow[2]={u2/(s2+u2),s2/(s2+u2)};
183  _flow=1+UseRandom::rnd2(flow[0],flow[1]);
184  _diagram=9+_flow;
185  return output;
186  }
187 
191  double qbarg2qbargME() const {
192  // scale
193  Energy2 u(uHat()),t(tHat()),s(sHat());
194  Energy4 u2(sqr(u)),s2(sqr(s)); // t2(sqr(t))
195  double flow[2]={u2/(s2+u2),s2/(s2+u2)};
196  _flow=1+UseRandom::rnd2(flow[0],flow[1]);
197  _diagram=12+_flow;
198  return (-4./9./s/u+1./t/t)*(s*s+u*u);
199  }
200 
204  double qq2qqME() const {
205  Energy2 u(uHat()),t(tHat());
206  Energy4 s2(sqr(sHat())),u2(sqr(u)),t2(sqr(t));
207  double output;
208  if(mePartonData()[0]->id()==mePartonData()[1]->id()) {
209  output = 0.5*(4./9.*((s2+u2)/t2+(s2+t2)/u2)
210  -8./27.*s2/u/t);
211  double flow[2]={(s2+u2)/t2,(s2+t2)/u2};
212  _flow=1+UseRandom::rnd2(flow[0],flow[1]);
213  }
214  else {
215  output = 4./9.*(s2+u2)/t2;
216  _flow=2;
217  }
218  _diagram = 15+_flow;
219  return output;
220  }
221 
225  double qbarqbar2qbarqbarME() const {
226  Energy2 u(uHat()),t(tHat());
227  Energy4 u2(sqr(u)),t2(sqr(t)),s2(sqr(sHat()));
228  double output;
229  if(mePartonData()[0]->id()==mePartonData()[1]->id()) {
230  output = 0.5*(4./9.*((s2+u2)/t2+(s2+t2)/u2)
231  -8./27.*s2/u/t);
232  double flow[2]={(s2+u2)/t2,(s2+t2)/u2};
233  _flow=1+UseRandom::rnd2(flow[0],flow[1]);
234  }
235  else {
236  output = 4./9.*(s2+u2)/t2;
237  _flow = 2;
238  }
239  _diagram = 17+_flow;
240  // final part of colour and spin factors
241  return output;
242  }
243 
247  double qqbar2qqbarME() const {
248  // type of process
249  bool diagon[2]={mePartonData()[0]->id()== -mePartonData()[1]->id(),
250  mePartonData()[0]->id()== mePartonData()[2]->id()};
251  // scale
252  Energy2 u(uHat()),t(tHat()),s(sHat());
253  Energy4 s2(sqr(s)),t2(sqr(t)),u2(sqr(u));
254  double output;
255  if(diagon[0]&&diagon[1]) {
256  output= (4./9.*((s2+u2)/t2+(u2+t2)/s2)
257  -8./27.*u2/s/t);
258  double flow[2]={(t2+u2)/s2,(s2+u2)/t2};
259  _flow=1+UseRandom::rnd2(flow[0],flow[1]);
260  }
261  else if(diagon[0]) {
262  output = (4./9.*(t2+u2)/s2);
263  _flow=1;
264  }
265  else {
266  output = (4./9.*(s2+u2)/t2);
267  _flow=2;
268  }
269  _diagram=19+_flow;
270  return output;
271  }
273 
274 protected:
275 
282  virtual IBPtr clone() const;
283 
288  virtual IBPtr fullclone() const;
290 
291 protected:
292 
298  virtual void doinit();
299 
301 
302 private:
303 
308  MEQCD2to2Fast & operator=(const MEQCD2to2Fast &) = delete;
309 
310 private:
311 
315  unsigned int _maxflavour;
316 
320  unsigned int _process;
321 
325  mutable unsigned int _flow;
326 
330  mutable unsigned int _diagram;
331 
336 
337 };
338 
339 }
340 
341 #endif /* HERWIG_MEQCD2to2Fast_H */
unsigned int _maxflavour
Maximum numbere of quark flavours to include.
unsigned int _process
Processes to include.
static int rnd2(double p0, double p1)
double qqbar2qqbarME() const
Matrix element for .
ThePEG::Ptr< InterfacedBase >::pointer IBPtr
constexpr auto sqr(const T &x) -> decltype(x *x)
double qg2qgME() const
Matrix element for .
MEQCD2to2Fast()
The default constructor.
Definition: MEQCD2to2Fast.h:38
double gg2ggME() const
Members to calculate the matrix elements.
bool _strictFlavourScheme
Exclude contributions with massive incominbg quarks.
double qq2qqME() const
Matrix element for .
double qbarqbar2qbarqbarME() const
Matrix element for .
unsigned int _flow
Colour flow.
The MEQCD2to2Fast class implements the matrix elements for QCD scattering processes using hard coded...
Definition: MEQCD2to2Fast.h:31
unsigned int _diagram
Diagram.
double qbarg2qbargME() const
Matrix elements for .
double qqbar2ggME() const
Matrix element for .
vector< DiagPtr > DiagramVector
-*- C++ -*-
static int rnd3(double p0, double p1, double p2)
double gg2qqbarME() const
Matrix element for .
The HwMEBase class serves a number of purposes.
Definition: HwMEBase.h:35