herwig is hosted by Hepforge, IPPP Durham
Herwig 7.3.0
DipoleSplittingKernel.h
1// -*- C++ -*-
2//
3// DipoleSplittingKernel.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_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#include "ThePEG/EventRecord/RhoDMatrix.h"
24#include "Herwig/Decay/DecayMatrixElement.h"
25#include "Herwig/Decay/TwoBodyDecayMatrixElement.h"
26
27namespace Herwig {
28
29using namespace ThePEG;
30
42
43public:
44
49
50public:
51
55 Ptr<AlphaSBase>::tptr alphaS() const { return theAlphaS; }
56
60 void alphaS(Ptr<AlphaSBase>::tptr ap) { theAlphaS = ap; }
61
65 Ptr<DipoleSplittingKinematics>::tptr splittingKinematics() const {
67 }
68
72 Ptr<DipoleMCCheck>::ptr mcCheck() const { return theMCCheck; }
73
77 void splittingKinematics(Ptr<DipoleSplittingKinematics>::tptr sp) {
79 }
80
84 Ptr<PDFRatio>::tptr pdfRatio() const { return thePDFRatio; }
85
89 void pdfRatio(Ptr<PDFRatio>::tptr sp) { thePDFRatio = sp; }
90
99 virtual int nDimAdditional() const { return 0; }
100
105
110
115
120
121public:
122
127 virtual bool canHandle(const DipoleIndex&) const = 0;
128
134 virtual bool canHandleEquivalent(const DipoleIndex& a,
135 const DipoleSplittingKernel& sk,
136 const DipoleIndex& b) const = 0;
137
142 virtual tcPDPtr emitter(const DipoleIndex&) const = 0;
143
148 virtual tcPDPtr emission(const DipoleIndex&) const = 0;
149
154 virtual tcPDPtr spectator(const DipoleIndex&) const = 0;
155
160 PDPtr flavour() const { return theFlavour; }
161
166 bool strictLargeN() const { return theStrictLargeN; }
167
168public:
169
174 virtual void startPresampling(const DipoleIndex&) {
175 presampling = true;
176 }
177
182 virtual void stopPresampling(const DipoleIndex&) {
183 presampling = false;
184 }
185
190 unsigned long presamplingPoints() const { return thePresamplingPoints; }
191
196 unsigned long maxtry() const { return theMaxtry; }
197
202 unsigned long freezeGrid() const { return theFreezeGrid; }
203
208 void freezeGrid(unsigned long n) { theFreezeGrid = n; }
209
213 void detuning(double d) { theDetuning = d; }
214
218 double detuning() const { return theDetuning; }
219
224 virtual double evaluate(const DipoleSplittingInfo&) const = 0;
225
230 virtual vector< pair<int, Complex> > generatePhi( const DipoleSplittingInfo& dInfo, const RhoDMatrix& rho) const = 0;
231
235 virtual DecayMEPtr matrixElement(const DipoleSplittingInfo& dInfo) const = 0;
236
240 void clearAlphaPDFCache() const {
241 theAlphaSCache.clear();
242 thePDFCache.clear();
243 }
244
249 virtual void accept(const DipoleSplittingInfo&, double, double, map<string,double>&) const;
250
255 virtual void veto(const DipoleSplittingInfo&, double, double, map<string,double>&) const;
256
264 virtual bool haveOverestimate(const DipoleSplittingInfo&) const { return false; }
265
270 virtual double overestimate(const DipoleSplittingInfo&) const { return -1.; }
271
278 virtual double invertOverestimateIntegral(const DipoleSplittingInfo&, double) const {
279 return -1.;
280 }
281
285 bool useThisKernel() const {
286 return theUseThisKernel;
287 }
288
289public:
290
295
300
305
310
314 unsigned int cmwScheme() const {return theCMWScheme;}
315
316protected:
317
322 Energy optScale = ZERO,
323 double rScaleFactor = 1.0,
324 double fScaleFactor = 1.0) const;
325
331
332public:
333
341
347 void persistentInput(PersistentIStream & is, int version);
349
356 static void Init();
357
358
359// If needed, insert declarations of virtual function defined in the
360// InterfacedBase class here (using ThePEG-interfaced-decl in Emacs).
361
362private:
363
367 Ptr<AlphaSBase>::ptr theAlphaS;
368
374
378 Ptr<DipoleSplittingKinematics>::ptr theSplittingKinematics;
379
384 Ptr<PDFRatio>::ptr thePDFRatio;
385
390 unsigned long thePresamplingPoints;
391
396 unsigned long theMaxtry;
397
402 static double theMaxPDFRatio;
403
408 unsigned long theFreezeGrid;
409
414
420
424 Ptr<DipoleMCCheck>::ptr theMCCheck;
425
431
436
441
446
451
457
462 unsigned int theCMWScheme=0;
463
467 mutable map<double,double> theAlphaSCache;
468
472 mutable map<double,double> thePDFCache;
473
478
482 bool theUseThisKernel = true;
483
484
485private:
486
492
498
499};
500
501}
502
503#include "ThePEG/Utilities/ClassTraits.h"
504
505namespace ThePEG {
506
511template <>
512struct BaseClassTrait<Herwig::DipoleSplittingKernel,1> {
514 typedef HandlerBase NthBase;
515};
516
519template <>
520struct ClassTraits<Herwig::DipoleSplittingKernel>
521 : public ClassTraitsBase<Herwig::DipoleSplittingKernel> {
523 static string className() { return "Herwig::DipoleSplittingKernel"; }
531 static string library() { return "HwDipoleShower.so"; }
532};
533
536}
537
538#endif /* HERWIG_DipoleSplittingKernel_H */
DipoleIndex is used to index splitting generators for a particular dipole.
DipoleSplittingInfo contains all parameters to generate a full dipole splitting.
DipoleSplittingKernel is the base class for all kernels used within the dipole shower.
bool strictLargeN() const
Return true, if this splitting kernel is supposed to work in a strict large-N limit,...
void persistentInput(PersistentIStream &is, int version)
Function used to read in object persistently.
void renormalizationScaleFactor(double f)
Set the renormalization scale factor.
virtual double overestimate(const DipoleSplittingInfo &) const
Return the overestimate to this splitting kernel for the given dipole splitting.
Ptr< AlphaSBase >::tptr alphaS() const
Return the alpha_s to be used.
Energy factorizationScaleFreeze() const
Get the freezing value for the factorization scale.
Ptr< PDFRatio >::ptr thePDFRatio
An optional PDF ratio object to be used when evaluating this kernel.
double theRenormalizationScaleFactor
The renormalization scale factor.
Energy theScreeningScale
An optional 'colour screening' scale for alternative intrinsic pt generation.
virtual tcPDPtr spectator(const DipoleIndex &) const =0
Return the spectator data after splitting, given a dipole index.
Ptr< DipoleMCCheck >::ptr theMCCheck
Pointer to a check histogram object.
bool presampling
True, if we are presampling.
PDPtr theFlavour
The flavour produced, if this cannot be determined from the dipole.
PDPtr flavour() const
Return the flavour produced, if this cannot be determined from the dipole.
map< double, double > thePDFCache
Cache for PDF evaluations.
unsigned long freezeGrid() const
Return the number of accepted points after which the grid should be frozen.
virtual bool canHandleEquivalent(const DipoleIndex &a, const DipoleSplittingKernel &sk, const DipoleIndex &b) const =0
Return true, if this splitting kernel is the same for the given index a, as the given splitting kerne...
unsigned long thePresamplingPoints
The number of points to presample this splitting generator.
Energy theFactorizationScaleFreeze
A freezing value for the factorization scale.
virtual double evaluate(const DipoleSplittingInfo &) const =0
Evaluate this splitting kernel for the given dipole splitting.
Ptr< DipoleSplittingKinematics >::tptr splittingKinematics() const
Return the splitting kinematics object.
void freezeGrid(unsigned long n)
Set the number of accepted points after which the grid should be frozen.
double detuning() const
Return the detuning factor applied to the sampling overestimate kernel.
void renormalizationScaleFreeze(Energy s)
Set the freezing value for the renormalization scale.
void detuning(double d)
Set a detuning factor to be applied to the sampling overestimate kernel.
virtual void accept(const DipoleSplittingInfo &, double, double, map< string, double > &) const
Update the variations vector at the given splitting using the indicated kernel and overestimate value...
void pdfRatio(Ptr< PDFRatio >::tptr sp)
Set the PDFRatio object.
unsigned long theMaxtry
The maximum number of trials to generate a splitting.
double theFactorizationScaleFactor
The factorization scale factor.
void splittingKinematics(Ptr< DipoleSplittingKinematics >::tptr sp)
Set the splitting kinematics object.
Ptr< DipoleSplittingKinematics >::ptr theSplittingKinematics
The splitting kinematics to be used.
unsigned long maxtry() const
Return the maximum number of trials to generate a splitting.
unsigned int theCMWScheme
Implementing CMW in the kernels.
unsigned int cmwScheme() const
Return the CMW sheme used by the kernel.
virtual int nDimAdditional() const
Return the number of additional parameter random variables needed to evaluate this kernel except the ...
double renormalizationScaleFactor() const
Get the renormalization scale factor.
void alphaS(Ptr< AlphaSBase >::tptr ap)
Set the alpha_s to be used.
virtual void startPresampling(const DipoleIndex &)
Inform this splitting kernel, that it is being presampled until a call to stopPresampling.
Energy renormalizationScaleFreeze() const
Get the freezing value for the renormalization scale.
virtual vector< pair< int, Complex > > generatePhi(const DipoleSplittingInfo &dInfo, const RhoDMatrix &rho) const =0
Evaluate rho_ii' V_ijk V*_i'jk / equivalent for initial-state splitting, required for generating spin...
void factorizationScaleFreeze(Energy s)
Set the freezing value for the factorization scale.
static AbstractClassDescription< DipoleSplittingKernel > initDipoleSplittingKernel
The static object used to initialize the description of this class.
virtual void stopPresampling(const DipoleIndex &)
Inform this splitting kernel, that it is not being presampled until a call to startPresampling.
virtual tcPDPtr emitter(const DipoleIndex &) const =0
Return the emitter data after splitting, given a dipole index.
void persistentOutput(PersistentOStream &os) const
Function used to write out object persistently.
static double theMaxPDFRatio
The maximum value for any pdf ratio.
void factorizationScaleFactor(double f)
Set the factorization scale factor.
virtual bool haveOverestimate(const DipoleSplittingInfo &) const
Return true, if this kernel is capable of delivering an overestimate to the kernel,...
map< double, double > theAlphaSCache
Cache for alphas evaluations.
virtual double invertOverestimateIntegral(const DipoleSplittingInfo &, double) const
Invert the integral over the overestimate w.r.t.
virtual bool canHandle(const DipoleIndex &) const =0
Return true, if this splitting kernel applies to the given dipole index.
void clearAlphaPDFCache() const
Clear the alphaPDF cache.
Ptr< AlphaSBase >::ptr theAlphaS
The alpha_s to be used.
DipoleSplittingKernel & operator=(const DipoleSplittingKernel &)=delete
The assignment operator is private and must never be called.
DipoleSplittingKernel()
The default constructor.
bool virtualitySplittingScale() const
Return true, if the virtuality of the splitting should be used as the argument of alphas rather than ...
Ptr< DipoleMCCheck >::ptr mcCheck() const
Return the mc check object.
double alphaPDF(const DipoleSplittingInfo &, Energy optScale=ZERO, double rScaleFactor=1.0, double fScaleFactor=1.0) const
Return the common factor of (alphas/2pi)*(pdf ratio)
double factorizationScaleFactor() const
Get the factorization scale factor.
virtual tcPDPtr emission(const DipoleIndex &) const =0
Return the emission data after splitting, given a dipole index.
unsigned long presamplingPoints() const
Return the number of points to presample this splitting generator.
Ptr< PDFRatio >::tptr pdfRatio() const
Return the PDFRatio object.
unsigned long theFreezeGrid
Return the number of accepted points after which the grid should be frozen.
virtual void veto(const DipoleSplittingInfo &, double, double, map< string, double > &) const
Update the variations vector at the given splitting using the indicated kernel and overestimate value...
double theDetuning
The detuning factor applied to the sampling overestimate kernel.
bool theUseThisKernel
True, if the kernel should be used.
bool theVirtualitySplittingScale
True, if the virtuality of the splitting should be used as the argument of alphas rather than the pt.
static void Init()
The standard Init function used to initialize the interfaces.
virtual DecayMEPtr matrixElement(const DipoleSplittingInfo &dInfo) const =0
Return the completely spin-unaveraged (i.e.
bool theStrictLargeN
True, if this splitting kernel is supposed to work in a strict large-N limit, i.e.
Energy theRenormalizationScaleFreeze
A freezing value for the renormalization scale.
-*- C++ -*-
constexpr ZeroUnit ZERO
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr
ThePEG::Ptr< ParticleData >::pointer PDPtr
static string className()
static string library()