MG5 generates too tiny cross-section

Asked by Samadrita Mukherjee

Hello,

We have written down our model file where a new scalar has interactions with SM quarks (both left and right handed.) In the Feynrules the interactions are written as,

    1/fa (del[ALP,mu]) (QLbar.Ga[mu].QL + uRbar.Ga[mu].uR + dRbar.Ga[mu].dR)

UFO files are generated successfully and other sanity checks like the hermiticity of the Lagrangian etc are also done.

When I import the UFO model in MG5 and try to generate the cross-section for a process

generate e p > e j ALP NP=1

it is giving a absurd small cross-section of the order of 10^{-31} pb.

However, with the interaction with

either QLbar.Ga[mu].QL or uRbar.Ga[mu].uR + dRbar.Ga[mu].dR

the cross-section is 10^{-7} pb which we expect it to be for the parameters we are running.

My question what is wrong with the cross-section in MG5, when we are just adding the left-handed doublet and the right-handed quarks together, why it is cancelling out something so strongly yielding a tiny cross-section?

P.S. I did the same exercise with QLbar.Ga[mu].QL - uRbar.Ga[mu].uR - dRbar.Ga[mu].dR, i.e. subtracting the right-handed quarks (it is wrong physics-wise the Axion interaction), and I got again 10^{-7} pb.

This is puzzling me a lot! Could you help in this regard?
If any card or log file information is required, I can provide.

Regards,
Samadrita

Question information

Language:
English Edit question
Status:
Answered
For:
MadGraph5_aMC@NLO Edit question
Assignee:
No assignee Edit question
Last query:
Last reply:
Revision history for this message
Olivier Mattelaer (olivier-mattelaer) said :
#1

Hi,

Sorry, I do not know the FeynRules convention and therefore can not really comment on the Lagrangian that you wrote.
If you provide me the UFO model (either by attaching it to a bug report and/or by sending it to me by email at <email address hidden> with a link to this issue) I can take a look at the UFO model and see if MG5aMC handle such model correctly by checking if the cancelation should occur according to the UFO model.
But I would bet that the likely place where something is wrong is either in your field definition of your Lagrangian and/or in some phase handling within FeynRules.

Cheers,

Olivier

Revision history for this message
Samadrita Mukherjee (samadritamisti) said :
#2

Dear Olivier,

Thanks for your prompt reply. I have sent you two UFO files in your email id, one with only left-handed quark doublets (QLbarQL_UFO) and the other containing QLbar.Ga[mu].QL + uRbar.Ga[mu].uR + dRbar.Ga[mu].dR.

The process I am checking is:

define e = e+ e-
generate e p > e j ax

where 'ax' is my new particle.

The respective cross-sections are :

1) Cross-section : 2.284e-07 +- 1.101e-09 pb
     Nb of events : 10000

2) Cross-section : 1.337e-31 +- 4.031e-34 pb
     Nb of events : 10000

I do not understand it because, there is no error message or warning anywhere, not in FeynRules or in MG5.

I wanted to switch on Pythia afterward, in the first case, while MG5 generates 10,000 points, after Pythia, I could only get 0.1% of these events and the Pythia shower jobs are being ended in 2 secs. However, for the second case, at least the showering is working fine.

I am attaching the log files herewith. Please find the UFO files in your email.

Thanks and Regards,
Samadrita

Revision history for this message
Olivier Mattelaer (olivier-mattelaer) said :
#3

Hi,

So for your " QLbar.Ga[mu].QL + uRbar.Ga[mu].uR + dRbar.Ga[mu].dR."
This correspond within the UFO to an interaction for c c~ ax
which is
\slash P_ax

While for the other this correspond to
\slash P_ax (1-\gamma_5)/2

Does this make sense for you given the Lagrangian?

So given those expression, an easy sanity check of MG5aMC is to look at the massless limit since in that case the factor "(1-\gamma)/2" should not have any impact on the process.

In that case:
\slash P_ax returns
     Cross-section : 1.125e-23 +- 4.089e-26 pb
     Nb of events : 10000

and \slash P_ax (1-\gamma_5)/2
     Cross-section : 6.547e-24 +- 1.754e-26 pb
     Nb of events : 10000

Which are both super-small even if not compatible with each other.

So I decide to test for a given phase-space point (with output standalone)

with \slash P_ax (1-\gamma_5)/2 (massless lepton/quark)

 n E px py pz m
 1 0.5000000E+03 0.0000000E+00 0.0000000E+00 0.5000000E+03 0.0000000E+00
 2 0.5000000E+03 0.0000000E+00 0.0000000E+00 -0.5000000E+03 0.0000000E+00
 3 0.4585788E+03 0.1694532E+03 0.3796537E+03 -0.1935025E+03 0.5394797E-05
 4 0.3640666E+03 -0.1832987E+02 -0.3477043E+03 0.1063496E+03 0.1009274E-04
 5 0.1773546E+03 -0.1511234E+03 -0.3194936E+02 0.8715287E+02 0.1000004E-02
 -----------------------------------------------------------------------------
 Matrix element = 7.1064889800383043E-039 GeV^ -2
 -----------------------------------------------------------------------------

with \slash P_ax (massless lepton/quark)
 -----------------------------------------------------------------------------
 Matrix element = 2.1561174109629372E-038 GeV^ -2
 -----------------------------------------------------------------------------

Which again does not agree.
So I decide to check the process via the check command "check e+ c > e+ c ax"
and the result for \slash P_ax
are
> Lorentz invariance results:
> Process Min element Max element Relative diff. Result
> e+ c > e+ c ax 5.0279934147e-39 1.0067466266e-38 5.0057012543e-01 Failed
> JAMP 0 2.0111973659e-38 4.0269865065e-38 5.0057012543e-01 Failed
> Summary: 0/1 passed, 1/1 failed
> Failed processes: e+ c > e+ c ax
> Process permutation results:
> Process Min element Max element Relative diff. Result
> e+ c > e+ c ax 3.7285211440e-38 3.7836224417e-38 1.4669926650e-02 Failed
> Summary: 0/1 passed, 1/1 failed
> Failed processes: Process: e+ c > e+ c ax

The permutation failure is a strong indication that this amplitude is just zero in the massless limit.
Same is happening in the second model \slash P_ax (1-\gamma_5)/2

> Lorentz invariance results:
> Process Min element Max element Relative diff. Result
> e+ c > e+ c ax 1.0354452192e-38 7.6740036477e-38 8.6507105460e-01 Failed
> JAMP 0 4.1417808767e-38 3.0696014591e-37 8.6507105460e-01 Failed
> Summary: 0/1 passed, 1/1 failed
> Failed processes: e+ c > e+ c ax
> Process permutation results:
> Process Min element Max element Relative diff. Result
> e+ c > e+ c ax 2.0571151139e-38 2.9387358771e-38 3.5294117647e-01 Failed
> Summary: 0/1 passed, 1/1 failed
> Failed processes: Process: e+ c > e+ c ax

Now If I return to the massive case and run the check command
I have the same behavior for the \slash P model (small matrix element which seems to be numerically zero)
> Lorentz invariance results:
> Process Min element Max element Relative diff. Result
> e+ c > e+ c ax 1.2910908546e-39 1.1687335644e-38 8.8953077981e-01 Failed
> JAMP 0 5.1643634185e-39 4.6749342575e-38 8.8953077981e-01 Failed
> Summary: 0/1 passed, 1/1 failed
> Failed processes: e+ c > e+ c ax
> Process permutation results:
> Process Min element Max element Relative diff. Result
> e+ c > e+ c ax 2.0307223757e-39 2.2144006811e-39 8.6536151195e-02 Failed
> Summary: 0/1 passed, 1/1 failed
> Failed processes: Process: e+ c > e+ c ax

But for the \slash P (1-gamma_5) model, I do have now a larger matrix-element and numerically stable.

> Lorentz invariance results:
> Process Min element Max element Relative diff. Result
> e+ c > e+ c ax 1.3426019546e-12 1.3426019546e-12 1.8049882115e-15 Passed
> Summary: 1/1 passed, 0/1 failed
> Process permutation results:
> Process Min element Max element Relative diff. Result
> e+ c > e+ c ax 1.8563731387e-13 1.8563731387e-13 0.0000000000e+00 Passed
> Summary: 1/1 passed, 0/1 failed

Now in order to look at the why this happens, it is interesting to look at the amplitude level.
For that, I have taken a look at the process "e+ c > e+ c ax / Z h "
In that case we only have two diagram.
Here is the value evaluated (in the massive case for the 16 helicitiy combination)

> AMP(1), AMP(2) (-1.64679070415800999E-013,-7.35022510959323217E-014) (1.64679070415800671E-013,7.35022510959321703E-014)
> AMP(1), AMP(2) (2.08820911961561582E-011,-4.08588685417426464E-011) (-2.08820911961561000E-011,4.08588685417425754E-011)
> AMP(1), AMP(2) (-3.10446217438834459E-024,-1.23736109850727720E-007) (6.64635935231933234E-036,1.23736109850727535E-007)
> AMP(1), AMP(2) (5.24362616888071501E-005,2.76427360643036656E-006) (-5.24362616888070349E-005,-2.76427360643034538E-006)
> AMP(1), AMP(2) (-1.64703860265815449E-011,3.95535361699618547E-011) (1.64703860265815255E-011,-3.95535361699617771E-011)
> AMP(1), AMP(2) (1.65826729693604588E-013,7.35627520574546213E-014) (-1.65826729693604210E-013,-7.35627520574544951E-014)
> AMP(1), AMP(2) (-1.54042334183659789E-004,4.81147546027791618E-005) (1.54042334183659437E-004,-4.81147546027790670E-005)
> AMP(1), AMP(2) (-1.40175670422713326E-008,8.72862064156581496E-008) (1.40175670422713111E-008,-8.72862064156579908E-008)
> AMP(1), AMP(2) (-1.40175670422713326E-008,-8.72862064156581496E-008) (1.40175670422713111E-008,8.72862064156580173E-008)
> AMP(1), AMP(2) (1.54042334183659789E-004,4.81147546027791618E-005) (-1.54042334183659437E-004,-4.81147546027790670E-005)
> AMP(1), AMP(2) (-1.65826729693604588E-013,7.35627520574546213E-014) (1.65826729693604235E-013,-7.35627520574544951E-014)
> AMP(1), AMP(2) (-1.64703860265815449E-011,-3.95535361699618547E-011) (1.64703860265815255E-011,3.95535361699617771E-011)
> AMP(1), AMP(2) (-5.24362616888071501E-005,2.76427360643036656E-006) (5.24362616888070349E-005,-2.76427360643034538E-006)
> AMP(1), AMP(2) (-3.10446217438834459E-024,1.23736109850727720E-007) (6.64635935231933234E-036,-1.23736109850727508E-007)
> AMP(1), AMP(2) (2.08820911961561549E-011,4.08588685417426464E-011) (-2.08820911961561000E-011,-4.08588685417425754E-011)
> AMP(1), AMP(2) (1.64679070415800999E-013,-7.35022510959323217E-014) (-1.64679070415800671E-013,7.35022510959321703E-014)
> Matrix element = 7.3440087780587226E-038 GeV^ -2

So one can clearly identify that the two diagram are equal up to a minus sign and therefore cancel each other.

If I hack the code to compute the current that connect to the photon, one can see that the sign is already there:

 current diag #1 (2.10581278740892859E-019,-4.32011051362256922E-003) (-5.02582269577534781E-004,2.64945451507959525E-005) (2.64945451507958813E-005,5.02582269577534781E-004) (-2.54946068659848850E-020,-6.79953785058962134E-004)
 curren diag #2 (-1.07831016198817747E-019,4.32011051362256055E-003) (5.02582269577533805E-004,-2.64945451507965183E-005) (-2.64945451507965183E-005,-5.02582269577533805E-004) (-1.05745802763202936E-018,6.79953785058962567E-004)

Which means that the cancelation does not depend of the virtuality of the photon propagator.
So in essence, the same cancelation does occur for the Z/Higgs for the exact same reason.

This is a non trivial cancellation and I can not prove it without the analytical computation.
Knowing that such cancelation occurs already for a e+ > e+ ax (and I did check it numerically)
Should make it textbook level, but I have to stop my investigation here.

So with those knowledge now, you should be able to either proof that
1) Your UFO model is correct or not (by computing the lorentz structure and check the one you expect)
2) Proof that if your UFO model is correct, the result is zero or not. (and potentially proof that MG5aMC is wrong)

At least during those investigation, I did not spot anything weird in the computation.
Within the code, the cancellation does indeed occur and I agree that the reason of why it is occuring is a bit obscur.
This also means that no simple "minus" sign and/or factor of two can be the source of such cancellation.
To me the only possible source of such cancelation is physics (assuming that the UFO model correspond to the physics that you are looking for.)

Obviously if you have the analytical comptutation showing that this should not cancel, I would be more than happy to take a look at the computation and compare it with the code.

Cheers,

Olivier

> On 21 Mar 2023, at 07:10, Samadrita Mukherjee <email address hidden> wrote:
>
> Question #705893 on MadGraph5_aMC@NLO changed:
> https://answers.launchpad.net/mg5amcnlo/+question/705893
>
> Status: Answered => Open
>
> Samadrita Mukherjee is still having a problem:
> Dear Olivier,
>
> Thanks for your prompt reply. I have sent you two UFO files in your
> email id, one with only left-handed quark doublets (QLbarQL_UFO) and the
> other containing QLbar.Ga[mu].QL + uRbar.Ga[mu].uR + dRbar.Ga[mu].dR.
>
> The process I am checking is:
>
> define e = e+ e-
> generate e p > e j ax
>
> where 'ax' is my new particle.
>
> The respective cross-sections are :
>
> 1) Cross-section : 2.284e-07 +- 1.101e-09 pb
> Nb of events : 10000
>
> 2) Cross-section : 1.337e-31 +- 4.031e-34 pb
> Nb of events : 10000
>
> I do not understand it because, there is no error message or warning
> anywhere, not in FeynRules or in MG5.
>
> I wanted to switch on Pythia afterward, in the first case, while MG5
> generates 10,000 points, after Pythia, I could only get 0.1% of these
> events and the Pythia shower jobs are being ended in 2 secs. However,
> for the second case, at least the showering is working fine.
>
> I am attaching the log files herewith. Please find the UFO files in your
> email.
>
> Thanks and Regards,
> Samadrita
>
> --
> You received this question notification because you are an answer
> contact for MadGraph5_aMC@NLO.

Revision history for this message
Olivier Mattelaer (olivier-mattelaer) said :
#4

Hi,

I discussed with a student this morning and he confirms to me that Pslash is indeed leading to some cancellation.
The origin seems to be related to parity.

Cheers,

Olivier

Revision history for this message
Samadrita Mukherjee (samadritamisti) said :
#5

Hi,

Thanks for the detailed reply. I still need to understand it thoroughly. Could you please tell me, in which files of the UFO, the interaction vertices can be seen as you wrote:

"So for your " QLbar.Ga[mu].QL + uRbar.Ga[mu].uR + dRbar.Ga[mu].dR."
This correspond within the UFO to an interaction for c c~ ax
which is
\slash P_ax

While for the other this correspond to
\slash P_ax (1-\gamma_5)/2"

I checked the 'vertices.py' file where in both cases it is the following:

V_77 = Vertex(name = 'V_77',
              particles = [ P.c__tilde__, P.c, P.ax ],
              color = [ 'Identity(1,2)' ],
              lorentz = [ L.FFS1 ],
              couplings = {(0,0):C.GC_6})

Revision history for this message
Olivier Mattelaer (olivier-mattelaer) said :
#6

So from that file, you get the information about the lorenz structure which is FFS1
this is defined in the lorentz.py file:

FFS1 = Lorentz(name = 'FFS1',
               spins = [ 2, 2, 1 ],
               structure = 'P(-1,3)*Gamma(-1,2,1)')

Cheers,

Olivier

Revision history for this message
Samadrita Mukherjee (samadritamisti) said :
#7

Dear Olivier,

With the help of the vertices you correctly pointed out, I have solved the issue at the FeynRules and independently at the UFO files. Now I am able to get the correct order of cross-section in MG5 for the process I am looking into.

Though, while I try to switch on Pythia, out of the 10,000 events that MG5 handles, Pythia only selects 0.1% of the events or even less. The Pythia jobs are finished within 1-2 seconds and in the HEPMC / root I don't have many events.

The terminal output and Log file of Pythia are herewith:

______________________________________________________________________________________________________

NFO: Update the dependent parameter of the param_card.dat
Generating 10000 events with run name run_01
survey run_01
INFO: compile directory
compile Source Directory
Using random number seed offset = 21
INFO: Running Survey
Creating Jobs
Working on SubProcesses
INFO: Compiling for process 1/4.
INFO: P1_emu_emuax
INFO: Compiling for process 2/4.
INFO: P1_emc_emcax
INFO: Compiling for process 3/4.
INFO: P1_emd_emdax
INFO: Compiling for process 4/4.
INFO: P1_ems_emsax
INFO: P1_emu_emuax
INFO: P1_emc_emcax
INFO: P1_emd_emdax
INFO: P1_ems_emsax
INFO: Idle: 1, Running: 1, Completed: 2 [ current time: 13h37 ]
INFO: Idle: 0, Running: 2, Completed: 2 [ current time: 13h37 ]
INFO: Idle: 0, Running: 0, Completed: 4 [ 0.52s ]
INFO: End survey
refine 10000
Creating Jobs
INFO: Refine results to 10000
INFO: Generating 10000.0 unweighted events.
sum of cpu time of last step: 4 seconds
INFO: Effective Luminosity 13200280667.271608 pb^-1
INFO: need to improve 2 channels
- Current estimate of cross-section: 9.090715798000001e-07 +- 1.3543097700426186e-08
    P1_emu_emuax
    P1_emc_emcax
    P1_emd_emdax
    P1_ems_emsax
INFO: Idle: 5, Running: 8, Completed: 0 [ current time: 13h37 ]
INFO: Idle: 0, Running: 0, Completed: 13 [ 1.4s ]
INFO: Combining runs
sum of cpu time of last step: 17 seconds
INFO: finish refine
refine 10000 --treshold=0.9
No need for second refine due to stability of cross-section
INFO: Combining Events
  === Results Summary for run: run_01 tag: tag_1 ===

     Cross-section : 9.139e-07 +- 4.181e-09 pb
     Nb of events : 10000

INFO: No version of lhapdf. Can not run systematics computation
store_events
INFO: Storing parton level results
INFO: End Parton
reweight -from_cards
decay_events -from_cards
INFO: Running Pythia8 [arXiv:1410.3012]
Splitting .lhe event file for PY8 parallelization...
Submitting Pythia8 jobs...
Pythia8 shower jobs: 1 Idle, 7 Running, 0 Done [1 seconds]
Pythia8 shower jobs: 0 Idle, 7 Running, 1 Done [2 seconds]
Pythia8 shower jobs: 0 Idle, 6 Running, 2 Done [2 seconds]
Pythia8 shower jobs: 0 Idle, 5 Running, 3 Done [2 seconds]
Pythia8 shower jobs: 0 Idle, 4 Running, 4 Done [2 seconds]
Pythia8 shower jobs: 0 Idle, 3 Running, 5 Done [2 seconds]
Pythia8 shower jobs: 0 Idle, 2 Running, 6 Done [2 seconds]
Pythia8 shower jobs: 0 Idle, 1 Running, 7 Done [2 seconds]
Pythia8 shower jobs: 0 Idle, 0 Running, 8 Done [2 seconds]
Merging results from the split PY8 runs...
INFO: Pythia8 shower finished after 2 seconds.
INFO: No delphes_card detected, so not run Delphes
  === Results Summary for run: run_01 tag: tag_1 ===

     Cross-section : 9.139e-07 +- 4.181e-09 pb
     Nb of events : 10000

INFO: storing files of previous run
INFO: Storing Pythia8 files of previous run
INFO: Done
quit
INFO:
more information in /home/samadrita/HEP-Packages/MG5_aMC_v3_1_1/ep_ejax_SM/index.html
______________________________________________________________________________________________________________

______________________________________________________________________________________________________________

  Processes, with strategy-dependent cross section info
  number xsec (pb) xerr (pb) xmax (pb)
       1 2.2767e-07 1.0438e-09 2.2767e-07

 -------- End LHA initialization information --------
 PYTHIA Error in ProcessContainer::constructProcess: setting mass failed
 PYTHIA Abort from Pythia::next: processLevel failed; giving up
 PYTHIA Abort from Pythia::next: reached end of Les Houches Events File
WARNING in MG5aMC_PY8_interface.cc: Reached end of LHEF after 1222 Normalisation will be decreased by-2.24%.

 *------- PYTHIA Event and Cross Section Statistics -------------------------------------------------------------*
 | |
 | Subprocess Code | Number of events | sigma +- delta |
 | | Tried Selected Accepted | (estimated) (mb) |
 | | | |
 |-----------------------------------------------------------------------------------------------------------------|
 | | | |
 | Les Houches User Process(es) 9999 | 1250 1250 1 | 3.643e-19 3.643e-19 |
 | ... whereof user classification code 1 | 1250 1250 1 | |
 | | | |
 | sum | 1250 1250 1 | 3.643e-19 3.643e-19 |
 | |
 *------- End PYTHIA Event and Cross Section Statistics ----------------------------------------------------------*

 *------- PYTHIA Error and Warning Messages Statistics ----------------------------------------------------------*
 | |
 | times message |
 | |
 | 249 Abort from Pythia::next: processLevel failed; giving up |
 | 1 Abort from Pythia::next: reached end of Les Houches Events File |
 | 1249 Error in ProcessContainer::constructProcess: setting mass failed |
 | 2 Info from SLHAinterface::initSLHA: No MODSEL found, keeping internal SUSY switched off |
 | 2 Info from SLHAinterface::initSLHA: importing MASS entries |
 | 2 Warning in SLHAinterface::initSLHA: ignoring DECAY tables |
 | 2 Warning in SLHAinterface::initSLHA: ignoring MASS entries |
 | 1 Warning in SLHAinterface::initSLHA: ignoring QNUMBERS |
 | 2 Warning in SLHAinterface::initSLHA: ignoring empty DECAY tables |
 | |
 *------- End PYTHIA Error and Warning Messages Statistics ------------------------------------------------------*

 Contribution of sample 0 to the inclusive cross section : 1.82136328e-19 +- 1.82136328e-19

Inclusive cross section: 1.82136328e-19 +- 1.82136328e-19 mb

Revision history for this message
Samadrita Mukherjee (samadritamisti) said :
#8

Dear Olivier,

In my last message, I enquired about very few events processed in Pythia for the process ep -> ej ax where as MG5 has as many events as I want to generate. Pythia only selects 0.1% of the events or even less. The Pythia jobs are finished within 1-2 seconds and in the HEPMC / root I don't have many events.

The terminal output and Log file of Pythia are attached in the earlier query.

I am using MG%-v3.1.1

Any idea why this is happening and how to solve it in order to get a proper HepMC output?

Best Regards,
Samadrita

Revision history for this message
Olivier Mattelaer (olivier-mattelaer) said :
#9

You might need some specific version of pythia8 and some dedicated tuning for DIS process.
Sorry, I have no experience on how to set pythia8 correctly for such type of process, I would suggest to read the pythia8 manual.

CHeers,

Olivier

Can you help with this problem?

Provide an answer of your own, or ask Samadrita Mukherjee for more information if necessary.

To post a message you must log in.