# cross-section difference with different calculation methods

Asked by Philipp Millet on 2018-02-21

Dear all,

First of all, I already had a look at [1] but it did not solve the problem yet.

I use MG to calculate cross-sections for a SUSY model imported in the UFO format. Decay widths and masses are calculated with Spheno and used in the param card.

I calculated cross sections in two ways:

1)
define smuL = se2 se2c
define mu = e2 e2bar
define fgq = u1 u1bar d1 d1bar
generate p p > smuL, (smuL > mu n1, (n1 > fgq fgq mu)) @0

2)
define smuL = se2 se2c
define mu = e2 e2bar
define fgq = u1 u1bar d1 d1bar
generate p p > mu n1 @0

The resulting cross-sections from method 1 is up to 30% larger than the one from method 2 * BR from the model (depending on the parameter point). If I use 'generate p p > smuL, (smuL > mu n1)' , I get consistent results to method 2, so I guess that the three body decay of the n1 is the culprit here.

Concerning the suggested solutions in [1]:

- cut_decays is False in both cases

- the bwcutoff parameter is set to 150000, so I assume that the processes should be identical

- the width of the n1 can be rather small (reaching 10^-10 at masses of ~10^3) but this should be still far enough away from the limit mentioned in [1]. I tested running with a much larger width for version 1 and got the same result (module ratio of the width)

Do you know what could be reason for this difference? Are these differences expected? Is maybe one of the method giving more accurate results?

If you need more input (param_cards, models, etc) please let me know. I'm using Madgraph 2.3.3.

Thanks a lot,
Philipp

## Question information

Language:
English Edit question
Status:
For:
Assignee:
No assignee Edit question
Last query:
2018-02-23
2018-02-23
 Olivier Mattelaer (olivier-mattelaer) said on 2018-02-21: #1

Hi Philipp,

Since you basically have check all the typical problem. This is going to be tough.

My first guess here is that MadGraph and Spheno will not agree on the value of the width.
What is the width computed by MadWidth? (just put Auto in the param_card)

If indeed those two does not agree, then we have found the culprint and can probably find which method is the best.

Cheers,

Olivier

> On 21 Feb 2018, at 22:02, Philipp Millet <email address hidden> wrote:
>
> New question #664785 on MadGraph5_aMC@NLO:
>
> Dear all,
>
> First of all, I already had a look at [1] but it did not solve the problem yet.
>
>
> I use MG to calculate cross-sections for a SUSY model imported in the UFO format. Decay widths and masses are calculated with Spheno and used in the param card.
>
> I calculated cross sections in two ways:
>
> 1)
> define smuL = se2 se2c
> define mu = e2 e2bar
> define fgq = u1 u1bar d1 d1bar
> generate p p > smuL, (smuL > mu n1, (n1 > fgq fgq mu)) @0
>
> 2)
> define smuL = se2 se2c
> define mu = e2 e2bar
> define fgq = u1 u1bar d1 d1bar
> generate p p > mu n1 @0
>
> The resulting cross-sections from method 1 is up to 30% larger than the one from method 2 * BR from the model (depending on the parameter point). If I use 'generate p p > smuL, (smuL > mu n1)' , I get consistent results to method 2, so I guess that the three body decay of the n1 is the culprit here.
>
> Concerning the suggested solutions in [1]:
>
> - cut_decays is False in both cases
>
> - the bwcutoff parameter is set to 150000, so I assume that the processes should be identical
>
> - the width of the n1 can be rather small (reaching 10^-10 at masses of ~10^3) but this should be still far enough away from the limit mentioned in [1]. I tested running with a much larger width for version 1 and got the same result (module ratio of the width)
>
> Do you know what could be reason for this difference? Are these differences expected? Is maybe one of the method giving more accurate results?
>
> If you need more input (param_cards, models, etc) please let me know. I'm using Madgraph 2.3.3.
>
> Thanks a lot,
> Philipp
>
>
>
>
>
> --

 Philipp Millet (millet-u) said on 2018-02-22: #2

Dear Olivier,

thanks for your fast reply. I tried to calculate the neutralino width by setting the width to Auto. However, the generation fails with:

Please note that the automatic computation of the width is
only valid in narrow-width approximation and at tree-level.
Command "launch" interrupted in sub-command:
"set max_npoint_for_channel 0" with error:
MadGraph5Error : Invalid restriction card (not same block)
set(['ye', 'umix', 'minpar', 'snumix', 'sphenolowenergy', 'fwcoef', 'msu2', 'gauge', 'nmix', 'uulmix', 'msd2', 'chargemix', 'rvlamlqd', 'phases', 'msl2', 'decay', 'higgslhc14', 'tu', 'selmix', 'flavorkitqfv', 'higgslhc13', 'higgslhc8', 'effhiggscouplings', 'rvtlle', 'td', 'rvlamudd', 'higgslhc7', 'usqmix', 'imfwcoef', 'uelmix', 'rvtlqd', 'dsqmix', 'msoft', 'msq2', 'higgsfcc100', 'yd', 'spheno', 'rvtudd', 'modsel', 'sminputs', 'yu', 'mse2', 'pseudoscalarmix', 'te', 'gaugegut', 'uurmix', 'udlmix', 'vmix', 'udrmix', 'mass', 'uermix', 'flavorkitlfv', 'rvlamlle', 'hmix', 'scalarmix']) != set(['imyu', 'umix', 'tu', 'snumix', 'imuelmix', 'imuurmix', 'dsqmix', 'imye', 'imyd', 'imte', 'nmix', 'uulmix', 'imrvlamlle', 'imrvtlle', 'effhiggscouplings', 'rvlamlqd', 'phases', 'decay', 'imudlmix', 'selmix', 'imsnumix', 'imuulmix', 'rvtlle', 'imuermix', 'td', 'imrvlamudd', 'imtd', 'usqmix', 'imselmix', 'uelmix', 'rvtlqd', 'imnmix', 'imudrmix', 'ye', 'imrvtlqd', 'imhmix', 'yd', 'imdsqmix', 'imtu', 'sminputs', 'yu', 'rvtudd', 'imvmix', 'imrvtudd', 'rvlamudd', 'imrvlamlqd', 'uurmix', 'udlmix', 'imumix', 'imusqmix', 'vmix', 'udrmix', 'imphases', 'mass', 'uermix', 'te', 'rvlamlle', 'hmix']).
Missing block: imyu,imuelmix,imuurmix,imvmix,imyd,imusqmix,imudlmix,imnmix,imuulmix,imuermix,imrvlamudd,imtd,imte,imselmix,imsnumix,imudrmix,imrvtlqd,imhmix,imtu,imrvtudd,imye,imrvlamlle,imrvlamlqd,imumix,imrvtlle,imphases,imdsqmix
Unknown block : sphenolowenergy,msoft,msu2,gauge,msd2,chargemix,msl2,higgslhc14,minpar,flavorkitqfv,higgslhc13,higgslhc8,higgslhc7,imfwcoef,spheno,fwcoef,higgsfcc100,modsel,mse2,pseudoscalarmix,gaugegut,msq2,flavorkitlfv,scalarmix

I looks like it is missing some blocks in the param_card. When generating events it just takes the default values for these blocks and seems to work fine. You can find an example param_card in [1].

 Olivier Mattelaer (olivier-mattelaer) said on 2018-02-22: #3

Hi,

Actually this is now always going to fail, even with normal generation,
since this is much too dangerous to use random value for missing block (this indicates that your model is not compatible with

To avoid problem with genuine ommision, I have add a command to add the missing block
"update missing" but that command is not available in 2.2.3.

Cheers,

Olivier

> On 22 Feb 2018, at 16:02, Philipp Millet <email address hidden> wrote:
>
> Question #664785 on MadGraph5_aMC@NLO changed:
>
>
> Philipp Millet is still having a problem:
> Dear Olivier,
>
> thanks for your fast reply. I tried to calculate the neutralino width by
> setting the width to Auto. However, the generation fails with:
>
>
> Please note that the automatic computation of the width is
> only valid in narrow-width approximation and at tree-level.
> Command "launch" interrupted in sub-command:
> "set max_npoint_for_channel 0" with error:
> MadGraph5Error : Invalid restriction card (not same block)
> set(['ye', 'umix', 'minpar', 'snumix', 'sphenolowenergy', 'fwcoef', 'msu2', 'gauge', 'nmix', 'uulmix', 'msd2', 'chargemix', 'rvlamlqd', 'phases', 'msl2', 'decay', 'higgslhc14', 'tu', 'selmix', 'flavorkitqfv', 'higgslhc13', 'higgslhc8', 'effhiggscouplings', 'rvtlle', 'td', 'rvlamudd', 'higgslhc7', 'usqmix', 'imfwcoef', 'uelmix', 'rvtlqd', 'dsqmix', 'msoft', 'msq2', 'higgsfcc100', 'yd', 'spheno', 'rvtudd', 'modsel', 'sminputs', 'yu', 'mse2', 'pseudoscalarmix', 'te', 'gaugegut', 'uurmix', 'udlmix', 'vmix', 'udrmix', 'mass', 'uermix', 'flavorkitlfv', 'rvlamlle', 'hmix', 'scalarmix']) != set(['imyu', 'umix', 'tu', 'snumix', 'imuelmix', 'imuurmix', 'dsqmix', 'imye', 'imyd', 'imte', 'nmix', 'uulmix', 'imrvlamlle', 'imrvtlle', 'effhiggscouplings', 'rvlamlqd', 'phases', 'decay', 'imudlmix', 'selmix', 'imsnumix', 'imuulmix', 'rvtlle', 'imuermix', 'td', 'imrvlamudd', 'imtd', 'usqmix', 'imselmix', 'uelmix', 'rvtlqd', 'imnmix', 'imudrmix', 'ye', 'imrvtlqd', 'imhmix', 'yd', 'imdsqmix', 'imtu', 'sminputs', 'yu', 'rvtudd', 'imvmix', 'imrvtudd', 'rvlamudd', 'imrvlamlqd', 'uurmix', 'udlmix', 'imumix', 'imusqmix', 'vmix', 'udrmix', 'imphases', 'mass', 'uermix', 'te', 'rvlamlle', 'hmix']).
> Missing block: imyu,imuelmix,imuurmix,imvmix,imyd,imusqmix,imudlmix,imnmix,imuulmix,imuermix,imrvlamudd,imtd,imte,imselmix,imsnumix,imudrmix,imrvtlqd,imhmix,imtu,imrvtudd,imye,imrvlamlle,imrvlamlqd,imumix,imrvtlle,imphases,imdsqmix
> Unknown block : sphenolowenergy,msoft,msu2,gauge,msd2,chargemix,msl2,higgslhc14,minpar,flavorkitqfv,higgslhc13,higgslhc8,higgslhc7,imfwcoef,spheno,fwcoef,higgsfcc100,modsel,mse2,pseudoscalarmix,gaugegut,msq2,flavorkitlfv,scalarmix
>
> I looks like it is missing some blocks in the param_card. When
> generating events it just takes the default values for these blocks and
> seems to work fine. You can find an example param_card in [1].
>
> [1] https://cernbox.cern.ch/index.php/s/atiWXieZzPbNwOO
>
> --

 Philipp Millet (millet-u) said on 2018-02-22: #4

Hi,

thanks again for your answer. I will check this with the author of the Spheno/UFO model.

In the meantime, how do I use the 'update missing' command? Do I have to specify it before/after the launch command in ./bin/madevent?

Thanks,
Philipp

 Philipp Millet (millet-u) said on 2018-02-22: #5

Ok, I think this is related to the same issue observed in [1]. I will have a look with a newer MG version.

 Olivier Mattelaer (olivier-mattelaer) said on 2018-02-22: #6

When you have the question:

Do you want to edit a card (press enter to bypass editing)?
/------------------------------------------------------------\
| 1. param : param_card.dat |
| 2. run : run_card.dat |
\------------------------------------------------------------/
you can also
- enter the path to a valid card or banner.
- use the 'set' command to modify a parameter directly.
The set option works only for param_card and run_card.
- call an external program (ASperGE/MadWidth/...).
Type 'help' for the list of available command
[0, done, 1, param, 2, run, 3, madanalysis5_parton, enter path][90s to answer]
>

One of the available command that you can type is "update missing".
This will update your param_card to contain all the missing block/entry that are missing.

Cheers,

Olivier

 Philipp Millet (millet-u) said on 2018-02-23: #7

Dear Olivier,

thanks again for your help! I updated to the newest MG version and used the update missing command to write the empty blocks. However, calculating the width of the n1 seems to take around ~30hours (extrapolated after the first 10h) and seems to use only one of the available cores (when generating events, the same madgraph setup uses all available cores). Is there a way to split the calculation onto multiple cores?

Thanks!
Philipp

 Philipp Millet (millet-u) said on 2018-02-23: #8

Using

" compute_widths n1 --body_decay=3.0025 --path=./output/param_card.dat --output=./output/param_card_decaywidth_bodydecay30025_precision_0p01.dat --precision_channel=0.001"

I managed to calculate the n1 width in a more reasonable time. The result is:

DECAY 1000022 5.992400e-10
# BR NDA ID1 ID2 ...
2.655363e-01 3 -1 1 14 # 1.59119972412e-10
2.653361e-01 3 -14 -1 1 # 1.59000004564e-10
2.351145e-01 3 -1 2 13 # 1.4089001298e-10
2.340131e-01 3 -13 -2 1 # 1.40230010044e-10

compared to Spheno values of:

DECAY 1000022 4.22235416E-10 # Chi_1
# BR NDA ID1 ID2
# BR NDA ID1 ID2 ID3
2.64578146E-01 3 1 -1 14 # BR(Chi_1 -> Fd_1 Fd_1^* Fv_2 )
2.64578146E-01 3 -1 1 -14 # BR(Chi_1 -> Fd_1^* Fd_1 Fv_2^* )
2.35421854E-01 3 1 -2 -13 # BR(Chi_1 -> Fd_1 Fu_1^* Fe_2^* )
2.35421854E-01 3 -1 2 13 # BR(Chi_1 -> Fd_1^* Fu_1 Fe_2 )

To come back to the original question: there is a difference between the widths (provided the above command I used is the good one (?)). The cross section for the two calculation methods are (bwcutoff 150000):

- generate p p > smuL, (smuL > mu n1, (n1 > fgq fgq mu)): 0.001364

- generate p p > mu n1 * BR(n1 -> q q mu): 0.000972202

Is it possible to tell which value one should use? How is the width difference affecting the outcome (given that MG should read the width from the slha anyway)?

Thanks a lot!
Philipp

 Olivier Mattelaer (olivier-mattelaer) said on 2018-02-23: #9

Hi,

Your command is valid but prevent you to include 4 body decay (which I guess is the part that was slowing you down since indeed the diagram generation part is single core).
At the same time, since Spheno does not have any 4 body decay and that we have the same channel it should be fine.

So clearly the problem is that the width does not match.
In this case:
generate p p > mu n1 @0
The cross-section does not depend of the width at all
while in the case of
generate p p > mu n1, (n1 > fgq fgq mu) @0

This is proportional to BR= Gamma_{n1 > fgq fgq mu}/ Gamma_tot.
Note that we do not use the BR information (since we do not use the NWA) and that we do a phase space integration.
But the final result will be proportional to 1/Gamma_tot (where Gamma_tot is your input in the param_card).

Now since we do integrate over the full matrix-element and do not use either Gamma_{n1 > fgq fgq mu} or BR, you will face the situation that
if you do not take the LO width in the param_card then you can have "weird" result as you observe here.

This being explain, the question is what is the correct width?
Do you know if Spheno is computing the LO width? if not then this explains the issue.

If it does, then you can cross-check those widths by integrating them with madevent:

generate 1000022 > 1 -1 14
add process 1000022 > 1 -1 -14
add process 1000022 > 1 -1 -13
add process 1000022 > 1 1 -13
output
launch -f

This is independent of the Auto method since the way to generate the matrix-element is different.
The phase-space integrator is on the other hand basically the same.

Cheers,

Olivier

> On 23 Feb 2018, at 16:07, Philipp Millet <email address hidden> wrote:
>
> Question #664785 on MadGraph5_aMC@NLO changed:
>
> Using
>
> " compute_widths n1 --body_decay=3.0025 --path=./output/param_card.dat
> --output=./output/param_card_decaywidth_bodydecay30025_precision_0p01.dat
> --precision_channel=0.001"
>
> I managed to calculate the n1 width in a more reasonable time. The
> result is:
>
> DECAY 1000022 5.992400e-10
> # BR NDA ID1 ID2 ...
> 2.655363e-01 3 -1 1 14 # 1.59119972412e-10
> 2.653361e-01 3 -14 -1 1 # 1.59000004564e-10
> 2.351145e-01 3 -1 2 13 # 1.4089001298e-10
> 2.340131e-01 3 -13 -2 1 # 1.40230010044e-10
>
> compared to Spheno values of:
>
> DECAY 1000022 4.22235416E-10 # Chi_1
> # BR NDA ID1 ID2
> # BR NDA ID1 ID2 ID3
> 2.64578146E-01 3 1 -1 14 # BR(Chi_1 -> Fd_1 Fd_1^* Fv_2 )
> 2.64578146E-01 3 -1 1 -14 # BR(Chi_1 -> Fd_1^* Fd_1 Fv_2^* )
> 2.35421854E-01 3 1 -2 -13 # BR(Chi_1 -> Fd_1 Fu_1^* Fe_2^* )
> 2.35421854E-01 3 -1 2 13 # BR(Chi_1 -> Fd_1^* Fu_1 Fe_2 )
>
> To come back to the original question: there is a difference between the
> widths (provided the above command I used is the good one (?)). The
> cross section for the two calculation methods are (bwcutoff 150000):
>
> - generate p p > smuL, (smuL > mu n1, (n1 > fgq fgq mu)): 0.001364
>
> - generate p p > mu n1 * BR(n1 -> q q mu): 0.000972202
>
> Is it possible to tell which value one should use? How is the width
> difference affecting the outcome (given that MG should read the width
> from the slha anyway)?
>
> Thanks a lot!
> Philipp
>
> --