Institut de recherche en mathématique et physiqueCentre de Cosmologie, Physique des Particules et Phénoménologie

MonteCarlo Simulation
Olivier Mattelaer

Kind of measurement
Peak

Shape

Rate

Background directly
measured from data.
Theory needed only forparameter extraction

Background SHAPE needed.
Flexible MC for both signal
and background validated andtuned to data

Relies on prediction for bothshape and normalization.
Complicated interplay of best
simulations and data

Mattelaer Olivier

Japan 2023

X

Theory side
Lagrangian

Feynman Rule

• This is Where
the new idea
are expressed• Sameinformation as
the Lagrangian
FeynRules

Mattelaer Olivier

Japan 2023

Cross-section• What is theprecision?
X

Monte-Carlo Physics

Our goal• Cross-section• Differential cross-section• Un-weighted events
Mattelaer Olivier

Japan 2023

2

Simulation of collider events

Simulation of collider events

Mattelaer Olivier

Japan 2023

3

What are the MC for?
ScalesTeV

GeV
Sherpa artist

MeV
Mattelaer Olivier

Japan 2023

4

What are the MC for?
ScalesTeV

GeV

1. High-Q2 Scattering

2. Parton Shower☞ where BSM physics lies

☞ process dependent☞ rst principles description☞ it can be systematically improved

MeV

fi

Mattelaer Olivier

3. Hadronization

4. Underlying EventJapan 2023

5

What are the MC for?
Scales

1. High-Q2 Scattering

2. Parton Shower

TeV

GeV

☞ QCD -”known physics”☞ universal/ process independent☞ rst principles description
MeV

fi

Mattelaer Olivier

3. Hadronization

4. Underlying EventJapan 2023

6

What are the MC for?
Scales

1. High-Q2 Scattering

2. Parton Shower

TeV

GeV

2

☞ low Q physics☞ universal/ process independent☞ model-based description

MeV
Mattelaer Olivier

3. Hadronization

4. Underlying EventJapan 2023

7

What are the MC for?
Scales

1. High-Q2 Scattering

2. Parton Shower

TeV
2☞ low Q physics☞ energy and process dependent☞ model-based description
GeV

MeV
Mattelaer Olivier

3. Hadronization

4. Underlying EventJapan 2023

8

What are the MC for?
Scales

TeV

GeV

MeV
Mattelaer Olivier

Japan 2023

9

To Remember• Multi-scale problem➡ New physics visible only at High scale➡ Problem split in different scale

Factorisation theorem

Mattelaer Olivier

Japan 2023

10

Master formula for the LHC

µF

`+

`

x1 E

x2 E

µF
p

p
long distancelong distancea,bMattelaer Olivier

dx1 dx2 d
Phase-space
integralˆab
FS fa (x1 , µF )fb (x2 , µF ) ⇥Parton density
functionsJapan 2023

X (ŝ, µF , µR )Parton-level crosssection11

Perturbative expansion
dˆ
⇥ab
•

X (ŝ, µF , µR )The parton-level cross section can be computed as aseries in perturbation theory, using the coupling
constant as an expansion parameter, schematically:⇤⇤ˆ = ⇤ Born 1 +LOpredictions•

Mattelaer Olivier

Parton-level cross sections2⇥⇤ (1) +

NLO
corrections

s

2⇥

⇥2

⇤ (2) +

NNLO
corrections

s

2⇥

⇥3

⇤ (3) + . . .

⌅N3LO or NNNLO
corrections

Including higher corrections improves predictions
and reduces theoretical uncertainties
Japan 2023

12

Improved predictions
d =

a,b

⇤
ˆ=⇤

dx1 dx2 fa (x1 , µF )fb (x2 , µF ) dˆ
⇥ab

Born

⇤

1+

s

2⇥

⇤

(1)

+

s

2⇥

⇥2

⇤

(2)

+

s

2⇥

⇥3

X (ŝ, µF , µR )⇤ (3) + . . .

⌅• Leading Order predictions candepend strongly on therenormalization and factorization
scales• Including higher order corrections
reduces the dependence on thesescalesMattelaer Olivier

Japan 2023

13

LO
LO computation (top quark pair)normalized toone

Same variable
At LO:- Large scale uncertainty- but mainly in the Normalisation- LO is good for shape
Mattelaer Olivier

Japan 2023

14

σX =

Going NNLO...?

!" 1
a,b

2
2
Q
Q
2
2
dx1 dx2 fa (x1 , µF )fb (x2 , µF ) × σ̂ab→X (x1 , x2 , αS (µ2R ), 2 , 2 )µF µR0= σ0 +state-of-the-art.
αS σ1 + αS2 σ2 +. . . are only
the current
There• NNLO isσ̂ab→X
a few results available: Higgs, Drell-Yan, ttbar
Why?
• Why do we need it?
● A NNLO computation gives control on the➡ control of the uncertainties in a
uncertainties of a perturbative
calculation.

calculation● It’s “mandatory” if NLO corrections are very large to➡It isof“mandatory”
check the behaviour
the perturbativeifseries
NLO corrections

are very large to check the behavior● It’s the best we have! It is needed for Standard Candlesof the perturbative series
and for really exploiting all the available information, forexample that of NNLO PDF’s.
➡It is needed for Standard Candles

and very precise tests of perturbation theory, exploiting all thee.g. for
CERN Academic Training available
Lectures - Mayinformation,
2012
47 determining NNLO PDF sets Fabio Maltoni
Wednesday 2 May 2012

Let’s focus on LO
Mattelaer Olivier

Japan 2023

X

Tevatron vs. the LHC• Tevatron: 2 TeV proton-antiproton collider➡

Most important: q-q annihilation (85% of t t )• LHC: 7-14 TeV proton-proton collider➡

Most important: g-g annihilation (90% of t t )———Mattelaer Olivier

Japan 2023

X

Hadron Collidersa,bdx1 dx2 dˆab
FS fa (x1 , µF )fb (x2 , µF ) ⇥X (ŝ, µF , µR )Cross sections at a collider depend on:• Coupling strength• Coupling to what?
b(light quarks, gluons, heavy quarks,
EW gauge bosons?)
Mass

W

•
Z
• Single production/pair productiontMattelaer Olivier

Japan 2023

X

Parton densities

n Kinematics
Ratio of Luminosity: LHC at 7 TeV vs Tevatron

.007
5

0.07
5

pdf’s measured in deep-inelastic scattering!

!! Power of collider can be
fully characterized by ratioof parton luminosities!! Ratio larger for gg than qq!! Due to steap rise of gluon
towards low x!! MX=100 GeV!! gg: R"10, e.g. Higgs!! qq: R"3, e.g. W and Z!! MX=800 GeV

ramaticallytowards
lowe.g.
x SUSY
gg: ŝ),
R"1000,At small x!!(small
gluon domination.

large x valence quarks
sections for At
LHC,!! qq:e.g.
R"20, e.g. Z’
Mattelaer Olivier

LHC formidable at large mass –For low mass, Tevatron backgrounds smaller

Japan 2023

11!

X

Hadron colliders
dx1 dx2 dˆab
FS fa (x1 , µF )fb (x2 , µF ) ⇥sity: LHC a,b
at 7 TeV vs Tevatron

X (ŝ, µF , µR )an beby ratioes
han qq
gluon

ggs
and Z
SUSY
Mattelaer Olivier

11!

Japan 2023

X

To Remember
dx1 dx2 da,bPhase-space
integralˆab
FS fa (x1 , µF )fb (x2 , µF ) ⇥Parton density
functions

X (ŝ, µF , µR )Parton-level crosssection• PDF: content of the proton➡ De

ne the physics/processes that will
dominate on your accelerator• LO: good for shape• NLO/NNLO: Reduce scale uncertainty• Computation are inclusive (+ any jet) due
to renormalization/factorization scalefiMattelaer Olivier

Japan 2023

15

Matrix-Element
Calculate a given process (e.g. gluino pair)s s~ > go go WEIGHTED=2

page 1/1• Determine the production mechanism
2

4

2

4

go

s~
s~

go

g

Easy
enough

sl
go

s

s
3

1
diagram 1

go1

3
diagram 2

QCD=2, QED=0

QCD=2, QED=0• Evaluate
the matrix-element!
1|M|2 dΦ(n)σ=➡Need Feynman Rules!
2s• Phase-Space! Integration1|M|2 dΦ(n)
σ=
2s
4

2

2

go

s~

sr

go

1

3

go

1

QCD=2, QED=0

diagram 4

32

QCD=2, QED=04(in general)

Now

go

s~

Mattelaer Oliviers

Hard
Tommorow
Very
Hard

sl

s

diagram 3

3

go

s~

sr

Japan 2023

16

Monte Carlo Integration

Mattelaer Olivier

Japan 2023

17

Matrix-Element
Calculate a given process (e.g. gluino pair)s s~ > go go WEIGHTED=2

page 1/1• Determine the production mechanism
2

4

2

4

go

s~
s~

go

g

Easy
enough

sl
go

s

s
3

1
diagram 1

go1

3
diagram 2

QCD=2, QED=0

QCD=2, QED=0• Evaluate
the matrix-element!
1|M|2 dΦ(n)σ=➡Need Feynman Rules!
2s• Phase-Space! Integration1|M|2 dΦ(n)
σ=
2s
4

2

2

go

s~

sr

go

1

3

go

1

QCD=2, QED=0

diagram 4

32

QCD=2, QED=04(in general)

Now

go

s~

Mattelaer Olivier

s

Hard
Very
Hard

sl

s

diagram 3

3

go

s~

sr

Japan 2023

X

Monte Carlo Integration
Calculations of cross section or decay widths involve
integrations over high-dimension phase space of very
peaked functions:

1
σ=
2s

!

Dim[Φ(n)] ∼ 3n

|M|2 dΦ(n)General and exible method is needed

Mattelaer Olivierfl

Not only integrating but also generates events

Japan 2023

18

Integration
! 1

π
I=
dx cos x
2
0

Z

(q 2

dq 2
M 2 + iM )2

dx C√
IN = 0.637 ± 0.307/ N√
IN = 0.637 ± 0.307/ N

simpson

MC

3

0,638

0,3

5

0,63670,8

200,636620,6

100

0,636619

0,65

1000

0,636619

0,636

Mattelaer Olivier

Z

Method of evaluationp• MonteCarlo• Trapezium• Simpson

Japan 2023

1/ N
1/N 2
1/N 4
19

Integration
! 1

π
I=
dx cos x
2
0

Z

(q 2

dq 2
M 2 + iM )2• MonteCarlo• Trapezium• Simpson

Mattelaer Olivier

1/ N
1/N 2
1/N 4

dx C

√
IN = 0.637 ± 0.307/ N√
IN = 0.637 ± 0.307/ N

Method of evaluationp

Z

More Dimensionp

1/ N
1/N 2/d
1/N 4/dJapan 2023

20

Integration
! 1

Z

π
I=
dx cos x
2
0

(q 2

Z

dq 2
M 2 + iM )2dx C

V = VN = 0√
IN = 0.637 ± 0.307/ N√
IN = 0.637 ± 0.307/ N

I=

! x2
x1

V = (x2 − x1 )1 !
IN = (x2 − x1 )
f (x)
N i=1
N

f (x)dx

! x2
x1

!
2 1
2
VN = (x2 − x1 )[f (x)]2 − IN
N i=1
N[f (x)]2 dx − I 2!
I = IN ± VN /N Can be minimized!
Mattelaer Olivier

Japan 2023

21

Importance Sampling√
IN = 0.637 ± 0.307/ N

! 1

π
I=
dx cos x
2
0

1.0

1.0

0.8

0.8

0.6

0.6

0.4

0.4

0.2

0.2

0.0
0.0

I=

Z 1
0

0.2

0.4

dx(1

0.6

0.8

1.0

0.0
0.0

0.1

0.2

0.3

0.4

0.5

Z 1!
⇡
⇡
π
cos
cos
x
ξ
x[ξ]
cos
2
2 2
2x
2
2
I = = dx(1 dξcx )cx )1−x[ξ](1 2cx2 )(1 cx2 ) 0 ξ1!1√
IN = 0.637 ± 0.307/ N√
IN = 0.637 ± 0.031/ N

fi

The Phase-Space parametrization is important to have an
ef cient computation!
Mattelaer Olivier

Japan 2023

22

Importance Sampling
Z

(q 2

dq 2
M 2 + iM )2⇠ = arctan

✓ 2
q

M2
M

◆
√
IN = 0.637 ± 0.307/ N

Why Importance Sampling?

Probability of usingthat point p(x)The change of variable ensure that the evaluation of
the function is done where the function is the largest!
Mattelaer Olivier

Japan 2023

23

Cut Impact• Events are generated according to our best
knowledge of the function (with importancesampling➡ Adding a cut needs to modi

phase-space integrator

ed the➡ Not possible for custom cut (hardcodedby the user)

Mattelaer Olivier

Run card cut

Japan 2023fi

No cutCustom cut

24

Cut Impact
No cutRun card cut

Custom cut

No cutRun card cut

Custom cut

DAB

Might miss the contribution and think it is just zero.
Mattelaer Olivier

Japan 2023

25

Importance Sampling
Key Point• Generate the random point in a distribution
which is close to the function to integrate.
• This is a change of variable, such that thefunction is atter in this new variable.
• Needs to know an approximate function.

Adaptative Monte-Carlofl

Mattelaer Olivierfl• Create an approximation of the function on
the ight!
Japan 2023

26

Adaptative Monte-Carlo• Create an approximation of the function on
the ight!
Algorithm

1. Creates bin such thateach of them have thesame contribution.
➡Many bins where thefunction is large
2. Use the approximate
for the importancesampling method.

flMattelaer Olivier

Japan 2023

27

VEGAS
More than one Dimension• VEGAS works only with 1(few) dimension➡memory problem
Solution• Use projection on the axis→

p(x)= p(x)•p(y)•p(z)…
• We need toensure thefactorization !
➡ Additionalchange ofvariable

Mattelaer Olivier

Japan 2023

28

choice of the phase-space parametrization has a strong impact on theMonte-Carlo
Integration
efficiency of the MC integration :efficiency of an adaptative MC integration :• The choice of the parameterisation has a
any
peak
is
alignedalong
a single
directionof the P-S
case
1 :strong
any
peak
is
alignedalong
asingledirectionof the P-S
impact on the ef ciency
parametrization
parametrization
y2 y2

y2 y2

Grid

y1 y1

y1 y1

→the
theadaptive
adaptative
Monte-Carlo
P-S
integration
veryefficientefficient→Monte-Carlo
P-S
integration
isisvery
The adaptiveMonte-Carlo
Technique
picks
point
in interesting areas
The technique is efficient

MadWeight – p. 12/29

Japan 2023
fi

Mattelaer Olivier

29

Monte-Carlo Integration

choiceof thephase-spaceparametrization
a strong
impacton thechoiceof thephase-spaceparametrization
has hasa strong
impacton thechoiceof theparametrization
has a strong
impact on
the on the
choiceof phase-space
the phase-spaceparametrization
has a strong
impactefficiencyof theintegration
hoiceof theparametrization
efficiency
ofphase-space
the
MC MC
integration: : has a strong impact on theefficiencyof the:efficiency
ofMC
the integration
MC integration:
ficiency of the MC integration :•theThechoiceof theparametrization
hasofinathe
someare
not casealignedalong
aachange
singledirection
solution
topeaks
previous:: performchangeofvariables
inP-S
order
solutiontotheprevious
caseperformaofvariablesorder
solution to the previous case : perform a change of variables in orderstrong
impactthedirection
ef ciency
some peaks are
not aligned
along ona singleof the P-S
to align
the
peaks
along
a
singledirection
of
the
P-S
parametrization
parametrization
to align
the peaks
alongalong
a single
directionof theP-Sparametrization
to
align
the
peaks
a
singledirection
of
the
P-S
parametrization
parametrization
y2

y2 y2

y2

y2

y2

yy11++y2y2 y2
y1 + y2Rotation
Grid
Grid
y1
y1

y1

y1

y1y1

y1 − y2

y1 − y2y1 −y1 y2→ the adaptive Monte-Carlo P-S integration is very efficient→theadaptive
Monte-CarloP-S integration
converges
slowly→theadaptive
Monte-Carlo
P-S
integrationis very
efficient→ theadaptive
Monte-Carlo
P-S
integrationisvery
efficientTheadaptive
Monte-Carlo
Techniques
picks
point
→ the adaptive Monte-Carlo P-S integration converges
slowlypointsininteresting
everywhereareas
The techniqueintegral converges
is efficient
slowlyMadWeight – p. 12/29
MadWeight – p. 12/29

MadWeight – p. 12/29

MadWeight – p. 12/29
MadWeight – p. 12

Japan 2023
fi

Mattelaer Olivier

30

Multi-channel
What do we do if there isno transformation thataligns all integrand peaks
to the chosen axes?
Vegas is bound to fail!
Solution: use different transformations = channels

p(x) =n!

αi pi (x)

with

i=1

n
!

αi = 1i=1with each pi(x) taking care of one “peak” at the time
Mattelaer Olivier

Japan 2023

X

Multi-channelp(x) =n!

αi pi (x)

i=1

with

n
!

αi = 1i=1

p1(x)Mattelaer Olivier

p2(x)Japan 2023

X

Multi-channelp(x) =n!

αi pi (x)

i=1
n
!

with

αi = 1

i=1

Then,
I=

!

f (x)dx =

n
"
i=1

αi

!

f (x)
pi (x)dxp(x)⇡1Mattelaer Olivier

Japan 2023

XExample: QCD 2 → 2
2

u~

2

4

u~

g

g

u
g

u

u

u u~ > g g QED=0
u u~ > g g QED=0

page 1/1
page 1/1

31
diagram 1

2
2
u~
u~

4
4

2
2

4
4

u~
u~

1

QCD=2

diagram 2

3

2

g
g

u~

g

g
g
g
g

u
u

u

g
g

u
u

u
u
3
3

1
1
diagram 1
diagram 1

1
1

QCD=2
QCD=2

1
1
/ =
ŝ
(p1 + p2 )2
u~
u~

3
3

g
g

u

3
3
diagram 2
diagram 2

2
2

g
g

QCD=2
QCD=2

1
1
/ =
(p1 p3 )2
t̂

g

1

4
diagram 3

QCD=2

1
1
/ =
û
(p1 p4 )2

u
u

u
u

g
g

Three very different pole structures contributing
to the same matrix element.
1
1

4
4

diagram 3
diagram 3

QCD=2
QCD=2

Diagrams made by MadGraph5

Mattelaer Olivier

Japan 2023

31

Multi-channel
Consider the integration of an amplitude |M|^2 at tree level which many
contributing diagrams. We would like to have a basis of functions,

f=

n
!

fi

with

fi ≥ 0 ,

∀ i,

i=1

such that:
1. we know how to integrate each one of them,
2. they describe all possible peaks,
giving us the combined integral

I=

!

! (Φ)
! =
dΦf

n !
"
i=1

n
"
!
f
(
Φ)
i
!
!
dΦ gi (Φ)=
Ii ,!
gi (Φ)
i=1

Does such a basis exist?
Mattelaer Olivier

Japan 2023

X

Single-Diagram-Enhanced
techniqueMulti-channel
based on single
diagrams*

Trick in MadEvent: Split the complexity
Z*Method used in MadGraph

Z P
Z
2
2
X
|M
|
|M
|
i
i
2
2
Pi
P
|Mtot |2 =
|M
|
=
|M
|
tot
tot
2
2
j |Mj |
j |Mj |i
– Any single diagram is “easy” to integrate (polestructures/suitable integration variables known
from the propagators)⇡1– Divide integration into pieces, based on diagrams– All other peaks taken care of by denominator sum

N Integral– Errors add in quadrature so no extra cost
– “Weight” functions already calculated during |M|2 calculation– Parallel in nature

Mattelaer Olivier

Japan 2023

32

Z

Z P
Z
2
2
X
|M
|
|M
|
i
i
2
2
Pi
P
|Mtot |2 =
|M
|
=
|M
|
tot
tot
2
2
j |Mj |
j |Mj |iterm of the above sum.
each term might not begauge invariant

Mattelaer Olivier

Japan 2023

33

To Remember• Phase-Space integration is dif cult
• We need to know the function➡ Be careful with cuts

• MadGraph split the integral in different
contribution linked to the Feynman Diagram➡Those are not the contribution of a given

Mattelaer Olivier

Japan 2023

fi

diagram34

Can we do Better?

• Importance sampling/VEGAS is learning a
Carlo integration
function
HOT TOPIC: Machine Learning➡ has a strong impact on thethe phase-space parametrizationof the MC integration :➡ Lot of work in progressn to the previous case : perform a change of variables in order(Last weeknumber)the peaks along a single direction of the P-S parametrization

multi-channel VEGAS
y1 + y2y2

y1 − y2

y1

adaptive Monte-Carlo P-S integration is very efficient(Preliminary)

W2j
MadWeight – p. 12/29

Variance reduce by a factor 3 (so convergence 9 times faster)
Event generation also three times more ef cientJapan 2023
fi

Mattelaer Olivier

35

Event Generation

Mattelaer Olivier

Japan 2023

36

What is the goal?
• Cross-section• But large theoretical uncertainty• Differential Cross-Section• Provided as sample of events• Sample size is problematic• Those events will need tohave full detector simulation

Mattelaer Olivier

Japan 2023

37

How to get sample?
• Monte-Carlo integration use
random points
#100

#100

• We can keep those• (Uncorrelated) sample• Points not distributed as thereal function• Need to keep track of theimportance of each point(weight)dσdO

OMattelaer Olivier• Typically a lot of event have
low information
Japan 2023

38

Do we need to keep small weight?
• Let’s put a minimum• Discard events below theminimum• NO! We loose cross-section/ bias
ourself
#100

#100
#0
threshold• Let’s put a minimum# 50• But keep 50% of the events below

2x• Multiply the weight of each eventby 2 (preserve cross-section)• We loose information

Mattelaer Olivier

fi• But we gain in le size

Japan 2023

39

Do we need to keep small weight?
• Let’s put a threshold• But keep 50% of the events below• Multiply the weight of each eventby 2 (preserve cross-section)• Let’s improve#100

Threshold

#8
# 50
2x
wmin

• Let’s make the thresholdproportional to the weight
100w%• Keep each event with w
thres
probability

wthres• If kept multiply his weight by w
• So the new weight is wthres
Mattelaer Olivier

Japan 2023

40

Unweighted events
Events distributed as in nature• All bins should event event
proportional to their cross-section(Up to Poisson distribution)• All events should have the same
weight• This correspond to the smallest lesize or maximum compression
dσ
dO

O

Japan 2023
fi

Mattelaer Olivier

41

Do we need to keep small weight?
• Let’s improve• Let’s make the thresholdproportional to the weight• So the new weight is wthres• Let’s all event have the same weight#100

threshold• So set wthres > max(w)• Maximal compression

Mattelaer Olivier

Japan 2023

42

Event generation
1 N
1 N f(xi)f(x)dx =f(xi) =
wthres
∑∑∫N i=1
N i=1 wthres
1. pick xi
f(x)2. calculate f(xi)
3. pick y ∈ [0, max( f )]4. Compare:
if y < f(xi) accept event,else reject it.

Mattelaer Olivier

Japan 2023

43

Event generation
Improved by combining with importance sampling:
1. pick x distributed as p(x)f(x)

2. calculate f(x) and p(x)3. pick 0<y<1
4. Compare:if f(x)>y p(x) accept event,else reject it.

Mattelaer Olivier

fimuch better ef ciency!!!
Japan 2023

44

Monte-Carlo Summary
Bad Point• Slow Convergence (especially in low number
of Dimension)• Need to know the function• Impact on cut
Good Point• Complex area of Integration• Easy error estimate• quick estimation of the integral• Possibility to have unweighted events
Mattelaer Olivier

Japan 2023

45