herwig is hosted by Hepforge, IPPP Durham
Herwig  7.1.5
DipoleSplittingKernel.h
1 // -*- C++ -*-
2 //
3 // DipoleSplittingKernel.h is a part of Herwig - A multi-purpose Monte Carlo event generator
4 // Copyright (C) 2002-2017 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_DipoleSplittingKernel_H
10 #define HERWIG_DipoleSplittingKernel_H
11 //
12 // This is the declaration of the DipoleSplittingKernel class.
13 //
14 
15 #include "ThePEG/Handlers/HandlerBase.h"
16 #include "ThePEG/StandardModel/AlphaSBase.h"
17 #include "ThePEG/PDF/PDF.h"
18 
19 #include "Herwig/Shower/Dipole/Utility/PDFRatio.h"
20 #include "Herwig/Shower/Dipole/Base/DipoleSplittingInfo.h"
21 #include "Herwig/Shower/Dipole/Kinematics/DipoleSplittingKinematics.h"
22 
23 namespace Herwig {
24 
25 using namespace ThePEG;
26 
38 
39 public:
40 
47 
51  virtual ~DipoleSplittingKernel();
53 
54 public:
55 
59  Ptr<AlphaSBase>::tptr alphaS() const { return theAlphaS; }
60 
64  void alphaS(Ptr<AlphaSBase>::tptr ap) { theAlphaS = ap; }
65 
69  Ptr<DipoleSplittingKinematics>::tptr splittingKinematics() const {
70  return theSplittingKinematics;
71  }
72 
76  Ptr<DipoleMCCheck>::ptr mcCheck() const { return theMCCheck; }
77 
81  void splittingKinematics(Ptr<DipoleSplittingKinematics>::tptr sp) {
82  theSplittingKinematics = sp;
83  }
84 
88  Ptr<PDFRatio>::tptr pdfRatio() const { return thePDFRatio; }
89 
93  void pdfRatio(Ptr<PDFRatio>::tptr sp) { thePDFRatio = sp; }
94 
103  virtual int nDimAdditional() const { return 0; }
104 
108  void renormalizationScaleFreeze(Energy s) { theRenormalizationScaleFreeze = s; }
109 
113  void factorizationScaleFreeze(Energy s) { theFactorizationScaleFreeze = s; }
114 
118  Energy renormalizationScaleFreeze() const { return theRenormalizationScaleFreeze; }
119 
123  Energy factorizationScaleFreeze() const { return theFactorizationScaleFreeze; }
124 
125 public:
126 
131  virtual bool canHandle(const DipoleIndex&) const = 0;
132 
138  virtual bool canHandleEquivalent(const DipoleIndex& a,
139  const DipoleSplittingKernel& sk,
140  const DipoleIndex& b) const = 0;
141 
146  virtual tcPDPtr emitter(const DipoleIndex&) const = 0;
147 
152  virtual tcPDPtr emission(const DipoleIndex&) const = 0;
153 
158  virtual tcPDPtr spectator(const DipoleIndex&) const = 0;
159 
164  PDPtr flavour() const { return theFlavour; }
165 
170  bool strictLargeN() const { return theStrictLargeN; }
171 
172 public:
173 
178  virtual void startPresampling(const DipoleIndex&) {
179  presampling = true;
180  }
181 
186  virtual void stopPresampling(const DipoleIndex&) {
187  presampling = false;
188  }
189 
194  unsigned long presamplingPoints() const { return thePresamplingPoints; }
195 
200  unsigned long maxtry() const { return theMaxtry; }
201 
206  unsigned long freezeGrid() const { return theFreezeGrid; }
207 
212  void freezeGrid(unsigned long n) { theFreezeGrid = n; }
213 
217  void detuning(double d) { theDetuning = d; }
218 
222  double detuning() const { return theDetuning; }
223 
228  virtual double evaluate(const DipoleSplittingInfo&) const = 0;
229 
233  void clearAlphaPDFCache() const {
234  theAlphaSCache.clear();
235  thePDFCache.clear();
236  }
237 
242  virtual void accept(const DipoleSplittingInfo&, double, double, map<string,double>&) const;
243 
248  virtual void veto(const DipoleSplittingInfo&, double, double, map<string,double>&) const;
249 
257  virtual bool haveOverestimate(const DipoleSplittingInfo&) const { return false; }
258 
263  virtual double overestimate(const DipoleSplittingInfo&) const { return -1.; }
264 
271  virtual double invertOverestimateIntegral(const DipoleSplittingInfo&, double) const {
272  return -1.;
273  }
274 
278  bool useThisKernel() const {
279  return theUseThisKernel;
280  }
281 
282 public:
283 
287  double factorizationScaleFactor() const { return theFactorizationScaleFactor; }
288 
292  void factorizationScaleFactor(double f) { theFactorizationScaleFactor = f; }
293 
297  double renormalizationScaleFactor() const { return theRenormalizationScaleFactor; }
298 
302  void renormalizationScaleFactor(double f) { theRenormalizationScaleFactor = f; }
303 
304 protected:
305 
309  double alphaPDF(const DipoleSplittingInfo&,
310  Energy optScale = ZERO,
311  double rScaleFactor = 1.0,
312  double fScaleFactor = 1.0) const;
313 
318  bool virtualitySplittingScale() const { return theVirtualitySplittingScale; }
319 
320 public:
321 
328  void persistentOutput(PersistentOStream & os) const;
329 
335  void persistentInput(PersistentIStream & is, int version);
337 
344  static void Init();
345 
346 
347 // If needed, insert declarations of virtual function defined in the
348 // InterfacedBase class here (using ThePEG-interfaced-decl in Emacs).
349 
350 private:
351 
356 
362 
366  Ptr<DipoleSplittingKinematics>::ptr theSplittingKinematics;
367 
372  Ptr<PDFRatio>::ptr thePDFRatio;
373 
378  unsigned long thePresamplingPoints;
379 
384  unsigned long theMaxtry;
385 
390  unsigned long theFreezeGrid;
391 
395  double theDetuning;
396 
402 
406  Ptr<DipoleMCCheck>::ptr theMCCheck;
407 
413 
418 
423 
428 
433 
439 
444  unsigned int theCMWScheme=0;
445 
449  mutable map<double,double> theAlphaSCache;
450 
454  mutable map<double,double> thePDFCache;
455 
460 
464  bool theUseThisKernel = true;
465 
466 
467 private:
468 
474 
479  DipoleSplittingKernel & operator=(const DipoleSplittingKernel &) = delete;
480 
481 };
482 
483 }
484 
485 #include "ThePEG/Utilities/ClassTraits.h"
486 
487 namespace ThePEG {
488 
493 template <>
494 struct BaseClassTrait<Herwig::DipoleSplittingKernel,1> {
496  typedef HandlerBase NthBase;
497 };
498 
501 template <>
502 struct ClassTraits<Herwig::DipoleSplittingKernel>
503  : public ClassTraitsBase<Herwig::DipoleSplittingKernel> {
505  static string className() { return "Herwig::DipoleSplittingKernel"; }
513  static string library() { return "HwDipoleShower.so"; }
514 };
515 
518 }
519 
520 #endif /* HERWIG_DipoleSplittingKernel_H */
unsigned long freezeGrid() const
Return the number of accepted points after which the grid should be frozen.
Ptr< DipoleMCCheck >::ptr mcCheck() const
Return the mc check object.
double theFactorizationScaleFactor
The factorization scale factor.
transient_pointer tptr
Energy renormalizationScaleFreeze() const
Get the freezing value for the renormalization scale.
static AbstractClassDescription< DipoleSplittingKernel > initDipoleSplittingKernel
The static object used to initialize the description of this class.
bool strictLargeN() const
Return true, if this splitting kernel is supposed to work in a strict large-N limit, i.e.
Ptr< DipoleSplittingKinematics >::ptr theSplittingKinematics
The splitting kinematics to be used.
void freezeGrid(unsigned long n)
Set the number of accepted points after which the grid should be frozen.
virtual void stopPresampling(const DipoleIndex &)
Inform this splitting kernel, that it is not being presampled until a call to startPresampling.
void factorizationScaleFactor(double f)
Set the factorization scale factor.
bool theStrictLargeN
True, if this splitting kernel is supposed to work in a strict large-N limit, i.e.
DipoleIndex is used to index splitting generators for a particular dipole.
void splittingKinematics(Ptr< DipoleSplittingKinematics >::tptr sp)
Set the splitting kinematics object.
bool virtualitySplittingScale() const
Return true, if the virtuality of the splitting should be used as the argument of alphas rather than ...
unsigned long theFreezeGrid
Return the number of accepted points after which the grid should be frozen.
Ptr< DipoleMCCheck >::ptr theMCCheck
Pointer to a check histogram object.
Ptr< AlphaSBase >::ptr theAlphaS
The alpha_s to be used.
unsigned long thePresamplingPoints
The number of points to presample this splitting generator.
void alphaS(Ptr< AlphaSBase >::tptr ap)
Set the alpha_s to be used.
Energy theRenormalizationScaleFreeze
A freezing value for the renormalization scale.
double theDetuning
The detuning factor applied to the sampling overestimate kernel.
unsigned long maxtry() const
Return the maximum number of trials to generate a splitting.
virtual double invertOverestimateIntegral(const DipoleSplittingInfo &, double) const
Invert the integral over the overestimate w.r.t.
map< double, double > theAlphaSCache
Cache for alphas evaluations.
bool theVirtualitySplittingScale
True, if the virtuality of the splitting should be used as the argument of alphas rather than the pt...
Ptr< PDFRatio >::tptr pdfRatio() const
Return the PDFRatio object.
void renormalizationScaleFreeze(Energy s)
Set the freezing value for the renormalization scale.
Ptr< PDFRatio >::ptr thePDFRatio
An optional PDF ratio object to be used when evaluating this kernel.
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr
bool presampling
True, if we are presampling.
ThePEG::Ptr< ParticleData >::pointer PDPtr
Ptr< AlphaSBase >::tptr alphaS() const
Return the alpha_s to be used.
void detuning(double d)
Set a detuning factor to be applied to the sampling overestimate kernel.
map< double, double > thePDFCache
Cache for PDF evaluations.
void clearAlphaPDFCache() const
Clear the alphaPDF cache.
unsigned long presamplingPoints() const
Return the number of points to presample this splitting generator.
void factorizationScaleFreeze(Energy s)
Set the freezing value for the factorization scale.
virtual int nDimAdditional() const
Return the number of additional parameter random variables needed to evaluate this kernel except the ...
void pdfRatio(Ptr< PDFRatio >::tptr sp)
Set the PDFRatio object.
Energy theFactorizationScaleFreeze
A freezing value for the factorization scale.
PDPtr theFlavour
The flavour produced, if this cannot be determined from the dipole.
virtual void startPresampling(const DipoleIndex &)
Inform this splitting kernel, that it is being presampled until a call to stopPresampling.
unsigned long theMaxtry
The maximum number of trials to generate a splitting.
-*- C++ -*-
void renormalizationScaleFactor(double f)
Set the renormalization scale factor.
double theRenormalizationScaleFactor
The renormalization scale factor.
double factorizationScaleFactor() const
Get the factorization scale factor.
DipoleSplittingInfo contains all parameters to generate a full dipole splitting.
double renormalizationScaleFactor() const
Get the renormalization scale factor.
Ptr< DipoleSplittingKinematics >::tptr splittingKinematics() const
Return the splitting kinematics object.
Energy factorizationScaleFreeze() const
Get the freezing value for the factorization scale.
virtual bool haveOverestimate(const DipoleSplittingInfo &) const
Return true, if this kernel is capable of delivering an overestimate to the kernel, and of inverting the integral over the overestimate w.r.t.
virtual double overestimate(const DipoleSplittingInfo &) const
Return the overestimate to this splitting kernel for the given dipole splitting.
PDPtr flavour() const
Return the flavour produced, if this cannot be determined from the dipole.
constexpr ZeroUnit ZERO
double detuning() const
Return the detuning factor applied to the sampling overestimate kernel.
DipoleSplittingKernel is the base class for all kernels used within the dipole shower.
Energy theScreeningScale
An optional &#39;colour screening&#39; scale for alternative intrinsic pt generation.