herwig
is hosted by
Hepforge
,
IPPP Durham
Herwig
7.2.1
Shower
Dipole
Kernels
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
27
namespace
Herwig
{
28
29
using namespace
ThePEG
;
30
41
class
DipoleSplittingKernel
:
public
HandlerBase
{
42
43
public
:
44
50
DipoleSplittingKernel
();
51
55
virtual
~
DipoleSplittingKernel
();
57
58
public
:
59
63
Ptr<AlphaSBase>::tptr
alphaS
()
const
{
return
theAlphaS; }
64
68
void
alphaS
(
Ptr<AlphaSBase>::tptr
ap) { theAlphaS = ap; }
69
73
Ptr<DipoleSplittingKinematics>::tptr
splittingKinematics
()
const
{
74
return
theSplittingKinematics;
75
}
76
80
Ptr<DipoleMCCheck>::ptr
mcCheck
()
const
{
return
theMCCheck; }
81
85
void
splittingKinematics
(Ptr<DipoleSplittingKinematics>::tptr sp) {
86
theSplittingKinematics = sp;
87
}
88
92
Ptr<PDFRatio>::tptr
pdfRatio
()
const
{
return
thePDFRatio; }
93
97
void
pdfRatio
(Ptr<PDFRatio>::tptr sp) { thePDFRatio = sp; }
98
107
virtual
int
nDimAdditional
()
const
{
return
0; }
108
112
void
renormalizationScaleFreeze
(Energy s) { theRenormalizationScaleFreeze = s; }
113
117
void
factorizationScaleFreeze
(Energy s) { theFactorizationScaleFreeze = s; }
118
122
Energy
renormalizationScaleFreeze
()
const
{
return
theRenormalizationScaleFreeze; }
123
127
Energy
factorizationScaleFreeze
()
const
{
return
theFactorizationScaleFreeze; }
128
129
public
:
130
135
virtual
bool
canHandle(
const
DipoleIndex
&)
const
= 0;
136
142
virtual
bool
canHandleEquivalent(
const
DipoleIndex
& a,
143
const
DipoleSplittingKernel
& sk,
144
const
DipoleIndex
& b)
const
= 0;
145
150
virtual
tcPDPtr
emitter(
const
DipoleIndex
&)
const
= 0;
151
156
virtual
tcPDPtr
emission(
const
DipoleIndex
&)
const
= 0;
157
162
virtual
tcPDPtr
spectator(
const
DipoleIndex
&)
const
= 0;
163
168
PDPtr
flavour
()
const
{
return
theFlavour; }
169
174
bool
strictLargeN
()
const
{
return
theStrictLargeN; }
175
176
public
:
177
182
virtual
void
startPresampling
(
const
DipoleIndex
&) {
183
presampling =
true
;
184
}
185
190
virtual
void
stopPresampling
(
const
DipoleIndex
&) {
191
presampling =
false
;
192
}
193
198
unsigned
long
presamplingPoints
()
const
{
return
thePresamplingPoints; }
199
204
unsigned
long
maxtry
()
const
{
return
theMaxtry; }
205
210
unsigned
long
freezeGrid
()
const
{
return
theFreezeGrid; }
211
216
void
freezeGrid
(
unsigned
long
n) { theFreezeGrid = n; }
217
221
void
detuning
(
double
d) { theDetuning = d; }
222
226
double
detuning
()
const
{
return
theDetuning; }
227
232
virtual
double
evaluate(
const
DipoleSplittingInfo
&)
const
= 0;
233
238
virtual
vector< pair<int, Complex> > generatePhi(
const
DipoleSplittingInfo
& dInfo,
const
RhoDMatrix
& rho)
const
= 0;
239
243
virtual
DecayMEPtr matrixElement(
const
DipoleSplittingInfo
& dInfo)
const
= 0;
244
248
void
clearAlphaPDFCache
()
const
{
249
theAlphaSCache.clear();
250
thePDFCache.clear();
251
}
252
257
virtual
void
accept(
const
DipoleSplittingInfo
&,
double
,
double
, map<string,double>&)
const
;
258
263
virtual
void
veto(
const
DipoleSplittingInfo
&,
double
,
double
, map<string,double>&)
const
;
264
272
virtual
bool
haveOverestimate
(
const
DipoleSplittingInfo
&)
const
{
return
false
; }
273
278
virtual
double
overestimate
(
const
DipoleSplittingInfo
&)
const
{
return
-1.; }
279
286
virtual
double
invertOverestimateIntegral
(
const
DipoleSplittingInfo
&,
double
)
const
{
287
return
-1.;
288
}
289
293
bool
useThisKernel
()
const
{
294
return
theUseThisKernel;
295
}
296
297
public
:
298
302
double
factorizationScaleFactor
()
const
{
return
theFactorizationScaleFactor; }
303
307
void
factorizationScaleFactor
(
double
f) { theFactorizationScaleFactor = f; }
308
312
double
renormalizationScaleFactor
()
const
{
return
theRenormalizationScaleFactor; }
313
317
void
renormalizationScaleFactor
(
double
f) { theRenormalizationScaleFactor = f; }
318
319
protected
:
320
324
double
alphaPDF(
const
DipoleSplittingInfo
&,
325
Energy optScale =
ZERO
,
326
double
rScaleFactor = 1.0,
327
double
fScaleFactor = 1.0)
const
;
328
333
bool
virtualitySplittingScale
()
const
{
return
theVirtualitySplittingScale; }
334
335
public
:
336
343
void
persistentOutput(
PersistentOStream
& os)
const
;
344
350
void
persistentInput(
PersistentIStream
& is,
int
version);
352
359
static
void
Init();
360
361
362
// If needed, insert declarations of virtual function defined in the
363
// InterfacedBase class here (using ThePEG-interfaced-decl in Emacs).
364
365
private
:
366
370
Ptr<AlphaSBase>::ptr
theAlphaS
;
371
376
Energy
theScreeningScale
;
377
381
Ptr<DipoleSplittingKinematics>::ptr
theSplittingKinematics
;
382
387
Ptr<PDFRatio>::ptr
thePDFRatio
;
388
393
unsigned
long
thePresamplingPoints
;
394
399
unsigned
long
theMaxtry
;
400
405
static
double
theMaxPDFRatio
;
406
411
unsigned
long
theFreezeGrid
;
412
416
double
theDetuning
;
417
422
PDPtr
theFlavour
;
423
427
Ptr<DipoleMCCheck>::ptr
theMCCheck
;
428
433
bool
theStrictLargeN
;
434
438
double
theFactorizationScaleFactor
;
439
443
double
theRenormalizationScaleFactor
;
444
448
Energy
theRenormalizationScaleFreeze
;
449
453
Energy
theFactorizationScaleFreeze
;
454
459
bool
theVirtualitySplittingScale
;
460
465
unsigned
int
theCMWScheme=0;
466
470
mutable
map<double,double>
theAlphaSCache
;
471
475
mutable
map<double,double>
thePDFCache
;
476
480
bool
presampling
;
481
485
bool
theUseThisKernel =
true
;
486
487
488
private
:
489
494
static
AbstractClassDescription<DipoleSplittingKernel>
initDipoleSplittingKernel
;
495
500
DipoleSplittingKernel
& operator=(
const
DipoleSplittingKernel
&) =
delete
;
501
502
};
503
504
}
505
506
#include "ThePEG/Utilities/ClassTraits.h"
507
508
namespace
ThePEG
{
509
514
template
<>
515
struct
BaseClassTrait
<Herwig::
DipoleSplittingKernel
,1> {
517
typedef
HandlerBase
NthBase;
518
};
519
522
template
<>
523
struct
ClassTraits
<Herwig::
DipoleSplittingKernel
>
524
:
public
ClassTraitsBase
<Herwig::DipoleSplittingKernel> {
526
static
string
className() {
return
"Herwig::DipoleSplittingKernel"
; }
534
static
string
library() {
return
"HwDipoleShower.so"
; }
535
};
536
539
}
540
541
#endif
/* HERWIG_DipoleSplittingKernel_H */
ThePEG::PersistentIStream
Herwig::DipoleSplittingKernel::freezeGrid
unsigned long freezeGrid() const
Return the number of accepted points after which the grid should be frozen.
Definition:
DipoleSplittingKernel.h:210
Herwig::DipoleSplittingKernel::mcCheck
Ptr< DipoleMCCheck >::ptr mcCheck() const
Return the mc check object.
Definition:
DipoleSplittingKernel.h:80
Herwig::DipoleSplittingKernel::theFactorizationScaleFactor
double theFactorizationScaleFactor
The factorization scale factor.
Definition:
DipoleSplittingKernel.h:438
Ptr< AlphaSBase >::tptr
transient_pointer tptr
Herwig::DipoleSplittingKernel::renormalizationScaleFreeze
Energy renormalizationScaleFreeze() const
Get the freezing value for the renormalization scale.
Definition:
DipoleSplittingKernel.h:122
Herwig::DipoleSplittingKernel::initDipoleSplittingKernel
static AbstractClassDescription< DipoleSplittingKernel > initDipoleSplittingKernel
The static object used to initialize the description of this class.
Definition:
DipoleSplittingKernel.h:494
Herwig::DipoleSplittingKernel::strictLargeN
bool strictLargeN() const
Return true, if this splitting kernel is supposed to work in a strict large-N limit, i.e.
Definition:
DipoleSplittingKernel.h:174
Herwig::DipoleSplittingKernel::theSplittingKinematics
Ptr< DipoleSplittingKinematics >::ptr theSplittingKinematics
The splitting kinematics to be used.
Definition:
DipoleSplittingKernel.h:381
Herwig::DipoleSplittingKernel::freezeGrid
void freezeGrid(unsigned long n)
Set the number of accepted points after which the grid should be frozen.
Definition:
DipoleSplittingKernel.h:216
ThePEG::PersistentOStream
ThePEG::AbstractClassDescription
Herwig::DipoleSplittingKernel::stopPresampling
virtual void stopPresampling(const DipoleIndex &)
Inform this splitting kernel, that it is not being presampled until a call to startPresampling.
Definition:
DipoleSplittingKernel.h:190
Herwig::DipoleSplittingKernel::factorizationScaleFactor
void factorizationScaleFactor(double f)
Set the factorization scale factor.
Definition:
DipoleSplittingKernel.h:307
Herwig::DipoleSplittingKernel::theStrictLargeN
bool theStrictLargeN
True, if this splitting kernel is supposed to work in a strict large-N limit, i.e.
Definition:
DipoleSplittingKernel.h:433
Herwig::DipoleIndex
DipoleIndex is used to index splitting generators for a particular dipole.
Definition:
DipoleSplittingInfo.h:35
Herwig::DipoleSplittingKernel::splittingKinematics
void splittingKinematics(Ptr< DipoleSplittingKinematics >::tptr sp)
Set the splitting kinematics object.
Definition:
DipoleSplittingKernel.h:85
Herwig::DipoleSplittingKernel::virtualitySplittingScale
bool virtualitySplittingScale() const
Return true, if the virtuality of the splitting should be used as the argument of alphas rather than ...
Definition:
DipoleSplittingKernel.h:333
Herwig::DipoleSplittingKernel::theFreezeGrid
unsigned long theFreezeGrid
Return the number of accepted points after which the grid should be frozen.
Definition:
DipoleSplittingKernel.h:411
ThePEG
Herwig::DipoleSplittingKernel::theMaxPDFRatio
static double theMaxPDFRatio
The maximum value for any pdf ratio.
Definition:
DipoleSplittingKernel.h:405
Herwig::DipoleSplittingKernel::theMCCheck
Ptr< DipoleMCCheck >::ptr theMCCheck
Pointer to a check histogram object.
Definition:
DipoleSplittingKernel.h:427
Herwig::DipoleSplittingKernel::theAlphaS
Ptr< AlphaSBase >::ptr theAlphaS
The alpha_s to be used.
Definition:
DipoleSplittingKernel.h:370
Herwig::DipoleSplittingKernel::thePresamplingPoints
unsigned long thePresamplingPoints
The number of points to presample this splitting generator.
Definition:
DipoleSplittingKernel.h:393
Herwig::DipoleSplittingKernel::alphaS
void alphaS(Ptr< AlphaSBase >::tptr ap)
Set the alpha_s to be used.
Definition:
DipoleSplittingKernel.h:68
Herwig::DipoleSplittingKernel::theRenormalizationScaleFreeze
Energy theRenormalizationScaleFreeze
A freezing value for the renormalization scale.
Definition:
DipoleSplittingKernel.h:448
Herwig::DipoleSplittingKernel::theDetuning
double theDetuning
The detuning factor applied to the sampling overestimate kernel.
Definition:
DipoleSplittingKernel.h:416
Herwig::DipoleSplittingKernel::maxtry
unsigned long maxtry() const
Return the maximum number of trials to generate a splitting.
Definition:
DipoleSplittingKernel.h:204
Herwig::DipoleSplittingKernel::invertOverestimateIntegral
virtual double invertOverestimateIntegral(const DipoleSplittingInfo &, double) const
Invert the integral over the overestimate w.r.t.
Definition:
DipoleSplittingKernel.h:286
Herwig::DipoleSplittingKernel::theAlphaSCache
map< double, double > theAlphaSCache
Cache for alphas evaluations.
Definition:
DipoleSplittingKernel.h:470
Herwig::DipoleSplittingKernel::theVirtualitySplittingScale
bool theVirtualitySplittingScale
True, if the virtuality of the splitting should be used as the argument of alphas rather than the pt...
Definition:
DipoleSplittingKernel.h:459
Herwig::DipoleSplittingKernel::pdfRatio
Ptr< PDFRatio >::tptr pdfRatio() const
Return the PDFRatio object.
Definition:
DipoleSplittingKernel.h:92
Herwig::DipoleSplittingKernel::renormalizationScaleFreeze
void renormalizationScaleFreeze(Energy s)
Set the freezing value for the renormalization scale.
Definition:
DipoleSplittingKernel.h:112
Herwig::DipoleSplittingKernel::thePDFRatio
Ptr< PDFRatio >::ptr thePDFRatio
An optional PDF ratio object to be used when evaluating this kernel.
Definition:
DipoleSplittingKernel.h:387
ThePEG::tcPDPtr
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr
Herwig::DipoleSplittingKernel::presampling
bool presampling
True, if we are presampling.
Definition:
DipoleSplittingKernel.h:480
ThePEG::PDPtr
ThePEG::Ptr< ParticleData >::pointer PDPtr
Herwig::DipoleSplittingKernel::alphaS
Ptr< AlphaSBase >::tptr alphaS() const
Return the alpha_s to be used.
Definition:
DipoleSplittingKernel.h:63
Herwig::DipoleSplittingKernel::detuning
void detuning(double d)
Set a detuning factor to be applied to the sampling overestimate kernel.
Definition:
DipoleSplittingKernel.h:221
Herwig::DipoleSplittingKernel::thePDFCache
map< double, double > thePDFCache
Cache for PDF evaluations.
Definition:
DipoleSplittingKernel.h:475
Herwig::DipoleSplittingKernel::clearAlphaPDFCache
void clearAlphaPDFCache() const
Clear the alphaPDF cache.
Definition:
DipoleSplittingKernel.h:248
Herwig::DipoleSplittingKernel::presamplingPoints
unsigned long presamplingPoints() const
Return the number of points to presample this splitting generator.
Definition:
DipoleSplittingKernel.h:198
ThePEG::HandlerBase
Herwig::DipoleSplittingKernel::factorizationScaleFreeze
void factorizationScaleFreeze(Energy s)
Set the freezing value for the factorization scale.
Definition:
DipoleSplittingKernel.h:117
Herwig::DipoleSplittingKernel::nDimAdditional
virtual int nDimAdditional() const
Return the number of additional parameter random variables needed to evaluate this kernel except the ...
Definition:
DipoleSplittingKernel.h:107
Herwig::DipoleSplittingKernel::pdfRatio
void pdfRatio(Ptr< PDFRatio >::tptr sp)
Set the PDFRatio object.
Definition:
DipoleSplittingKernel.h:97
Herwig::DipoleSplittingKernel::theFactorizationScaleFreeze
Energy theFactorizationScaleFreeze
A freezing value for the factorization scale.
Definition:
DipoleSplittingKernel.h:453
Herwig::DipoleSplittingKernel::theFlavour
PDPtr theFlavour
The flavour produced, if this cannot be determined from the dipole.
Definition:
DipoleSplittingKernel.h:422
Herwig::DipoleSplittingKernel::startPresampling
virtual void startPresampling(const DipoleIndex &)
Inform this splitting kernel, that it is being presampled until a call to stopPresampling.
Definition:
DipoleSplittingKernel.h:182
Herwig::DipoleSplittingKernel::theMaxtry
unsigned long theMaxtry
The maximum number of trials to generate a splitting.
Definition:
DipoleSplittingKernel.h:399
Herwig
-*- C++ -*-
Definition:
BasicConsistency.h:17
Herwig::DipoleSplittingKernel::renormalizationScaleFactor
void renormalizationScaleFactor(double f)
Set the renormalization scale factor.
Definition:
DipoleSplittingKernel.h:317
ThePEG::ClassTraits
Herwig::DipoleSplittingKernel::theRenormalizationScaleFactor
double theRenormalizationScaleFactor
The renormalization scale factor.
Definition:
DipoleSplittingKernel.h:443
Herwig::DipoleSplittingKernel::factorizationScaleFactor
double factorizationScaleFactor() const
Get the factorization scale factor.
Definition:
DipoleSplittingKernel.h:302
Herwig::DipoleSplittingInfo
DipoleSplittingInfo contains all parameters to generate a full dipole splitting.
Definition:
DipoleSplittingInfo.h:209
Herwig::DipoleSplittingKernel::renormalizationScaleFactor
double renormalizationScaleFactor() const
Get the renormalization scale factor.
Definition:
DipoleSplittingKernel.h:312
Herwig::DipoleSplittingKernel::splittingKinematics
Ptr< DipoleSplittingKinematics >::tptr splittingKinematics() const
Return the splitting kinematics object.
Definition:
DipoleSplittingKernel.h:73
ThePEG::RhoDMatrix
Herwig::DipoleSplittingKernel::factorizationScaleFreeze
Energy factorizationScaleFreeze() const
Get the freezing value for the factorization scale.
Definition:
DipoleSplittingKernel.h:127
Herwig::DipoleSplittingKernel::haveOverestimate
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.
Definition:
DipoleSplittingKernel.h:272
Herwig::DipoleSplittingKernel::overestimate
virtual double overestimate(const DipoleSplittingInfo &) const
Return the overestimate to this splitting kernel for the given dipole splitting.
Definition:
DipoleSplittingKernel.h:278
ThePEG::BaseClassTrait
Herwig::DipoleSplittingKernel::flavour
PDPtr flavour() const
Return the flavour produced, if this cannot be determined from the dipole.
Definition:
DipoleSplittingKernel.h:168
ZERO
constexpr ZeroUnit ZERO
Herwig::DipoleSplittingKernel::detuning
double detuning() const
Return the detuning factor applied to the sampling overestimate kernel.
Definition:
DipoleSplittingKernel.h:226
Herwig::DipoleSplittingKernel::useThisKernel
bool useThisKernel() const
Definition:
DipoleSplittingKernel.h:293
Herwig::DipoleSplittingKernel
DipoleSplittingKernel is the base class for all kernels used within the dipole shower.
Definition:
DipoleSplittingKernel.h:41
Ptr< AlphaSBase >::ptr
pointer ptr
ThePEG::ClassTraitsBase
Herwig::DipoleSplittingKernel::theScreeningScale
Energy theScreeningScale
An optional 'colour screening' scale for alternative intrinsic pt generation.
Definition:
DipoleSplittingKernel.h:376
Generated on Sat Apr 11 2020 14:50:29 for Herwig by
1.8.13