Merge lp:~mg5core1/mg5amcnlo/2.6.4 into lp:mg5amcnlo

Proposed by Olivier Mattelaer on 2018-11-09
Status: Merged
Merged at revision: 280
Proposed branch: lp:~mg5core1/mg5amcnlo/2.6.4
Merge into: lp:mg5amcnlo
Diff against target: 5979 lines (+775/-3545)
68 files modified
MadSpin/decay.py (+12/-3)
MadSpin/interface_madspin.py (+23/-8)
MadSpin/src/driver.f (+5/-5)
Template/LO/Source/run.inc (+4/-1)
Template/LO/Source/transpole.f (+20/-1)
Template/LO/SubProcesses/genps.f (+2/-2)
Template/LO/SubProcesses/myamp.f (+39/-112)
Template/LO/SubProcesses/reweight.f (+21/-7)
Template/NLO/MCatNLO/Makefile_MadFKS (+1/-1)
Template/NLO/SubProcesses/check_poles.f (+2/-2)
Template/NLO/SubProcesses/combine_plots_FO.sh (+56/-110)
Template/NLO/SubProcesses/driver_mintFO.f (+1/-1)
Template/NLO/SubProcesses/driver_mintMC.f (+1/-1)
Template/NLO/SubProcesses/makefile (+0/-6)
Template/NLO/SubProcesses/mint-integrator2.f (+3/-3)
Template/NLO/SubProcesses/read40.for (+1/-1)
Template/NLO/SubProcesses/split40.f (+0/-51)
UpdateNotes.txt (+17/-0)
VERSION (+2/-2)
aloha/aloha_writers.py (+3/-1)
madgraph/core/base_objects.py (+8/-2)
madgraph/core/helas_objects.py (+7/-0)
madgraph/interface/amcatnlo_interface.py (+4/-0)
madgraph/interface/amcatnlo_run_interface.py (+16/-2)
madgraph/interface/common_run_interface.py (+103/-13)
madgraph/interface/extended_cmd.py (+8/-0)
madgraph/interface/loop_interface.py (+13/-5)
madgraph/interface/madevent_interface.py (+16/-8)
madgraph/interface/madgraph_interface.py (+17/-5)
madgraph/interface/reweight_interface.py (+5/-4)
madgraph/iolibs/export_v4.py (+67/-5)
madgraph/iolibs/helas_call_writers.py (+23/-0)
madgraph/iolibs/template_files/matrix_loop_induced_madevent_group.inc (+1/-0)
madgraph/iolibs/template_files/matrix_madevent_group_v4.inc (+18/-0)
madgraph/iolibs/template_files/matrix_madevent_v4.inc (+16/-1)
madgraph/iolibs/template_files/matrix_madweight_group_v4.inc (+1/-0)
madgraph/madevent/sum_html.py (+1/-1)
madgraph/various/banner.py (+10/-6)
madgraph/various/lhe_parser.py (+13/-8)
madgraph/various/systematics.py (+1/-0)
models/check_param_card.py (+9/-5)
models/import_ufo.py (+41/-9)
tests/IOTests.py (+11/-5)
tests/acceptance_tests/test_cmd_madloop.py (+2/-2)
tests/acceptance_tests/test_export_fks.py (+1/-1)
tests/acceptance_tests/test_madspin.py (+54/-2)
tests/acceptance_tests/test_model_equivalence.py (+3/-4)
tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%GOLEM_interface.f (+0/-748)
tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%TIR_interface.f (+3/-11)
tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_matrix.f (+1/-1)
tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%GOLEM_interface.f (+0/-748)
tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%TIR_interface.f (+3/-11)
tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%loop_matrix.f (+1/-1)
tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_madevent_group/matrix1.f (+29/-6)
tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_madevent_nogroup/matrix.f (+20/-3)
tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%GOLEM_interface.f (+0/-823)
tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%TIR_interface.f (+3/-11)
tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%loop_matrix.f (+1/-1)
tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%GOLEM_interface.f (+0/-755)
tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%TIR_interface.f (+3/-11)
tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%loop_matrix.f (+1/-1)
tests/test_manager.py (+1/-1)
tests/unit_tests/iolibs/test_export_v4.py (+7/-0)
tests/unit_tests/madspin/test_madspin.py (+7/-0)
vendor/CutTools/src/makefile (+10/-4)
vendor/CutTools/src/qcdloop/aaxex.f (+1/-1)
vendor/IREGI/src/qcdloop/ff/aaxex.f (+1/-1)
vendor/IREGI/src/qcdloop/ff/makefile (+1/-1)
To merge this branch: bzr merge lp:~mg5core1/mg5amcnlo/2.6.4
Reviewer Review Type Date Requested Status
MadTeam 2018-11-09 Pending
Review via email: mp+358542@code.launchpad.net

Commit message

2.6.4 is now ready. I will procede to merging it.

To post a comment you must log in.

Preview Diff

[H/L] Next/Prev Comment, [J/K] Next/Prev File, [N/P] Next/Prev Hunk
=== modified file 'MadSpin/decay.py'
--- MadSpin/decay.py 2018-03-22 14:00:54 +0000
+++ MadSpin/decay.py 2018-11-09 10:32:12 +0000
@@ -3963,14 +3963,23 @@
3963 elif mode=='full':3963 elif mode=='full':
3964 stdin_text="5 0 0 0 \n" # before closing, write down the seed 3964 stdin_text="5 0 0 0 \n" # before closing, write down the seed
3965 external = self.calculator[('full',path)]3965 external = self.calculator[('full',path)]
3966 external.stdin.write(stdin_text)3966 try:
3967 external.stdin.write(stdin_text)
3968 except Exception:
3969 continue
3967 ranmar_state=external.stdout.readline()3970 ranmar_state=external.stdout.readline()
3968 ranmar_file=pjoin(path,'ranmar_state.dat')3971 ranmar_file=pjoin(path,'ranmar_state.dat')
3969 ranmar=open(ranmar_file, 'w')3972 ranmar=open(ranmar_file, 'w')
3970 ranmar.write(ranmar_state)3973 ranmar.write(ranmar_state)
3971 ranmar.close()3974 ranmar.close()
3972 external.stdin.close()3975 try:
3973 external.stdout.close()3976 external.stdin.close()
3977 except Exception:
3978 continue
3979 try:
3980 external.stdout.close()
3981 except Exception:
3982 continue
3974 external.terminate()3983 external.terminate()
3975 del external3984 del external
3976 else:3985 else:
39773986
=== modified file 'MadSpin/interface_madspin.py'
--- MadSpin/interface_madspin.py 2018-06-22 06:52:02 +0000
+++ MadSpin/interface_madspin.py 2018-11-09 10:32:12 +0000
@@ -199,10 +199,11 @@
199199
200 self.inputfile = inputfile200 self.inputfile = inputfile
201 if self.options['spinmode'] == 'none' and \201 if self.options['spinmode'] == 'none' and \
202 (self.options['input_format'] != 'lhe' or (self.options['input_format'] == 'auto' and '.lhe' in inputfile[:-5])): 202 (self.options['input_format'] not in ['lhe','auto'] or
203 (self.options['input_format'] == 'auto' and '.lhe' not in inputfile[-7:])):
203 self.banner = banner.Banner()204 self.banner = banner.Banner()
204 self.setup_for_pure_decay()205 self.setup_for_pure_decay()
205 return 206 return
206 207
207 if inputfile.endswith('.gz'):208 if inputfile.endswith('.gz'):
208 misc.gunzip(inputfile)209 misc.gunzip(inputfile)
@@ -984,8 +985,11 @@
984 #misc.sprint(i, particle.pdg, particle.pid)985 #misc.sprint(i, particle.pdg, particle.pid)
985 #misc.sprint(self.final_state, evt_decayfile)986 #misc.sprint(self.final_state, evt_decayfile)
986 # check if we need to decay the particle 987 # check if we need to decay the particle
987 if not (particle.pdg in self.final_state or particle.pdg in evt_decayfile):988 if self.final_state and particle.pdg not in self.final_state:
988 continue # nothing to do for this particle989 continue # nothing to do for this particle
990 if particle.pdg not in evt_decayfile:
991 continue # nothing to do for this particle
992
989 # check how the decay need to be done993 # check how the decay need to be done
990 nb_decay = len(evt_decayfile[particle.pdg])994 nb_decay = len(evt_decayfile[particle.pdg])
991 if nb_decay == 0:995 if nb_decay == 0:
@@ -1004,11 +1008,18 @@
1004 cumul = 01008 cumul = 0
1005 for j,events in evt_decayfile[particle.pdg].items():1009 for j,events in evt_decayfile[particle.pdg].items():
1006 cumul += events.cross1010 cumul += events.cross
1007 if r < cumul:1011 if r <= cumul:
1008 decay_file = events1012 decay_file = events
1009 decay_file_nb = j1013 decay_file_nb = j
1010 else:
1011 break1014 break
1015 else:
1016 # security for numerical accuracy issue... (unlikely but better safe)
1017 if (cumul-tot)/tot < 1e-5:
1018 decay_file = events
1019 decay_file_nb = j
1020 else:
1021 misc.sprint(j,cumul, events.cross, tot, (tot-cumul)/tot)
1022 raise Exception
1012 1023
1013 if self.options['new_wgt'] == 'BR':1024 if self.options['new_wgt'] == 'BR':
1014 tot_width = float(self.banner.get('param', 'decay', abs(pdg)).value)1025 tot_width = float(self.banner.get('param', 'decay', abs(pdg)).value)
@@ -1056,8 +1067,13 @@
1056 if len(hepmc_output) == 0:1067 if len(hepmc_output) == 0:
1057 hepmc_output.append(lhe_parser.Particle(event=hepmc_output))1068 hepmc_output.append(lhe_parser.Particle(event=hepmc_output))
1058 hepmc_output[0].color2 = 01069 hepmc_output[0].color2 = 0
1070 hepmc_output[0].status = -1
1071 hepmc_output.nexternal+=1
1059 decayed_particle = lhe_parser.Particle(particle, hepmc_output)1072 decayed_particle = lhe_parser.Particle(particle, hepmc_output)
1073 decayed_particle.mother1 = hepmc_output[0]
1074 decayed_particle.mother2 = hepmc_output[0]
1060 hepmc_output.append(decayed_particle)1075 hepmc_output.append(decayed_particle)
1076 hepmc_output.nexternal+=1
1061 decayed_particle.add_decay(decay)1077 decayed_particle.add_decay(decay)
1062 # change the weight associate to the event1078 # change the weight associate to the event
1063 if self.options['new_wgt'] == 'cross-section':1079 if self.options['new_wgt'] == 'cross-section':
@@ -1075,7 +1091,6 @@
1075 else:1091 else:
1076 hepmc_output.wgt = event.wgt1092 hepmc_output.wgt = event.wgt
1077 hepmc_output.nexternal = len(hepmc_output) # the append does not update nexternal1093 hepmc_output.nexternal = len(hepmc_output) # the append does not update nexternal
1078 hepmc_output.assign_mother()
1079 output_lhe.write(str(hepmc_output))1094 output_lhe.write(str(hepmc_output))
1080 else:1095 else:
1081 if counter==0:1096 if counter==0:
@@ -1425,7 +1440,7 @@
1425 if self.options['fixed_order']:1440 if self.options['fixed_order']:
1426 production, counterevt= production[0], production[1:]1441 production, counterevt= production[0], production[1:]
1427 if curr_event and curr_event % 1000 == 0 and float(str(curr_event)[1:]) ==0:1442 if curr_event and curr_event % 1000 == 0 and float(str(curr_event)[1:]) ==0:
1428 print "decaying event number %s. Efficiency: %s [%s s]" % (curr_event, 1/self.efficiency, time.time()-start)1443 logger.info("decaying event number %s. Efficiency: %s [%s s]" % (curr_event, 1/self.efficiency, time.time()-start))
1429 while 1:1444 while 1:
1430 nb_try +=11445 nb_try +=1
1431 decays = self.get_decay_from_file(production, evt_decayfile, nb_event-curr_event)1446 decays = self.get_decay_from_file(production, evt_decayfile, nb_event-curr_event)
14321447
=== modified file 'MadSpin/src/driver.f'
--- MadSpin/src/driver.f 2017-08-03 18:49:25 +0000
+++ MadSpin/src/driver.f 2018-11-09 10:32:12 +0000
@@ -340,10 +340,10 @@
340 if (jac.lt.0d0) then340 if (jac.lt.0d0) then
341 counter2=counter2+1 341 counter2=counter2+1
342 counter3=counter3+1 342 counter3=counter3+1
343c if (counter3.gt.500) then343 if (counter3.gt.500*8*100) then
344c write(*,*) "500_pts_failed_stop_executation"344 write(*,*) "500_pts_failed_stop_executation"
345c stop345 stop
346c endif346 endif
347 if (counter2.ge.8) then ! use another topology to generate PS points347 if (counter2.ge.8) then ! use another topology to generate PS points
348 do k=1,n_max_cg348 do k=1,n_max_cg
349 amp2(k)=0d0349 amp2(k)=0d0
@@ -960,7 +960,7 @@
960 include 'nexternal_prod.inc'960 include 'nexternal_prod.inc'
961 !include 'run.inc'961 !include 'run.inc'
962c arguments962c arguments
963 double precision jac,x(36),p(0:3,nexternal), p_prod(0:3,nexternal)963 double precision jac,x(36),p(0:3,nexternal), p_prod(0:3,nexternal_prod)
964 integer itree(2,-nexternal:-1)964 integer itree(2,-nexternal:-1)
965 double precision qmass(-nexternal:0),qwidth(-nexternal:0)965 double precision qmass(-nexternal:0),qwidth(-nexternal:0)
966c common966c common
967967
=== modified file 'Template/LO/Source/run.inc'
--- Template/LO/Source/run.inc 2018-02-01 22:30:56 +0000
+++ Template/LO/Source/run.inc 2018-11-09 10:32:12 +0000
@@ -91,4 +91,7 @@
91 double precision mxxmin4pdg(0:25) 91 double precision mxxmin4pdg(0:25)
92 logical mxxpart_antipart(1:25)92 logical mxxpart_antipart(1:25)
93 common/TO_PDG_SPECIFIC_CUT/pdg_cut, ptmin4pdg,ptmax4pdg, Emin4pdg, Emax4pdg, etamin4pdg,93 common/TO_PDG_SPECIFIC_CUT/pdg_cut, ptmin4pdg,ptmax4pdg, Emin4pdg, Emax4pdg, etamin4pdg,
94 &etamax4pdg, mxxmin4pdg,mxxpart_antipart
95\ No newline at end of file94\ No newline at end of file
95 &etamax4pdg, mxxmin4pdg,mxxpart_antipart
96
97 double precision small_width_treatment
98 common/narrow_width/small_width_treatment
96\ No newline at end of file99\ No newline at end of file
97100
=== modified file 'Template/LO/Source/transpole.f'
--- Template/LO/Source/transpole.f 2016-06-15 22:46:42 +0000
+++ Template/LO/Source/transpole.f 2018-11-09 10:32:12 +0000
@@ -28,13 +28,24 @@
28 double precision z,zmin,zmax,xmin,xmax,ez28 double precision z,zmin,zmax,xmin,xmax,ez
29 double precision pole1,width1,x,xc29 double precision pole1,width1,x,xc
30 double precision a,b30 double precision a,b
31c
32c small width treatment
33c
34 double precision small_width_treatment
35 common/narrow_width/small_width_treatment
31c-----36c-----
32c Begin Code37c Begin Code
33c-----38c-----
34 pole=pole139 pole=pole1
35 width=width140 width=width1
41
36 x = x142 x = x1
37 if (pole .gt. 0d0) then43 if (pole .gt. 0d0) then
44 if (width.lt.pole*small_width_treatment)then
45 width = pole * small_width_treatment
46 jac = jac * width/width1
47 endif
48
38 zmin = atan((-pole)/width)/width49 zmin = atan((-pole)/width)/width
39 zmax = atan((1d0-pole)/width)/width50 zmax = atan((1d0-pole)/width)/width
40 if (x .gt. del .and. x .lt. 1d0-del) then51 if (x .gt. del .and. x .lt. 1d0-del) then
@@ -175,7 +186,11 @@
175c186c
176 double precision pole1,width1,y1,jac187 double precision pole1,width1,y1,jac
177 real*8 x188 real*8 x
178189c
190c small width treatment
191c
192 double precision small_width_treatment
193 common/narrow_width/small_width_treatment
179c194c
180c Local195c Local
181c196c
@@ -191,6 +206,10 @@
191 width=width1206 width=width1
192 y = y1207 y = y1
193 if (pole .gt. 0d0) then !BW 208 if (pole .gt. 0d0) then !BW
209 if (width.lt.pole*small_width_treatment)then
210 width = pole * small_width_treatment
211 jac = jac * width/width1
212 endif
194 zmin = atan((-pole)/width)/width213 zmin = atan((-pole)/width)/width
195 zmax = atan((1d0-pole)/width)/width214 zmax = atan((1d0-pole)/width)/width
196 z = atan((y-pole)/width)/width215 z = atan((y-pole)/width)/width
197216
=== modified file 'Template/LO/SubProcesses/genps.f'
--- Template/LO/SubProcesses/genps.f 2018-06-20 11:54:42 +0000
+++ Template/LO/SubProcesses/genps.f 2018-11-09 10:32:12 +0000
@@ -41,7 +41,7 @@
41 integer mincfig,maxcfig !Range of configurations41 integer mincfig,maxcfig !Range of configurations
42 integer invar42 integer invar
43 double precision wgt !(input and output)43 double precision wgt !(input and output)
44 double precision x(maxdim),p(maxdim) !x,p (output) [p(0:3,nexternal)]44 double precision x(*),p(*) !x,p (output) [p(0:3,nexternal)]
45c45c
46c Local46c Local
47c47c
@@ -568,7 +568,7 @@
568 if (abs(lpp(1)) .eq. 3) m1 = 0.000511d0568 if (abs(lpp(1)) .eq. 3) m1 = 0.000511d0
569 if (abs(lpp(2)) .eq. 3) m2 = 0.000511d0569 if (abs(lpp(2)) .eq. 3) m2 = 0.000511d0
570 if (mass_ion(1).ge.0d0) m1 = mass_ion(1)570 if (mass_ion(1).ge.0d0) m1 = mass_ion(1)
571 if (mass_ion(2).ge.0d0) m1 = mass_ion(2)571 if (mass_ion(2).ge.0d0) m2 = mass_ion(2)
572 if(ebeam(1).lt.m1.and.lpp(1).ne.9) ebeam(1)=m1572 if(ebeam(1).lt.m1.and.lpp(1).ne.9) ebeam(1)=m1
573 if(ebeam(2).lt.m2.and.lpp(2).ne.9) ebeam(2)=m2573 if(ebeam(2).lt.m2.and.lpp(2).ne.9) ebeam(2)=m2
574 pi1(0)=ebeam(1)574 pi1(0)=ebeam(1)
575575
=== modified file 'Template/LO/SubProcesses/myamp.f'
--- Template/LO/SubProcesses/myamp.f 2016-12-17 16:26:56 +0000
+++ Template/LO/SubProcesses/myamp.f 2018-11-09 10:32:12 +0000
@@ -1,91 +1,3 @@
1 double precision function testamp(p)
2c*****************************************************************************
3c Approximates matrix element by propagators
4c*****************************************************************************
5 implicit none
6c
7c Constants
8c
9 include 'genps.inc'
10 include 'maxconfigs.inc'
11 include 'nexternal.inc'
12 double precision zero
13 parameter (zero = 0d0)
14c
15c Arguments
16c
17 double precision p(0:3,nexternal)
18c integer iconfig
19c
20c Local
21c
22 double precision xp(0:3,-nexternal:nexternal)
23 double precision mpole(-nexternal:0),shat,tsgn
24 integer i,j,iconfig
25
26 double precision prmass(-nexternal:0,lmaxconfigs)
27 double precision prwidth(-nexternal:0,lmaxconfigs)
28 integer pow(-nexternal:0,lmaxconfigs)
29 logical first_time
30c
31c Global
32c
33 integer iforest(2,-max_branch:-1,lmaxconfigs)
34 common/to_forest/ iforest
35 integer mapconfig(0:lmaxconfigs), this_config
36 common/to_mconfigs/mapconfig, this_config
37
38 include 'coupl.inc'
39c
40c External
41c
42 double precision dot
43
44 save prmass,prwidth,pow
45 data first_time /.true./
46c-----
47c Begin Code
48c-----
49 iconfig = this_config
50 if (first_time) then
51c include 'props.inc'
52 first_time=.false.
53 endif
54
55 do i=1,nexternal
56 mpole(-i)=0d0
57 do j=0,3
58 xp(j,i)=p(j,i)
59 enddo
60 enddo
61c mpole(-3) = 174**2
62c shat = dot(p(0,1),p(0,2))/(1800)**2
63 shat = dot(p(0,1),p(0,2))/(10)**2
64c shat = 1d0
65 testamp = 1d0
66 tsgn = +1d0
67 do i=-1,-(nexternal-3),-1 !Find all the propagotors
68 if (iforest(1,i,iconfig) .eq. 1) tsgn=-1d0
69 do j=0,3
70 xp(j,i) = xp(j,iforest(1,i,iconfig))
71 $ +tsgn*xp(j,iforest(2,i,iconfig))
72 enddo
73 if (prwidth(i,iconfig) .ne. 0d0 .and. .false.) then
74 testamp=testamp/((dot(xp(0,i),xp(0,i))
75 $ -prmass(i,iconfig)**2)**2
76 $ -(prmass(i,iconfig)*prwidth(i,iconfig))**2)
77 else
78 testamp = testamp/((dot(xp(0,i),xp(0,i)) -
79 $ prmass(i,iconfig)**2)
80 $ **(pow(i,iconfig)))
81 endif
82 testamp=testamp*shat**(pow(i,iconfig))
83c write(*,*) i,iconfig,pow(i,iconfig),prmass(i,iconfig)
84 enddo
85c testamp = 1d0/dot(xp(0,-1),xp(0,-1))
86 testamp=abs(testamp)
87c testamp = 1d0
88 end
891
90 logical function cut_bw(p)2 logical function cut_bw(p)
91c*****************************************************************************3c*****************************************************************************
@@ -114,6 +26,7 @@
11426
115 double precision prmass(-nexternal:0,lmaxconfigs)27 double precision prmass(-nexternal:0,lmaxconfigs)
116 double precision prwidth(-nexternal:0,lmaxconfigs)28 double precision prwidth(-nexternal:0,lmaxconfigs)
29 double precision prwidth_tmp(-nexternal:0,lmaxconfigs)
117 integer pow(-nexternal:0,lmaxconfigs)30 integer pow(-nexternal:0,lmaxconfigs)
118 logical first_time, onshell31 logical first_time, onshell
119 double precision xmass32 double precision xmass
@@ -152,13 +65,14 @@
152c65c
153 double precision dot66 double precision dot
15467
155 save prmass,prwidth,pow68 save prmass,prwidth,pow,prwidth_tmp
156 data first_time /.true./69 data first_time /.true./
157c-----70c-----
158c Begin Code71c Begin Code
159c----- 72c-----
160 cut_bw = .false. !Default is we passed the cut73 cut_bw = .false. !Default is we passed the cut
161 iconfig = this_config74 iconfig = this_config
75
162 if (first_time) then76 if (first_time) then
163 include 'props.inc'77 include 'props.inc'
164 nbw = 078 nbw = 0
@@ -209,10 +123,12 @@
209c write(*,*) 'xmass',xmass,prmass(i,iconfig)123c write(*,*) 'xmass',xmass,prmass(i,iconfig)
210c124c
211c Here we set if the BW is "on-shell" for LesHouches125c Here we set if the BW is "on-shell" for LesHouches
212c126c
127 prwidth_tmp(i,iconfig) = max(prwidth(i,iconfig), prmass(i,iconfig)*small_width_treatment)
128
213 onshell = (abs(xmass - prmass(i,iconfig)) .lt.129 onshell = (abs(xmass - prmass(i,iconfig)) .lt.
214 $ bwcutoff*prwidth(i,iconfig).and.130 $ bwcutoff*prwidth_tmp(i,iconfig).and.
215 $ (prwidth(i,iconfig)/prmass(i,iconfig).lt.0.1d0.or.131 $ (prwidth_tmp(i,iconfig)/prmass(i,iconfig).lt.0.1d0.or.
216 $ gForceBW(i,iconfig).eq.1))132 $ gForceBW(i,iconfig).eq.1))
217 if(onshell)then133 if(onshell)then
218c Remove on-shell forbidden s-channels (gForceBW=2) (JA 2/10/11)134c Remove on-shell forbidden s-channels (gForceBW=2) (JA 2/10/11)
@@ -264,10 +180,10 @@
264c For decay-chain syntax use BWcutoff here too (22/12/14)180c For decay-chain syntax use BWcutoff here too (22/12/14)
265 if (gForceBW(i, iconfig).eq.1) then181 if (gForceBW(i, iconfig).eq.1) then
266 onshell = (abs(xmass - prmass(i,iconfig)) .lt.182 onshell = (abs(xmass - prmass(i,iconfig)) .lt.
267 $ bwcutoff*prwidth(i,iconfig))183 $ bwcutoff*prwidth_tmp(i,iconfig))
268 else184 else
269 onshell = (abs(xmass - prmass(i,iconfig)) .lt.185 onshell = (abs(xmass - prmass(i,iconfig)) .lt.
270 $ 5d0*prwidth(i,iconfig))186 $ 5d0*prwidth_tmp(i,iconfig))
271 endif187 endif
272188
273 if (onshell .and. (lbw(nbw).eq. 2) .or.189 if (onshell .and. (lbw(nbw).eq. 2) .or.
@@ -313,6 +229,7 @@
313229
314 double precision prmass(-nexternal:0,lmaxconfigs)230 double precision prmass(-nexternal:0,lmaxconfigs)
315 double precision prwidth(-nexternal:0,lmaxconfigs)231 double precision prwidth(-nexternal:0,lmaxconfigs)
232 double precision prwidth_tmp(-nexternal:0,lmaxconfigs)
316 integer pow(-nexternal:0,lmaxconfigs)233 integer pow(-nexternal:0,lmaxconfigs)
317234
318 integer idup(nexternal,maxproc,maxsproc)235 integer idup(nexternal,maxproc,maxsproc)
@@ -382,11 +299,21 @@
382299
383c-----300c-----
384c Begin Code301c Begin Code
385c----- 302c-----
303 iconfig = this_config
304c needs to be initialise to avoid segfault
305 do i = -nexternal,-1
306 prwidth(i,iconfig) = 0
307 prmass(i,iconfig) =0
308 enddo
386 include 'props.inc'309 include 'props.inc'
387c etmin = 10310c etmin = 10
388 nt = 0311 nt = 0
389 iconfig = this_config312 do i = -nexternal,-1
313 prwidth_tmp(i,iconfig) = max(prwidth(i,iconfig), prmass(i,iconfig)*small_width_treatment)
314 enddo
315
316
390 mtot = 0d0317 mtot = 0d0
391 etot = 0d0 !Total energy needed318 etot = 0d0 !Total energy needed
392 spmass = 0d0 !Keep track of BW masses for shat319 spmass = 0d0 !Keep track of BW masses for shat
@@ -423,7 +350,7 @@
423c Look for identical particles to map radiation processes350c Look for identical particles to map radiation processes
424 call idenparts(iden_part, iforest(1,-max_branch,iconfig),351 call idenparts(iden_part, iforest(1,-max_branch,iconfig),
425 $ sprop(1,-max_branch,iconfig), gForceBW(-max_branch,iconfig),352 $ sprop(1,-max_branch,iconfig), gForceBW(-max_branch,iconfig),
426 $ prwidth(-nexternal,iconfig))353 $ prwidth_tmp(-nexternal,iconfig))
427354
428c Start loop over propagators355c Start loop over propagators
429 do i=-1,-(nexternal-3),-1356 do i=-1,-(nexternal-3),-1
@@ -446,15 +373,15 @@
446 xm(i)=max(xm(i),max(xqcutij(l1,l2),0d0))373 xm(i)=max(xm(i),max(xqcutij(l1,l2),0d0))
447 endif374 endif
448c write(*,*) 'iconfig,i',iconfig,i375c write(*,*) 'iconfig,i',iconfig,i
449c write(*,*) prwidth(i,iconfig),prmass(i,iconfig)376c write(*,*) prwidth_tmp(i,iconfig),prmass(i,iconfig)
450 if (prwidth(i,iconfig) .gt. 0 ) then377 if (prwidth_tmp(i,iconfig) .gt. 0 ) then
451 nbw=nbw+1378 nbw=nbw+1
452c JA 6/8/2011 Set xe(i) for resonances379c JA 6/8/2011 Set xe(i) for resonances
453 if (gforcebw(i,iconfig).eq.1) then380 if (gforcebw(i,iconfig).eq.1) then
454 xm(i) = max(xm(i), prmass(i,iconfig)-bwcutoff*prwidth(i,iconfig))381 xm(i) = max(xm(i), prmass(i,iconfig)-bwcutoff*prwidth_tmp(i,iconfig))
455 bwcut_for_PS(i) = bwcutoff382 bwcut_for_PS(i) = bwcutoff
456 else if (lbw(nbw).eq.1) then383 else if (lbw(nbw).eq.1) then
457 xm(i) = max(xm(i), prmass(i,iconfig)-5d0*prwidth(i,iconfig))384 xm(i) = max(xm(i), prmass(i,iconfig)-5d0*prwidth_tmp(i,iconfig))
458 bwcut_for_PS(i) = 5d0385 bwcut_for_PS(i) = 5d0
459 else386 else
460 bwcut_for_PS(i) = 5d0387 bwcut_for_PS(i) = 5d0
@@ -465,19 +392,19 @@
465c Either: required onshell and daughter masses too large392c Either: required onshell and daughter masses too large
466c Or: forced and daughter masses too large393c Or: forced and daughter masses too large
467c Or: required offshell and forced394c Or: required offshell and forced
468 if(prwidth(i,iconfig) .gt. 0.and.395 if(prwidth_tmp(i,iconfig) .gt. 0.and.
469 $ (lbw(nbw).eq.1.and.396 $ (lbw(nbw).eq.1.and.
470 $ (prmass(i,iconfig)+bwcut_for_PS(i)*prwidth(i,iconfig).lt.xm(i)397 $ (prmass(i,iconfig)+bwcut_for_PS(i)*prwidth_tmp(i,iconfig).lt.xm(i)
471 $ .or.prmass(i,iconfig)-bwcut_for_PS(i)*prwidth(i,iconfig).gt.dsqrt(stot))398 $ .or.prmass(i,iconfig)-bwcut_for_PS(i)*prwidth_tmp(i,iconfig).gt.dsqrt(stot))
472 $ .or.gforcebw(i,iconfig).eq.1.and.399 $ .or.gforcebw(i,iconfig).eq.1.and.
473 $ prmass(i,iconfig)+bwcutoff*prwidth(i,iconfig).lt.xm(i)400 $ prmass(i,iconfig)+bwcutoff*prwidth_tmp(i,iconfig).lt.xm(i)
474 $ .or.lbw(nbw).eq.2.and.gforcebw(i,iconfig).eq.1))401 $ .or.lbw(nbw).eq.2.and.gforcebw(i,iconfig).eq.1))
475 $ then402 $ then
476c Write results.dat and quit403c Write results.dat and quit
477 call write_null_results()404 call write_null_results()
478 stop405 stop
479 endif406 endif
480 if (prwidth(i,iconfig) .gt. 0 .and. lbw(nbw) .le. 1) then !B.W.407 if (prwidth_tmp(i,iconfig) .gt. 0 .and. lbw(nbw) .le. 1) then !B.W.
481 if (i .eq. -(nexternal-(nincoming+1))) then !This is s-hat408 if (i .eq. -(nexternal-(nincoming+1))) then !This is s-hat
482 j = 3*(nexternal-2)-4+1 !set i to ndim+1409 j = 3*(nexternal-2)-4+1 !set i to ndim+1
483c-----410c-----
@@ -489,14 +416,14 @@
489 $ .or. lbw(nbw).eq.1) then416 $ .or. lbw(nbw).eq.1) then
490 write(*,*) 'Setting PDF BW',j,nbw,prmass(i,iconfig)417 write(*,*) 'Setting PDF BW',j,nbw,prmass(i,iconfig)
491 spole(j)=prmass(i,iconfig)*prmass(i,iconfig)/stot418 spole(j)=prmass(i,iconfig)*prmass(i,iconfig)/stot
492 swidth(j) = prwidth(i,iconfig)*prmass(i,iconfig)/stot419 swidth(j) = prwidth(i,iconfig)*prmass(i,iconfig)/stot ! keep the real width here (important for the jacobian)
493 endif420 endif
494 else if((prmass(i,iconfig)+bwcut_for_PS(i)*prwidth(i,iconfig)).ge.xm(i)421 else if((prmass(i,iconfig)+bwcut_for_PS(i)*prwidth_tmp(i,iconfig)).ge.xm(i)
495 $ .and. iden_part(i).eq.0 .or. lbw(nbw).eq.1) then422 $ .and. iden_part(i).eq.0 .or. lbw(nbw).eq.1) then
496c JA 02/13 Only allow BW if xm below M+5*Gamma423c JA 02/13 Only allow BW if xm below M+5*Gamma
497 write(*,*) 'Setting BW',i,nbw,prmass(i,iconfig)424 write(*,*) 'Setting BW',i,nbw,prmass(i,iconfig)
498 spole(-i)=prmass(i,iconfig)*prmass(i,iconfig)/stot425 spole(-i)=prmass(i,iconfig)*prmass(i,iconfig)/stot
499 swidth(-i) = prwidth(i,iconfig)*prmass(i,iconfig)/stot426 swidth(-i) = prwidth(i,iconfig)*prmass(i,iconfig)/stot ! keep the real width here (important for the jacobian)
500 endif427 endif
501c JA 4/1/2011 Set grid in case there is no BW (radiation process)428c JA 4/1/2011 Set grid in case there is no BW (radiation process)
502 if (swidth(-i) .eq. 0d0 .and.429 if (swidth(-i) .eq. 0d0 .and.
@@ -509,7 +436,7 @@
509c Set spmass for BWs436c Set spmass for BWs
510 if (swidth(-i) .ne. 0d0)437 if (swidth(-i) .ne. 0d0)
511 $ spmass=spmass-xm(i) +438 $ spmass=spmass-xm(i) +
512 $ max(xm(i),prmass(i,iconfig)-bwcut_for_PS(i)*prwidth(i,iconfig))439 $ max(xm(i),prmass(i,iconfig)-bwcut_for_PS(i)*prwidth_tmp(i,iconfig))
513 else !1/x^pow440 else !1/x^pow
514 a=prmass(i,iconfig)**2/stot441 a=prmass(i,iconfig)**2/stot
515c JA 4/1/2011 always set grid442c JA 4/1/2011 always set grid
@@ -530,7 +457,7 @@
530 endif457 endif
531 endif458 endif
532 if (xo.eq.0d0) xo=1d0/stot459 if (xo.eq.0d0) xo=1d0/stot
533c if (prwidth(i, iconfig) .eq. 0d0.or.iden_part(i).gt.0) then 460c if (prwidth_tmp(i, iconfig) .eq. 0d0.or.iden_part(i).gt.0) then
534 call setgrid(-i,xo,a,1)461 call setgrid(-i,xo,a,1)
535c else 462c else
536c write(*,*) 'Using flat grid for BW',i,nbw,463c write(*,*) 'Using flat grid for BW',i,nbw,
@@ -610,7 +537,7 @@
610 spole(i)= -2.0d0 ! 1/s pole537 spole(i)= -2.0d0 ! 1/s pole
611 write(*,*) "Transforming s_hat 1/s ",i,xo, smin, stot538 write(*,*) "Transforming s_hat 1/s ",i,xo, smin, stot
612 else539 else
613 write(*,*) "Transforming s_hat BW ",spole(i),swidth(i)540 write(*,*) "Transforming s_hat BW ",spole(i), max(swidth(i), spole(i)*small_width_treatment)
614 endif541 endif
615 endif542 endif
616543
617544
=== modified file 'Template/LO/SubProcesses/reweight.f'
--- Template/LO/SubProcesses/reweight.f 2018-06-08 14:42:29 +0000
+++ Template/LO/SubProcesses/reweight.f 2018-11-09 10:32:12 +0000
@@ -6,6 +6,7 @@
6c ----6c ----
7c p p > t t~ (up to 2jet)7c p p > t t~ (up to 2jet)
8c p p > w+ (up to 3 jet)8c p p > w+ (up to 3 jet)
9c p p > j j w+ (ordering seems important)
9c p p > z t t~ j j (no MLM needed)10c p p > z t t~ j j (no MLM needed)
10c11c
11c12c
@@ -243,13 +244,15 @@
243 integer iddgluon, iddother, idgluon, idother244 integer iddgluon, iddother, idgluon, idother
244 logical isqcd245 logical isqcd
245 external isqcd246 external isqcd
247 integer get_color
248 external get_color
246249
247 idmo=ipdg(imo)250 idmo=ipdg(imo)
248 idda1=ipdg(ida1)251 idda1=ipdg(ida1)
249 idda2=ipdg(ida2)252 idda2=ipdg(ida2)
250253
251 if (btest(mlevel,4)) then254 if (btest(mlevel,4)) then
252 write(*,*) ' updating ipart for: ',ida1,ida2,' -> ',imo255 write(*,*) 'updating ipart for: ',ida1,ida2,' -> ',imo
253 endif256 endif
254257
255 if (btest(mlevel,4)) then258 if (btest(mlevel,4)) then
@@ -291,7 +294,6 @@
291 $ ' (',ipdg(imo),')'294 $ ' (',ipdg(imo),')'
292 return295 return
293 endif 296 endif
294
295c FS clustering297c FS clustering
296c Transmit parton PDG code for parton vertex298c Transmit parton PDG code for parton vertex
297 if(isjet(idmo)) then299 if(isjet(idmo)) then
@@ -370,14 +372,26 @@
370c quark -> gluon-quark or Z-quark or h-quark or W-quark372c quark -> gluon-quark or Z-quark or h-quark or W-quark
371 ipart(1,imo)=ipart(1,ida2)373 ipart(1,imo)=ipart(1,ida2)
372 ipart(2,imo)=0374 ipart(2,imo)=0
373 else375 else if(iabs(get_color(idmo)).eq.3.and.iabs(get_color(idda1)).eq.3.and.get_color(idda2).eq.1) then
376c exotic q > q' Scalar
377 ipart(1,imo)=ipart(1,ida1)
378 ipart(2,imo)=0
379 else if(iabs(get_color(idmo)).eq.3.and.iabs(get_color(idda2)).eq.3.and.get_color(idda1).eq.1) then
380c exotic q > Scalar q'
381 ipart(1,imo)=ipart(1,ida2)
382 ipart(2,imo)=0
383 else if (get_color(idmo).eq.1) then
374c Color singlet384c Color singlet
375 ipart(1,imo)=ipart(1,ida1)385 ipart(1,imo)=ipart(1,ida1)
376 ipart(2,imo)=ipart(1,ida2)386 ipart(2,imo)=ipart(1,ida2)
387 else
388 write(*,*) idmo,'>', idda1, idda2, 'color', get_color(idmo),'>', get_color(idda1), get_color(idda2)
389 write(*,*) "failed for ipartupdate. Please retry without MLM/default dynamical scale"
390 stop 3
377 endif391 endif
378 392
379 if (btest(mlevel,4)) then393 if (btest(mlevel,4)) then
380 write(*,*) ' -> ',(ipart(i,imo),i=1,2),' (',ipdg(imo),')'394 write(*,*) 'XY -> ',(ipart(i,imo),i=1,2),' (',ipdg(imo),')'
381 endif395 endif
382396
383 return397 return
@@ -695,6 +709,7 @@
695c increasecode gives whether we should increase jcode at next vertex709c increasecode gives whether we should increase jcode at next vertex
696 increasecode=.false.710 increasecode=.false.
697 do n=1,nexternal-2711 do n=1,nexternal-2
712c write(*,*) 'QCD jet status (before n= ',n,'):',(iqjets(i),i=3,nexternal)
698 do i=1,2 ! index of the child in the interaction713 do i=1,2 ! index of the child in the interaction
699 do j=1,2 ! j index of the beam714 do j=1,2 ! j index of the beam
700 if(idacl(n,i).eq.ibeam(j))then715 if(idacl(n,i).eq.ibeam(j))then
@@ -769,7 +784,6 @@
769 pdgm = ipdgcl(imocl(n),igraphs(1),iproc)784 pdgm = ipdgcl(imocl(n),igraphs(1),iproc)
770 pdgid1 = ipdgcl(idacl(n,1),igraphs(1),iproc)785 pdgid1 = ipdgcl(idacl(n,1),igraphs(1),iproc)
771 pdgid2 = ipdgcl(idacl(n,2),igraphs(1),iproc)786 pdgid2 = ipdgcl(idacl(n,2),igraphs(1),iproc)
772
773 if (.not.isqcd(pdgm).and..not.isqcd(pdgid1).and..not.isqcd(pdgid2)) then787 if (.not.isqcd(pdgm).and..not.isqcd(pdgid1).and..not.isqcd(pdgid2)) then
774 ! this is to avoid to do weird stuff for w+ w- z (or h h h) 788 ! this is to avoid to do weird stuff for w+ w- z (or h h h)
775 ! this fix an issue for qq_zttxqq G1594.08789 ! this fix an issue for qq_zttxqq G1594.08
@@ -914,8 +928,8 @@
914 endif928 endif
915 enddo929 enddo
916 njetstore(iconfig)=njets930 njetstore(iconfig)=njets
917 if (btest(mlevel,4))931 if (btest(mlevel,4))
918 $ write(*,*) 'Storing jets: ',(iqjetstore(i,iconfig),i=1,njets)932 $ write(*,*) 'Storing jets: ',(iqjetstore(i,iconfig),i=1,njets)
919c Recluster without requiring chcluster933c Recluster without requiring chcluster
920 goto 100934 goto 100
921 else935 else
922936
=== modified file 'Template/NLO/MCatNLO/Makefile_MadFKS'
--- Template/NLO/MCatNLO/Makefile_MadFKS 2016-02-24 00:19:50 +0000
+++ Template/NLO/MCatNLO/Makefile_MadFKS 2018-11-09 10:32:12 +0000
@@ -1,7 +1,7 @@
1-include ../../Source/make_opts1-include ../../Source/make_opts
22
3DEBUG=3DEBUG=
4FF=$(FC) $(FFLAGS) $(DEBUG)4FF=$(FC) $(FFLAGS) $(DEBUG) -std=legacy
5CPP=$(CXX) $(CXXFLAGS) $(DEBUG) $(INCLOPTION)5CPP=$(CXX) $(CXXFLAGS) $(DEBUG) $(INCLOPTION)
6CC=$(CXX) $(CFLAGS) $(DEBUG) $(INCLOPTION)6CC=$(CXX) $(CFLAGS) $(DEBUG) $(INCLOPTION)
77
88
=== modified file 'Template/NLO/SubProcesses/check_poles.f'
--- Template/NLO/SubProcesses/check_poles.f 2017-06-18 20:58:19 +0000
+++ Template/NLO/SubProcesses/check_poles.f 2018-11-09 10:32:12 +0000
@@ -19,7 +19,7 @@
19 include 'genps.inc'19 include 'genps.inc'
20 include 'nexternal.inc'20 include 'nexternal.inc'
21 include 'nFKSconfigs.inc'21 include 'nFKSconfigs.inc'
22 double precision p(0:3, nexternal), prambo(0:3, nexternal)22 double precision p(0:3, nexternal), prambo(0:3,100)
23 double precision p_born(0:3,nexternal-1)23 double precision p_born(0:3,nexternal-1)
24 common/pborn/p_born24 common/pborn/p_born
25 double precision pswgt25 double precision pswgt
@@ -49,7 +49,7 @@
49 include 'coupl.inc'49 include 'coupl.inc'
50 include 'q_es.inc'50 include 'q_es.inc'
51 integer nsqso 51 integer nsqso
52 double precision pmass(nexternal), pmass_rambo(nexternal)52 double precision pmass(nexternal), pmass_rambo(100)
53 integer nfail53 integer nfail
54 logical first_time54 logical first_time
55 data first_time/.TRUE./55 data first_time/.TRUE./
5656
=== modified file 'Template/NLO/SubProcesses/combine_plots_FO.sh'
--- Template/NLO/SubProcesses/combine_plots_FO.sh 2013-05-25 12:52:41 +0000
+++ Template/NLO/SubProcesses/combine_plots_FO.sh 2018-11-09 10:32:12 +0000
@@ -1,6 +1,4 @@
1#!/bin/bash1#!/bin/bash
2TD=td_mac_intel
3PSOPEN=open
42
5# $? is the value of last executed command. A call to this function3# $? is the value of last executed command. A call to this function
6# after a failure will cause the program to quit the script4# after a failure will cause the program to quit the script
@@ -15,126 +13,74 @@
15fi13fi
16}14}
1715
18thisdir=`pwd`16function combine {
19if [ -f $thisdir/MADatNLO.top ]17counter=0
20then18i=1
21rm -f $thisdir/MADatNLO*.top
22fi
23make read40
24EXENAME=$thisdir/read40
25echo -n "" > dir19echo -n "" > dir
26counterp=020for plot_file in "${@:2}" ; do
27for p in P* ; do21 counter=$(($counter + 1))
28 cd $thisdir22 echo $plot_file >> dir
29 cd $p23 if [[ $(($counter % 40)) == 0 ]]; then
30 rm -f MADatNLO*24 echo -n "" > dir2
31 rm -f *read40*.out25 echo "40" > dir2
32 echo -n "" > dir
33 counter=0
34 i=1
35 for g in $* ; do
36 if [ ! -f $g/MADatNLO*.top ]
37 then
38 echo "No topdrawer file in $p/$g"
39 else
40 counter=$(($counter + 1))
41 echo $g"/MADatNLO.top" >> dir
42 if [[ $(($counter % 40)) == 0 ]]; then
43 echo -n "" > dir2
44 echo "40" > dir2
45 cat dir >> dir2
46 $EXENAME < dir2
47 teststatus Failure in step 1
48 rm -f dir
49 rm -f dir2
50 mv fort.88 MADatNLO$i\.top
51 mv read40.out "S1read40_"$i\.out
52 i=$(($i + 1))
53 echo -n "" > dir
54 fi
55 fi
56 done
57 if [[ $(($counter % 40)) != 0 ]] ; then
58 echo $(($counter % 40)) > dir2
59 cat dir >> dir226 cat dir >> dir2
60 while [[ $(($counter % 40)) -ne 0 ]]; do27 ./read40 < dir2
61 counter=$[$counter-1]28 teststatus Failure in step 1
62 done
63 $EXENAME < dir2
64 teststatus Failure in step 2
65 rm -f dir29 rm -f dir
66 rm -f dir230 rm -f dir2
67 mv fort.88 MADatNLO$i\.top 31 mv fort.88 $1_$i\.top
68 mv read40.out "S2read40_"$i\.out32 mv read40.out "S1read40_"$1_$i\.out
69 fi33 i=$(($i + 1))
70 34 echo -n "" > dir
71 counter=035 fi
72 for m in MADatNLO*.top ; do36done
73 counter=$[$counter+1]37if [[ $(($counter % 40)) != 0 ]] ; then
74 echo $m >> dir38 echo $(($counter % 40)) > dir2
75 done39 cat dir >> dir2
76 if [[ $counter -gt 1 ]] ; then40 ./read40 < dir2
77 echo $counter > dir241 teststatus Failure in step 2
78 cat dir >> dir2
79 $EXENAME < dir2
80 teststatus Failure in step 3
81 rm -f dir2
82 mv fort.88 MADatNLO.top
83 mv read40.out S3read40.out
84 else
85 mv MADatNLO1.top MADatNLO.top
86 fi
87 rm -f dir
88 cd $thisdir
89 counterp=$(($counterp + 1))
90 echo $p"/MADatNLO.top" >> dir
91 echo $counterp $p "done"
92done
93
94cd $thisdir
95make split40
96./split40
97
98i=0
99#pdir* are created by split
100for p in pdir* ; do
101 i=$(($i+1))
102 $EXENAME < $p
103 teststatus Failure in step 4
104 mv fort.88 MADatNLO$i\.top
105 mv read40.out "S4read40_"$i\.out
106 rm $p
107 echo $i $p "done"
108done
109
110if [[ $i -gt 1 ]] ; then
111 echo $i > dir2
112 for p in MADatNLO*.top ; do
113 echo $p >>dir2
114 done
115 $EXENAME < dir2
116 teststatus Failure in step 5
117 rm -f dir42 rm -f dir
118 rm -f dir243 rm -f dir2
119 mv fort.88 MADatNLO.top 44 mv fort.88 $1_$i\.top
120 mv read40.out S5read40.out45 mv read40.out "S2read40_"$1_$i\.out
46fi
47}
48
49make read40
50rm -f MADatNLO*.top
51
52if [ "$#" -gt 64000 ]; then
53 combine "MADatNLO_A" "$@"
54 combine "MADatNLO_B" MADatNLO_A*.top
55 combine "MADatNLO_C" MADatNLO_B*.top
56 combine "MADatNLO_D" MADatNLO_C*.top
57 mv MADatNLO_D_1.top MADatNLO.top
58elif [ "$#" -gt 1600 ]; then
59 combine "MADatNLO_A" "$@"
60 combine "MADatNLO_B" MADatNLO_A*.top
61 combine "MADatNLO_C" MADatNLO_B*.top
62 mv MADatNLO_C_1.top MADatNLO.top
63elif [ "$#" -gt 40 ]; then
64 combine "MADatNLO_A" "$@"
65 combine "MADatNLO_B" MADatNLO_A*.top
66 mv MADatNLO_B_1.top MADatNLO.top
67elif [ "$#" -gt 1 ]; then
68 combine "MADatNLO_A" "$@"
69 mv MADatNLO_A_1.top MADatNLO.top
121else70else
122 mv MADatNLO1.top MADatNLO.top71 cp $1 MADatNLO.top
123fi72fi
12473
125if [ -f $thisdir/read40.errors ]74
75if [ -f read40.errors ]
126then76then
127rm -f $thisdir/read40.errors77rm -f read40.errors
128fi78fi
129find . -name S"*"read40"*".out -exec fgrep -il "READ40 ERROR READ40 ERROR" {} \; > $thisdir/read40.errors79
13080find . -name S"*"read40"*".out -exec fgrep -il "READ40 ERROR READ40 ERROR" {} \; > read40.errors
131if [ ! -s $thisdir/read40.errors ] ; then81
82if [ ! -s read40.errors ] ; then
132 echo 'no errors found'83 echo 'no errors found'
133else84else
134 echo 'ERRORS! see read40.errors for details'85 echo 'ERRORS! see read40.errors for details'
135fi86fi
136
137#./NLO_Born3.py
138#$TD MADatNLO_combined.top "POSTSCR,ORIENT=3"
139#$PSOPEN MADatNLO_combined.ps
140##td2pdf MADatNLO_con
14187
=== modified file 'Template/NLO/SubProcesses/driver_mintFO.f'
--- Template/NLO/SubProcesses/driver_mintFO.f 2017-09-21 13:47:31 +0000
+++ Template/NLO/SubProcesses/driver_mintFO.f 2018-11-09 10:32:12 +0000
@@ -48,7 +48,7 @@
48 include "mint.inc"48 include "mint.inc"
49 integer nhits_in_grids(maxchannels)49 integer nhits_in_grids(maxchannels)
50 real* 8 xgrid(0:nintervals,ndimmax,maxchannels),ymax(nintervals50 real* 8 xgrid(0:nintervals,ndimmax,maxchannels),ymax(nintervals
51 $ ,ndimmax,maxchannels),ymax_virt(maxchannels),ans(nintegrals51 $ ,ndimmax,maxchannels),ymax_virt(0:maxchannels),ans(nintegrals
52 $ ,0:maxchannels),unc(nintegrals,0:maxchannels),chi2(nintegrals52 $ ,0:maxchannels),unc(nintegrals,0:maxchannels),chi2(nintegrals
53 $ ,0:maxchannels),x(ndimmax),itmax_fl53 $ ,0:maxchannels),x(ndimmax),itmax_fl
54 integer ixi_i,iphi_i,iy_ij,vn54 integer ixi_i,iphi_i,iy_ij,vn
5555
=== modified file 'Template/NLO/SubProcesses/driver_mintMC.f'
--- Template/NLO/SubProcesses/driver_mintMC.f 2017-09-21 13:47:31 +0000
+++ Template/NLO/SubProcesses/driver_mintMC.f 2018-11-09 10:32:12 +0000
@@ -55,7 +55,7 @@
55 common /event_normalisation/event_norm55 common /event_normalisation/event_norm
56c For MINT:56c For MINT:
57 real* 8 xgrid(0:nintervals,ndimmax,maxchannels),ymax(nintervals57 real* 8 xgrid(0:nintervals,ndimmax,maxchannels),ymax(nintervals
58 $ ,ndimmax,maxchannels),ymax_virt(maxchannels),ans(nintegrals58 $ ,ndimmax,maxchannels),ymax_virt(0:maxchannels),ans(nintegrals
59 $ ,0:maxchannels),unc(nintegrals,0:maxchannels),chi2(nintegrals59 $ ,0:maxchannels),unc(nintegrals,0:maxchannels),chi2(nintegrals
60 $ ,0:maxchannels),x(ndimmax)60 $ ,0:maxchannels),x(ndimmax)
61 integer ixi_i,iphi_i,iy_ij,vn,nhits_in_grids(maxchannels)61 integer ixi_i,iphi_i,iy_ij,vn,nhits_in_grids(maxchannels)
6262
=== modified file 'Template/NLO/SubProcesses/makefile'
--- Template/NLO/SubProcesses/makefile 2017-06-15 16:22:23 +0000
+++ Template/NLO/SubProcesses/makefile 2018-11-09 10:32:12 +0000
@@ -7,9 +7,6 @@
7# Files for the read40 combiner of top drawer files7# Files for the read40 combiner of top drawer files
8READ40=read40.o8READ40=read40.o
99
10# Files for the split40 splitter of top drawer files
11SPLIT40=split40.o
12
13# Files to collect all the events in the separate integration channels into a single file10# Files to collect all the events in the separate integration channels into a single file
14COLLECT_EVENTS=collect_events.o handling_lhe_events.o fill_MC_mshell.o11COLLECT_EVENTS=collect_events.o handling_lhe_events.o fill_MC_mshell.o
1512
@@ -21,9 +18,6 @@
21read40: $(READ40)18read40: $(READ40)
22 $(FC) $(LDFLAGS) -o read40 $(READ40)19 $(FC) $(LDFLAGS) -o read40 $(READ40)
2320
24split40: $(SPLIT40)
25 $(FC) $(LDFLAGS) -o split40 $(SPLIT40)
26
27collect_events: $(COLLECT_EVENTS) $(LIBS)21collect_events: $(COLLECT_EVENTS) $(LIBS)
28 $(FC) -o collect_events $(COLLECT_EVENTS) $(LINKLIBS) $(LDFLAGS)22 $(FC) -o collect_events $(COLLECT_EVENTS) $(LINKLIBS) $(LDFLAGS)
29 rm handling_lhe_events.o23 rm handling_lhe_events.o
3024
=== modified file 'Template/NLO/SubProcesses/mint-integrator2.f'
--- Template/NLO/SubProcesses/mint-integrator2.f 2017-05-09 14:07:18 +0000
+++ Template/NLO/SubProcesses/mint-integrator2.f 2018-11-09 10:32:12 +0000
@@ -77,7 +77,7 @@
77 $ ,ans(nintegrals,0:maxchannels),unc(nintegrals,0:maxchannels)77 $ ,ans(nintegrals,0:maxchannels),unc(nintegrals,0:maxchannels)
78 $ ,ans3(nintegrals,3),unc3(nintegrals,3),ans_l3(nintegrals)78 $ ,ans3(nintegrals,3),unc3(nintegrals,3),ans_l3(nintegrals)
79 $ ,unc_l3(nintegrals),chi2_l3(nintegrals),vol_chan79 $ ,unc_l3(nintegrals),chi2_l3(nintegrals),vol_chan
80 real * 8 xint_virt(maxchannels),ymax_virt(maxchannels)80 real * 8 xint_virt(maxchannels),ymax_virt(0:maxchannels)
81 $ ,ans_chan(0:maxchannels)81 $ ,ans_chan(0:maxchannels)
82 real * 8 x(ndimmax),vol82 real * 8 x(ndimmax),vol
83 real * 8 xacc(0:nintervals,ndimmax,maxchannels)83 real * 8 xacc(0:nintervals,ndimmax,maxchannels)
@@ -985,8 +985,8 @@
985 integer ndim,imode985 integer ndim,imode
986 include "mint.inc"986 include "mint.inc"
987 real * 8 fun,xgrid(0:nintervals,ndimmax,maxchannels)987 real * 8 fun,xgrid(0:nintervals,ndimmax,maxchannels)
988 $ ,ymax(nintervals,ndimmax,maxchannels),ymax_virt(maxchannels)988 $ ,ymax(nintervals,ndimmax,maxchannels)
989 $ ,x(ndimmax)989 $ ,ymax_virt(0:maxchannels),x(ndimmax)
990 real * 8 dx(ndimmax),xx(ndimmax),vol_chan,dummy990 real * 8 dx(ndimmax),xx(ndimmax),vol_chan,dummy
991 integer icell(ndimmax),ncell(ndimmax),ncell_virt991 integer icell(ndimmax),ncell(ndimmax),ncell_virt
992 integer ifold(ndimmax),kfold(ndimmax)992 integer ifold(ndimmax),kfold(ndimmax)
993993
=== modified file 'Template/NLO/SubProcesses/read40.for'
--- Template/NLO/SubProcesses/read40.for 2012-08-28 21:06:34 +0000
+++ Template/NLO/SubProcesses/read40.for 2018-11-09 10:32:12 +0000
@@ -5,7 +5,7 @@
5 real * 8 xfacts(1:40),step(1:40)5 real * 8 xfacts(1:40),step(1:40)
6 real * 8 xv(1:40,1:1000),yv(1:40,1:1000),er(1:40,1:1000)6 real * 8 xv(1:40,1:1000),yv(1:40,1:1000),er(1:40,1:1000)
7 real * 8 xa(1:1000),ya(1:1000),ea(1:1000)7 real * 8 xa(1:1000),ya(1:1000),ea(1:1000)
8 character * 70 file,file1,filename(1:40)8 character * 700 file,file1,filename(1:40)
9 character * 70 line,vartype,vartype29 character * 70 line,vartype,vartype2
10 character * 70 strscale,strlimx,strlimy,strtit10 character * 70 strscale,strlimx,strlimy,strtit
11 character * 11 chref1,chref711 character * 11 chref1,chref7
1212
=== removed file 'Template/NLO/SubProcesses/split40.f'
--- Template/NLO/SubProcesses/split40.f 2012-08-28 21:06:34 +0000
+++ Template/NLO/SubProcesses/split40.f 1970-01-01 00:00:00 +0000
@@ -1,51 +0,0 @@
1 PROGRAM SPLIT
2 IMPLICIT NONE
3 character*80 buff
4 integer number,i,j
5
6 open (unit=1, file="dir", status="old")
7
8 number=0
9 do
10 read(1,*,end=99,err=99) buff
11 number=number+1
12 enddo
13 99 continue
14 rewind(1)
15 i=1
16 do while (number.gt.40)
17 open(unit=2,file='pdir'//char(48+i),status='unknown')
18 write (2,'(i3)') 40
19 do j=1,40
20 read(1,'(a)',end=98,err=98) buff
21 write (2,'(a)') buff
22 enddo
23c do j=1,9
24c write (2,'(i1)') 1
25c enddo
26 close(2)
27 i=i+1
28 number=number-40
29 enddo
30 if (number.gt.0) then
31 open(unit=2,file='pdir'//char(48+i),status='unknown')
32 if (number.gt.99) then
33 write (2,'(i3)') number
34 elseif(number.gt.9)then
35 write (2,'(i2)') number
36 else
37 write (2,'(i1)') number
38 endif
39 do j=1,number
40 read(1,'(a)',end=98,err=98) buff
41 write (2,'(a)') buff
42 enddo
43c do j=1,number
44c write (2,'(i1)') 1
45c enddo
46 close(2)
47 endif
48
49
50 98 return
51 end
520
=== modified file 'UpdateNotes.txt'
--- UpdateNotes.txt 2018-06-22 06:52:02 +0000
+++ UpdateNotes.txt 2018-11-09 10:32:12 +0000
@@ -1,5 +1,22 @@
1Update notes for MadGraph5_aMC@NLO (in reverse time order)1Update notes for MadGraph5_aMC@NLO (in reverse time order)
22
32.6.4 (//)
4 OM: add specific treatement for small width (at LO only and not for loop-induced)
5 if the width is smaller than 1e-6 times the mass, a fake width (at that value) is used for the
6 numerical evaluation of the matrix-element. S-channel resonances are re-scaled according to
7 narrow-width approximation to return the correct total cross-section (the distribution of events
8 will on the other hand follow the new width).
9 The parameter '1e-6' can be changed by adding to (LO) run_card the parameter: "small_width_treatment"
10 OM: add a new command "install looptools" to trigger the question that is automatically trigger
11 the first time a loop computation is needed.
12 RF: Fixed a bug when using TopDrawer plots for f(N)LO runs, where the combination of the plots could lead
13 to completely wrong histograms/distributions in case of high-precision runs.
14 OM: Fix some MLM crash for some processes (in particular BSM processes with W').
15 OM: Fix a bug in the reweighting due to the new lhe format (the one avoiding some issue with py8)
16 OM: Fix a behavior for negative mass, the width was set to negative in the param_card automatically
17 making the Parton-shower (and other code) to crash since this does not follow the convention.
18 OM: Change compiler flag to support Mojave.
19
32.6.3.2 (22/06/18)202.6.3.2 (22/06/18)
4 OM: Fix a bug in auto-width when mass are below QCD scale.21 OM: Fix a bug in auto-width when mass are below QCD scale.
5 OM: Fix a bug for g b initial state where the mass in the lhe file was not always correctly assigned22 OM: Fix a bug for g b initial state where the mass in the lhe file was not always correctly assigned
623
=== modified file 'VERSION'
--- VERSION 2018-06-22 06:52:02 +0000
+++ VERSION 2018-11-09 10:32:12 +0000
@@ -1,5 +1,5 @@
1version = 2.6.3.21version = 2.6.4
2date = 2018-06-222date = 2018-11-30
33
44
55
66
=== modified file 'aloha/aloha_writers.py'
--- aloha/aloha_writers.py 2018-05-23 14:15:47 +0000
+++ aloha/aloha_writers.py 2018-11-09 10:32:12 +0000
@@ -523,7 +523,7 @@
523 if 'MP' in self.tag:523 if 'MP' in self.tag:
524 out.write(' complex*32 CI\n')524 out.write(' complex*32 CI\n')
525 if KERNEL.has_pi:525 if KERNEL.has_pi:
526 out.write(' double*16 PI\n')526 out.write(' REAL ( KIND = 16 ) PI\n')
527 else:527 else:
528 out.write(' complex*16 CI\n')528 out.write(' complex*16 CI\n')
529 if KERNEL.has_pi:529 if KERNEL.has_pi:
@@ -673,6 +673,8 @@
673 self.has_model_parameter = True673 self.has_model_parameter = True
674 if name.lower() in ['pi', 'as', 'mu_r', 'aewm1','g']:674 if name.lower() in ['pi', 'as', 'mu_r', 'aewm1','g']:
675 return name675 return name
676 if name.startswith(aloha.aloha_prefix):
677 return name
676 return '%s%s' % (aloha.aloha_prefix, name)678 return '%s%s' % (aloha.aloha_prefix, name)
677 679
678 if '_' in name:680 if '_' in name:
679681
=== modified file 'madgraph/core/base_objects.py'
--- madgraph/core/base_objects.py 2018-04-23 20:23:22 +0000
+++ madgraph/core/base_objects.py 2018-11-09 10:32:12 +0000
@@ -1487,7 +1487,7 @@
1487 def change_parameter_name_with_prefix(self, prefix='mdl_'):1487 def change_parameter_name_with_prefix(self, prefix='mdl_'):
1488 """ Change all model parameter by a given prefix.1488 """ Change all model parameter by a given prefix.
1489 Modify the parameter if some of them are identical up to the case"""1489 Modify the parameter if some of them are identical up to the case"""
14901490
1491 lower_dict={}1491 lower_dict={}
1492 duplicate = set()1492 duplicate = set()
1493 keys = self.get('parameters').keys()1493 keys = self.get('parameters').keys()
@@ -1568,7 +1568,13 @@
1568 for key in self['couplings'].keys():1568 for key in self['couplings'].keys():
1569 for coup in self['couplings'][key]:1569 for coup in self['couplings'][key]:
1570 coup.expr = rep_pattern.sub(replace, coup.expr)1570 coup.expr = rep_pattern.sub(replace, coup.expr)
1571 1571
1572 # change form-factor
1573 ff = [l.formfactors for l in self['lorentz'] if hasattr(l, 'formfactors')]
1574 ff = set(sum(ff,[])) # here we have the list of ff used in the model
1575 for f in ff:
1576 f.value = rep_pattern.sub(replace, f.value)
1577
1572 # change mass/width1578 # change mass/width
1573 for part in self['particles']:1579 for part in self['particles']:
1574 if str(part.get('mass')) in one_change:1580 if str(part.get('mass')) in one_change:
15751581
=== modified file 'madgraph/core/helas_objects.py'
--- madgraph/core/helas_objects.py 2017-08-10 13:09:28 +0000
+++ madgraph/core/helas_objects.py 2018-11-09 10:32:12 +0000
@@ -4569,6 +4569,13 @@
4569 return sum([d.get('wavefunctions') for d in \4569 return sum([d.get('wavefunctions') for d in \
4570 self.get('diagrams')], [])4570 self.get('diagrams')], [])
45714571
4572
4573 def get_all_mass_widths(self):
4574 """Gives a list of all widths used by this ME (from propagator)"""
4575
4576 return set([(d.get('mass'),d.get('width')) for d in self.get_all_wavefunctions()])
4577
4578
4572 def get_all_amplitudes(self):4579 def get_all_amplitudes(self):
4573 """Gives a list of all amplitudes for this ME"""4580 """Gives a list of all amplitudes for this ME"""
45744581
45754582
=== modified file 'madgraph/interface/amcatnlo_interface.py'
--- madgraph/interface/amcatnlo_interface.py 2018-06-08 13:35:22 +0000
+++ madgraph/interface/amcatnlo_interface.py 2018-11-09 10:32:12 +0000
@@ -445,7 +445,11 @@
445 proc_type=self.extract_process_type(line)445 proc_type=self.extract_process_type(line)
446 if proc_type[1] not in ['real', 'LOonly']:446 if proc_type[1] not in ['real', 'LOonly']:
447 run_interface.check_compiler(self.options, block=False)447 run_interface.check_compiler(self.options, block=False)
448 #validate_model will reset self._generate_info; to avoid
449 #this store it
450 geninfo = self._generate_info
448 self.validate_model(proc_type[1], coupling_type=proc_type[2])451 self.validate_model(proc_type[1], coupling_type=proc_type[2])
452 self._generate_info = geninfo
449453
450 #now generate the amplitudes as usual454 #now generate the amplitudes as usual
451 #self.options['group_subprocesses'] = 'False'455 #self.options['group_subprocesses'] = 'False'
452456
=== modified file 'madgraph/interface/amcatnlo_run_interface.py'
--- madgraph/interface/amcatnlo_run_interface.py 2018-06-14 15:24:33 +0000
+++ madgraph/interface/amcatnlo_run_interface.py 2018-11-09 10:32:12 +0000
@@ -1053,7 +1053,10 @@
1053 1053
1054 if self.last_mode in ['LO', 'NLO']:1054 if self.last_mode in ['LO', 'NLO']:
1055 self.switch['fixed_order'] = 'ON'1055 self.switch['fixed_order'] = 'ON'
1056 self.switch['fixed_order'] = 'OFF'1056 if self.proc_characteristics['ninitial'] == 1:
1057 self.switch['fixed_order'] = 'ON'
1058 else:
1059 self.switch['fixed_order'] = 'OFF'
10571060
1058 def color_for_fixed_order(self, switch_value):1061 def color_for_fixed_order(self, switch_value):
1059 1062
@@ -1158,6 +1161,11 @@
1158 if self.last_mode in ['LO', 'NLO', 'noshower', 'noshowerLO']:1161 if self.last_mode in ['LO', 'NLO', 'noshower', 'noshowerLO']:
1159 self.switch['shower'] = 'OFF'1162 self.switch['shower'] = 'OFF'
1160 return 1163 return
1164
1165 if self.proc_characteristics['ninitial'] == 1:
1166 self.switch['shower'] = 'OFF'
1167 return
1168
1161 1169
1162 if os.path.exists(pjoin(self.me_dir, 'Cards', 'shower_card.dat')):1170 if os.path.exists(pjoin(self.me_dir, 'Cards', 'shower_card.dat')):
1163 self.switch['shower'] = self.run_card['parton_shower'] 1171 self.switch['shower'] = self.run_card['parton_shower']
@@ -2603,7 +2611,13 @@
2603 devnull = open(os.devnull, 'w') 2611 devnull = open(os.devnull, 'w')
2604 2612
2605 if self.analyse_card['fo_analysis_format'].lower() == 'topdrawer':2613 if self.analyse_card['fo_analysis_format'].lower() == 'topdrawer':
2606 misc.call(['./combine_plots_FO.sh'] + folder_name, \2614 topfiles = []
2615 for job in jobs:
2616 if job['dirname'].endswith('.top'):
2617 topfiles.append(job['dirname'])
2618 else:
2619 topfiles.append(pjoin(job['dirname'],'MADatNLO.top'))
2620 misc.call(['./combine_plots_FO.sh'] + topfiles, \
2607 stdout=devnull, 2621 stdout=devnull,
2608 cwd=pjoin(self.me_dir, 'SubProcesses'))2622 cwd=pjoin(self.me_dir, 'SubProcesses'))
2609 files.cp(pjoin(self.me_dir, 'SubProcesses', 'MADatNLO.top'),2623 files.cp(pjoin(self.me_dir, 'SubProcesses', 'MADatNLO.top'),
26102624
=== modified file 'madgraph/interface/common_run_interface.py'
--- madgraph/interface/common_run_interface.py 2018-06-22 06:52:02 +0000
+++ madgraph/interface/common_run_interface.py 2018-11-09 10:32:12 +0000
@@ -963,13 +963,15 @@
963 else:963 else:
964 return None964 return None
965965
966 def ask_edit_cards(self, cards, mode='fixed', plot=True, first_cmd=None):966 def ask_edit_cards(self, cards, mode='fixed', plot=True, first_cmd=None, from_banner=None,
967 banner=None):
967 """ """968 """ """
968 if not self.options['madanalysis_path']:969 if not self.options['madanalysis_path']:
969 plot = False970 plot = False
970971
971 self.ask_edit_card_static(cards, mode, plot, self.options['timeout'],972 self.ask_edit_card_static(cards, mode, plot, self.options['timeout'],
972 self.ask, first_cmd=first_cmd)973 self.ask, first_cmd=first_cmd, from_banner=from_banner,
974 banner=banner)
973 975
974 for c in cards:976 for c in cards:
975 if not os.path.isabs(c):977 if not os.path.isabs(c):
@@ -3240,7 +3242,19 @@
3240 - Check that no width are too small (raise a warning if this is the case)3242 - Check that no width are too small (raise a warning if this is the case)
3241 3) if dependent is on True check for dependent parameter (automatic for scan)"""3243 3) if dependent is on True check for dependent parameter (automatic for scan)"""
3242 3244
3243 return self.static_check_param_card(path, self, run=run, dependent=dependent)3245 self.static_check_param_card(path, self, run=run, dependent=dependent)
3246
3247 card = param_card_mod.ParamCard(path)
3248 for param in card['decay']:
3249 width = param.value
3250 if width == 0:
3251 continue
3252 try:
3253 mass = card['mass'].get(param.lhacode).value
3254 except Exception:
3255 continue
3256
3257
3244 3258
3245 @staticmethod3259 @staticmethod
3246 def static_check_param_card(path, interface, run=True, dependent=False, 3260 def static_check_param_card(path, interface, run=True, dependent=False,
@@ -4138,6 +4152,21 @@
4138 return CommonRunCmd.install_lhapdf_pdfset_static(lhapdf_config, pdfsets_dir, 4152 return CommonRunCmd.install_lhapdf_pdfset_static(lhapdf_config, pdfsets_dir,
4139 filename.replace('.LHgrid',''), 4153 filename.replace('.LHgrid',''),
4140 lhapdf_version, alternate_path)4154 lhapdf_version, alternate_path)
4155 elif lhapdf_version.startswith('6.'):
4156 # try to do a simple wget
4157 wwwpath = "http://www.hepforge.org/archive/lhapdf/pdfsets/%s/%s.tar.gz"
4158 wwwpath %= ('.'.join(lhapdf_version.split('.')[:2]), filename)
4159 misc.wget(wwwpath, pjoin(pdfsets_dir, '%s.tar.gz' %filename))
4160 misc.call(['tar', '-xzpvf', '%s.tar.gz' %filename],
4161 cwd=pdfsets_dir)
4162 if os.path.exists(pjoin(pdfsets_dir, filename)) or \
4163 os.path.isdir(pjoin(pdfsets_dir, filename)):
4164 logger.info('%s successfully downloaded and stored in %s' \
4165 % (filename, pdfsets_dir))
4166 else:
4167 raise MadGraph5Error, \
4168 'Could not download %s into %s. Please try to install it manually.' \
4169 % (filename, pdfsets_dir)
4141 4170
4142 else:4171 else:
4143 raise MadGraph5Error, \4172 raise MadGraph5Error, \
@@ -4342,7 +4371,7 @@
43424371
4343 4372
4344 4373
4345 def __init__(self, question, cards=[], mode='auto', *args, **opt):4374 def __init__(self, question, cards=[], from_banner=None, banner=None, mode='auto', *args, **opt):
43464375
43474376
4348 self.load_default() 4377 self.load_default()
@@ -4365,7 +4394,8 @@
4365 self.all_vars = set()4394 self.all_vars = set()
4366 self.modified_card = set() #set of cards not in sync with filesystem4395 self.modified_card = set() #set of cards not in sync with filesystem
4367 # need to sync them before editing/leaving4396 # need to sync them before editing/leaving
43684397 self.init_from_banner(from_banner, banner)
4398
4369 #update default path by custom one if specify in cards4399 #update default path by custom one if specify in cards
4370 for card in cards:4400 for card in cards:
4371 if os.path.exists(card) and os.path.sep in cards:4401 if os.path.exists(card) and os.path.sep in cards:
@@ -4379,11 +4409,31 @@
4379 new_conflict = self.all_vars.intersection(new_vars)4409 new_conflict = self.all_vars.intersection(new_vars)
4380 self.conflict.union(new_conflict)4410 self.conflict.union(new_conflict)
4381 self.all_vars.union(new_vars)4411 self.all_vars.union(new_vars)
4412
4413
4414 def init_from_banner(self, from_banner, banner):
4415 """ defined card that need to be initialized from the banner file
4416 from_banner should be a list of card to load from the banner object
4417 """
4418
4419 if from_banner is None:
4420 self.from_banner = {}
4421 return
4422 misc.sprint(from_banner)
4423 self.from_banner = {}
4424 for card in from_banner:
4425 self.from_banner[card] = banner.charge_card(card)
4426 return self.from_banner
4427
43824428
4383 def get_path(self, name, cards):4429 def get_path(self, name, cards):
4384 """initialise the path if requested"""4430 """initialise the path if requested"""
43854431
4386 defname = '%s_default' % name4432 defname = '%s_default' % name
4433
4434 if name in self.from_banner:
4435 return self.from_banner[name]
4436
4387 if isinstance(cards, list):4437 if isinstance(cards, list):
4388 if name in cards:4438 if name in cards:
4389 return True4439 return True
@@ -4417,7 +4467,13 @@
4417 self.pname2block = {}4467 self.pname2block = {}
4418 self.restricted_value = {}4468 self.restricted_value = {}
4419 self.param_card = {}4469 self.param_card = {}
4420 if not self.get_path('param', cards):4470
4471 is_valid_path = self.get_path('param', cards)
4472 if not is_valid_path:
4473 self.param_consistency = False
4474 return []
4475 if isinstance(is_valid_path, param_card_mod.ParamCard):
4476 self.param_card = is_valid_path
4421 self.param_consistency = False4477 self.param_consistency = False
4422 return []4478 return []
44234479
@@ -4441,10 +4497,15 @@
4441 return self.pname2block.keys()4497 return self.pname2block.keys()
4442 4498
4443 def init_run(self, cards):4499 def init_run(self, cards):
4444 4500
4445 self.run_set = []4501 self.run_set = []
4446 if not self.get_path('run', cards):4502 is_valid_path = self.get_path('run', cards)
4447 return []4503 if not is_valid_path:
4504 return []
4505 if isinstance(is_valid_path, banner_mod.RunCard):
4506 self.run_card = is_valid_path
4507 return []
4508
4448 4509
4449 try:4510 try:
4450 self.run_card = banner_mod.RunCard(self.paths['run'], consistency='warning')4511 self.run_card = banner_mod.RunCard(self.paths['run'], consistency='warning')
@@ -4455,6 +4516,7 @@
4455 except IOError:4516 except IOError:
4456 run_card_def = {}4517 run_card_def = {}
44574518
4519
4458 if run_card_def:4520 if run_card_def:
4459 if self.run_card:4521 if self.run_card:
4460 self.run_set = run_card_def.keys() + self.run_card.hidden_param4522 self.run_set = run_card_def.keys() + self.run_card.hidden_param
@@ -5129,7 +5191,6 @@
5129 logger.warning(str(e))5191 logger.warning(str(e))
5130 return5192 return
51315193
5132
5133 start = 05194 start = 0
5134 if len(args) < 2:5195 if len(args) < 2:
5135 logger.warning('Invalid set command %s (need two arguments)' % line)5196 logger.warning('Invalid set command %s (need two arguments)' % line)
@@ -5206,6 +5267,7 @@
5206 5267
5207 if args[0] in ['run_card', 'param_card', 'MadWeight_card', 'shower_card',5268 if args[0] in ['run_card', 'param_card', 'MadWeight_card', 'shower_card',
5208 'delphes_card','madanalysis5_hadron_card','madanalysis5_parton_card']:5269 'delphes_card','madanalysis5_hadron_card','madanalysis5_parton_card']:
5270
5209 if args[1] == 'default':5271 if args[1] == 'default':
5210 logger.info('replace %s by the default card' % args[0],'$MG:BOLD')5272 logger.info('replace %s by the default card' % args[0],'$MG:BOLD')
5211 files.cp(self.paths['%s_default' %args[0][:-5]], self.paths[args[0][:-5]])5273 files.cp(self.paths['%s_default' %args[0][:-5]], self.paths[args[0][:-5]])
@@ -5262,7 +5324,10 @@
52625324
5263 #### RUN CARD5325 #### RUN CARD
5264 if args[start] in [l.lower() for l in self.run_card.keys()] and card in ['', 'run_card']:5326 if args[start] in [l.lower() for l in self.run_card.keys()] and card in ['', 'run_card']:
5327
5265 if args[start] not in self.run_set:5328 if args[start] not in self.run_set:
5329 if card in self.from_banner or 'run' in self.from_banner:
5330 raise Exception, "change not allowed for this card: event already generated!"
5266 args[start] = [l for l in self.run_set if l.lower() == args[start]][0]5331 args[start] = [l for l in self.run_set if l.lower() == args[start]][0]
52675332
5268 if args[start] in self.conflict and card == '':5333 if args[start] in self.conflict and card == '':
@@ -5279,9 +5344,10 @@
5279 logger.info('remove information %s from the run_card' % args[start],'$MG:BOLD')5344 logger.info('remove information %s from the run_card' % args[start],'$MG:BOLD')
5280 del self.run_card[args[start]]5345 del self.run_card[args[start]]
5281 else:5346 else:
5282 if args[0].startswith('sys_') or \5347 lower_name = args[0].lower()
5283 args[0] in self.run_card.list_parameter or \5348 if lower_name.startswith('sys_') or \
5284 args[0] in self.run_card.dict_parameter:5349 lower_name in self.run_card.list_parameter or \
5350 lower_name in self.run_card.dict_parameter:
5285 val = ' '.join(args[start+1:])5351 val = ' '.join(args[start+1:])
5286 val = val.split('#')[0]5352 val = val.split('#')[0]
5287 else:5353 else:
@@ -5626,6 +5692,30 @@
5626 self.run_card['mass_ion1'] != self.run_card['mass_ion2']):5692 self.run_card['mass_ion1'] != self.run_card['mass_ion2']):
5627 raise Exception, "Heavy ion profile for both beam are different but the symmetry used forbids it. \n Please generate your process with \"set group_subprocesses False\"."5693 raise Exception, "Heavy ion profile for both beam are different but the symmetry used forbids it. \n Please generate your process with \"set group_subprocesses False\"."
5628 5694
5695 # check the status of small width status from LO
5696 for param in self.param_card['decay']:
5697 width = param.value
5698 if width == 0 or isinstance(width,str):
5699 continue
5700 try:
5701 mass = self.param_card['mass'].get(param.lhacode).value
5702 except Exception:
5703 continue
5704 if isinstance(mass,str):
5705 continue
5706
5707 if mass:
5708 if abs(width/mass) < self.run_card['small_width_treatment']:
5709 logger.warning("Particle %s will use a fake width ( %s instead of %s ).\n" +
5710 "Cross-section will be rescaled according to NWA if needed." +
5711 "To force exact treatment reduce the value of 'small_width_treatment' parameter of the run_card",
5712 param.lhacode[0], mass*self.run_card['small_width_treatment'], width)
5713 elif abs(width/mass) < 1e-12:
5714 logger.error('The width of particle %s is too small for an s-channel resonance (%s). If you have this particle in an s-channel, this is likely to create numerical instabilities .', param.lhacode[0], width)
5715 if CommonRunCmd.sleep_for_error:
5716 time.sleep(5)
5717 CommonRunCmd.sleep_for_error = False
5718
56295719
5630 ########################################################################5720 ########################################################################
5631 # NLO specific check5721 # NLO specific check
56325722
=== modified file 'madgraph/interface/extended_cmd.py'
--- madgraph/interface/extended_cmd.py 2018-05-16 08:43:26 +0000
+++ madgraph/interface/extended_cmd.py 2018-11-09 10:32:12 +0000
@@ -2618,6 +2618,13 @@
26182618
2619 def postcmd(self, stop, line):2619 def postcmd(self, stop, line):
2620 2620
2621 # for diamond class arch where both branch defines the postcmd
2622 # set it up to be in coop mode
2623 try:
2624 out = super(ControlSwitch,self).postcmd(stop, line)
2625 except AttributeError:
2626 pass
2627
2621 line = line.strip()2628 line = line.strip()
2622 if ';' in line:2629 if ';' in line:
2623 line= [l for l in line.split(';') if l][-1] 2630 line= [l for l in line.split(';') if l][-1]
@@ -2676,6 +2683,7 @@
2676 def check_consistency(self, key, value):2683 def check_consistency(self, key, value):
2677 """check the consistency of the new flag with the old ones"""2684 """check the consistency of the new flag with the old ones"""
2678 2685
2686
2679 if key in self.last_changed:2687 if key in self.last_changed:
2680 self.last_changed.remove(key)2688 self.last_changed.remove(key)
2681 self.last_changed.append(key)2689 self.last_changed.append(key)
26822690
=== modified file 'madgraph/interface/loop_interface.py'
--- madgraph/interface/loop_interface.py 2018-03-22 14:00:54 +0000
+++ madgraph/interface/loop_interface.py 2018-11-09 10:32:12 +0000
@@ -21,6 +21,7 @@
21import time21import time
22import logging22import logging
23import re23import re
24import sys
2425
25import madgraph26import madgraph
26from madgraph import MG4DIR, MG5DIR, MadGraph5Error27from madgraph import MG4DIR, MG5DIR, MadGraph5Error
@@ -494,14 +495,19 @@
494 aloha.mp_precision = aloha_original_quad_mode495 aloha.mp_precision = aloha_original_quad_mode
495496
496497
497 def install_reduction_library(self):498 def install_reduction_library(self, force=False):
498 """Code to install the reduction library if needed"""499 """Code to install the reduction library if needed"""
499 500
500 opt = self.options501 opt = self.options
501 502
502 # Check if first time:503 # Check if first time:
503 if (opt['ninja'] is None) or (os.path.isfile(pjoin(MG5DIR, opt['ninja'],'libninja.a'))): 504 if not force and ((opt['ninja'] is None) or (os.path.isfile(pjoin(MG5DIR, opt['ninja'],'libninja.a')))):
504 return505 return
506
507 # do not trigger the question for tests
508 if 'test_manager.py' in sys.argv[0]:
509 from unittest.case import SkipTest
510 raise SkipTest
505 511
506 logger.info("First output using loop matrix-elements has been detected. Now asking for loop reduction:", '$MG:BOLD')512 logger.info("First output using loop matrix-elements has been detected. Now asking for loop reduction:", '$MG:BOLD')
507 to_install = self.ask('install', '0', ask_class=AskLoopInstaller, timeout=300, 513 to_install = self.ask('install', '0', ask_class=AskLoopInstaller, timeout=300,
@@ -950,9 +956,11 @@
950 if os.path.exists(pjoin(install_dir1, 'collier')):956 if os.path.exists(pjoin(install_dir1, 'collier')):
951 self.code['collier'] = pjoin(install_dir1, 'collier')957 self.code['collier'] = pjoin(install_dir1, 'collier')
952 if os.path.exists(pjoin(install_dir2, 'PJFry','bin','qd-config')):958 if os.path.exists(pjoin(install_dir2, 'PJFry','bin','qd-config')):
953 self.code['collier'] = pjoin(install_dir2, 'PJFry')959 self.code['pjfry'] = pjoin(install_dir2, 'PJFry')
954 if os.path.exists(pjoin(install_dir2, 'golem95')):960 if os.path.exists(pjoin(install_dir2, 'golem95')):
955 self.code['collier'] = pjoin(install_dir2, 'golem95')961 self.code['glem'] = pjoin(install_dir2, 'golem95')
962 if os.path.exists(pjoin(install_dir1, 'ninja')):
963 self.code['ninja'] = pjoin(install_dir2, 'ninja','lib')
956 964
957 # 1. create the question965 # 1. create the question
958 question, allowed_answer = self.create_question(first=True)966 question, allowed_answer = self.create_question(first=True)
959967
=== modified file 'madgraph/interface/madevent_interface.py'
--- madgraph/interface/madevent_interface.py 2018-06-08 13:35:22 +0000
+++ madgraph/interface/madevent_interface.py 2018-11-09 10:32:12 +0000
@@ -4142,9 +4142,14 @@
4142 self.configure_directory(html_opening =False)4142 self.configure_directory(html_opening =False)
4143 self.check_pythia8(args) 4143 self.check_pythia8(args)
41444144
4145 # Update the banner with the pythia card
4146 if not self.banner or len(self.banner) <=1:
4147 # Here the level keyword 'pythia' must not be changed to 'pythia8'.
4148 self.banner = banner_mod.recover_banner(self.results, 'pythia')
4149
4145 # the args are modify and the last arg is always the mode 4150 # the args are modify and the last arg is always the mode
4146 if not no_default:4151 if not no_default:
4147 self.ask_pythia_run_configuration(args[-1], pythia_version=8)4152 self.ask_pythia_run_configuration(args[-1], pythia_version=8, banner=self.banner)
41484153
4149 if self.options['automatic_html_opening']:4154 if self.options['automatic_html_opening']:
4150 misc.open_file(os.path.join(self.me_dir, 'crossx.html'))4155 misc.open_file(os.path.join(self.me_dir, 'crossx.html'))
@@ -4157,10 +4162,7 @@
4157 #"The normalisation of the hepmc output file will be wrong (i.e. non-standard).\n"+\4162 #"The normalisation of the hepmc output file will be wrong (i.e. non-standard).\n"+\
4158 #"Please use 'event_norm = average' in the run_card to avoid this problem.")4163 #"Please use 'event_norm = average' in the run_card to avoid this problem.")
41594164
4160 # Update the banner with the pythia card4165
4161 if not self.banner or len(self.banner) <=1:
4162 # Here the level keyword 'pythia' must not be changed to 'pythia8'.
4163 self.banner = banner_mod.recover_banner(self.results, 'pythia')
4164 4166
4165 if not self.options['mg5amc_py8_interface_path'] or not \4167 if not self.options['mg5amc_py8_interface_path'] or not \
4166 os.path.exists(pjoin(self.options['mg5amc_py8_interface_path'],4168 os.path.exists(pjoin(self.options['mg5amc_py8_interface_path'],
@@ -6126,7 +6128,7 @@
6126 return switch6128 return switch
6127 6129
6128 ############################################################################6130 ############################################################################
6129 def ask_pythia_run_configuration(self, mode=None, pythia_version=6):6131 def ask_pythia_run_configuration(self, mode=None, pythia_version=6, banner=None):
6130 """Ask the question when launching pythia"""6132 """Ask the question when launching pythia"""
6131 6133
6132 pythia_suffix = '' if pythia_version==6 else '%d'%pythia_version6134 pythia_suffix = '' if pythia_version==6 else '%d'%pythia_version
@@ -6187,10 +6189,16 @@
6187 if self.force:6189 if self.force:
6188 return mode6190 return mode
6189 6191
6192 if not banner:
6193 banner = self.banner
6194
6190 if auto:6195 if auto:
6191 self.ask_edit_cards(cards, mode='auto', plot=(pythia_version==6))6196 self.ask_edit_cards(cards, from_banner=['param', 'run'],
6197 mode='auto', plot=(pythia_version==6), banner=banner
6198 )
6192 else:6199 else:
6193 self.ask_edit_cards(cards, plot=(pythia_version==6))6200 self.ask_edit_cards(cards, from_banner=['param', 'run'],
6201 plot=(pythia_version==6), banner=banner)
61946202
6195 return mode6203 return mode
6196 6204
61976205
=== modified file 'madgraph/interface/madgraph_interface.py'
--- madgraph/interface/madgraph_interface.py 2018-06-22 06:52:02 +0000
+++ madgraph/interface/madgraph_interface.py 2018-11-09 10:32:12 +0000
@@ -2782,7 +2782,8 @@
2782 'gauge','lorentz', 'brs', 'cms']2782 'gauge','lorentz', 'brs', 'cms']
2783 _import_formats = ['model_v4', 'model', 'proc_v4', 'command', 'banner']2783 _import_formats = ['model_v4', 'model', 'proc_v4', 'command', 'banner']
2784 _install_opts = ['Delphes', 'MadAnalysis4', 'ExRootAnalysis',2784 _install_opts = ['Delphes', 'MadAnalysis4', 'ExRootAnalysis',
2785 'update', 'Golem95', 'PJFry', 'QCDLoop', 'maddm', 'maddump']2785 'update', 'Golem95', 'PJFry', 'QCDLoop', 'maddm', 'maddump',
2786 'looptools']
2786 2787
2787 # The targets below are installed using the HEPToolsInstaller.py script2788 # The targets below are installed using the HEPToolsInstaller.py script
2788 _advanced_install_opts = ['pythia8','zlib','boost','lhapdf6','lhapdf5','collier',2789 _advanced_install_opts = ['pythia8','zlib','boost','lhapdf6','lhapdf5','collier',
@@ -5840,7 +5841,7 @@
5840 function will overwrite any existing installation of the tool without 5841 function will overwrite any existing installation of the tool without
5841 warnings.5842 warnings.
5842 """5843 """
5843 5844
5844 # Make sure to avoid any border effect on custom_additional_options5845 # Make sure to avoid any border effect on custom_additional_options
5845 add_options = list(additional_options)5846 add_options = list(additional_options)
5846 5847
@@ -5857,6 +5858,10 @@
5857 if args[0] == 'update':5858 if args[0] == 'update':
5858 self.install_update(['update']+install_options['update_options'],wget=program)5859 self.install_update(['update']+install_options['update_options'],wget=program)
5859 return5860 return
5861 elif args[0] == 'looptools':
5862 self.install_reduction_library(force=True)
5863 return
5864
58605865
5861 plugin = self.install_plugin5866 plugin = self.install_plugin
5862 5867
@@ -5903,12 +5908,20 @@
5903 data = urllib.urlopen(cluster_path)5908 data = urllib.urlopen(cluster_path)
5904 except Exception:5909 except Exception:
5905 continue5910 continue
5911 if data.getcode() != 200:
5912 continue
5913
5906 break5914 break
5915
5907 else:5916 else:
5908 raise MadGraph5Error, '''Impossible to connect any of us servers.5917 raise MadGraph5Error, '''Impossible to connect any of us servers.
5909 Please check your internet connection or retry later'''5918 Please check your internet connection or retry later'''
5910 for line in data:5919 for wwwline in data:
5911 split = line.split()5920 split = wwwline.split()
5921 if len(split)!=2:
5922 if '--source' not in line:
5923 source = {0:'uiuc',1:'ucl'}[index]
5924 return self.do_install(line+' --source='+source, paths=paths, additional_options=additional_options)
5912 path[split[0]] = split[1]5925 path[split[0]] = split[1]
59135926
5914################################################################################5927################################################################################
@@ -7037,7 +7050,6 @@
7037 param_card = check_param_card.ParamCard(out_path.getvalue().split('\n'))7050 param_card = check_param_card.ParamCard(out_path.getvalue().split('\n'))
7038 7051
7039 for (block, lhacode) in put_to_one:7052 for (block, lhacode) in put_to_one:
7040 misc.sprint(block, lhacode)
7041 try:7053 try:
7042 param_card[block].get(lhacode).value = 17054 param_card[block].get(lhacode).value = 1
7043 except:7055 except:
70447056
=== modified file 'madgraph/interface/reweight_interface.py'
--- madgraph/interface/reweight_interface.py 2018-06-15 15:10:56 +0000
+++ madgraph/interface/reweight_interface.py 2018-11-09 10:32:12 +0000
@@ -233,8 +233,8 @@
233 233
234 if '=' in order:234 if '=' in order:
235 # get the type NLO QCD/QED/...235 # get the type NLO QCD/QED/...
236 order = order.split('=',1)[1]236 order = order.split('=',1)[1].strip()
237 237
238 # define the list of particles that are needed for the radiation238 # define the list of particles that are needed for the radiation
239 pert = fks_common.find_pert_particles_interactions(model,239 pert = fks_common.find_pert_particles_interactions(model,
240 pert_order = order)['soft_particles']240 pert_order = order)['soft_particles']
@@ -746,7 +746,8 @@
746 # Find new tag in the banner and add information if needed746 # Find new tag in the banner and add information if needed
747 if 'initrwgt' in self.banner and self.output_type == 'default': 747 if 'initrwgt' in self.banner and self.output_type == 'default':
748 if 'name=\'mg_reweighting\'' in self.banner['initrwgt']:748 if 'name=\'mg_reweighting\'' in self.banner['initrwgt']:
749 blockpat = re.compile(r'''<weightgroup name=\'mg_reweighting\'\s*>(?P<text>.*?)</weightgroup>''', re.I+re.M+re.S)749 blockpat = re.compile(r'''<weightgroup name=\'mg_reweighting\'\s*weight_name_strategy=\'includeIdInWeightName\'>(?P<text>.*?)</weightgroup>''', re.I+re.M+re.S)
750 misc.sprint(blockpat, self.banner['initrwgt'])
750 before, content, after = blockpat.split(self.banner['initrwgt'])751 before, content, after = blockpat.split(self.banner['initrwgt'])
751 header_rwgt_other = before + after752 header_rwgt_other = before + after
752 pattern = re.compile('<weight id=\'(?:rwgt_(?P<id>\d+)|(?P<id2>[_\w]+))(?P<rwgttype>\s*|_\w+)\'>(?P<info>.*?)</weight>', re.S+re.I+re.M)753 pattern = re.compile('<weight id=\'(?:rwgt_(?P<id>\d+)|(?P<id2>[_\w]+))(?P<rwgttype>\s*|_\w+)\'>(?P<info>.*?)</weight>', re.S+re.I+re.M)
@@ -1787,7 +1788,7 @@
1787 for i in range(len(pdg)):1788 for i in range(len(pdg)):
1788 if pdg[i] == oldpdg[i]:1789 if pdg[i] == oldpdg[i]:
1789 continue1790 continue
1790 if not self.model or not getattr(self.model, 'get_mass'):1791 if not self.model or not hasattr(self.model, 'get_mass'):
1791 continue1792 continue
1792 if self.model.get_mass(int(pdg[i])) == self.model.get_mass(int(oldpdg[i])):1793 if self.model.get_mass(int(pdg[i])) == self.model.get_mass(int(oldpdg[i])):
1793 continue1794 continue
17941795
=== modified file 'madgraph/iolibs/export_v4.py'
--- madgraph/iolibs/export_v4.py 2018-06-07 13:31:42 +0000
+++ madgraph/iolibs/export_v4.py 2018-11-09 10:32:12 +0000
@@ -12,6 +12,7 @@
12# For more information, visit madgraph.phys.ucl.ac.be and amcatnlo.web.cern.ch12# For more information, visit madgraph.phys.ucl.ac.be and amcatnlo.web.cern.ch
13#13#
14################################################################################14################################################################################
15from madgraph.iolibs.helas_call_writers import HelasCallWriter
15"""Methods and classes to export matrix elements to v4 format."""16"""Methods and classes to export matrix elements to v4 format."""
1617
17import copy18import copy
@@ -93,12 +94,28 @@
93 # - None, madgraph do nothing for initialisation94 # - None, madgraph do nothing for initialisation
94 exporter = 'v4'95 exporter = 'v4'
95 # language of the output 'v4' for Fortran output96 # language of the output 'v4' for Fortran output
96 # 'cpp' for C++ output 97 # 'cpp' for C++ output
97 98
98 99
99 def __init__(self, dir_path = "", opt=None):100 def __init__(self, dir_path = "", opt=None):
100 # cmd_options is a dictionary with all the optional argurment passed at output time 101 # cmd_options is a dictionary with all the optional argurment passed at output time
101 return102
103 # Activate some monkey patching for the helas call writer.
104 helas_call_writers.HelasCallWriter.customize_argument_for_all_other_helas_object = \
105 self.helas_call_writer_custom
106
107
108 # helper function for customise helas writter
109 @staticmethod
110 def custom_helas_call(call, arg):
111 """static method to customise the way aloha function call are written
112 call is the default template for the call
113 arg are the dictionary used for the call
114 """
115 return call, arg
116
117 helas_call_writer_custom = lambda x,y,z: x.custom_helas_call(y,z)
118
102119
103 def copy_template(self, model):120 def copy_template(self, model):
104 return121 return
@@ -160,6 +177,9 @@
160 #place holder to pass information to the run_interface177 #place holder to pass information to the run_interface
161 self.proc_characteristic = banner_mod.ProcCharacteristic()178 self.proc_characteristic = banner_mod.ProcCharacteristic()
162 179
180 # call mother class
181 super(ProcessExporterFortran,self).__init__(dir_path, opt)
182
163 183
164 #===========================================================================184 #===========================================================================
165 # process exporter fortran switch between group and not grouped185 # process exporter fortran switch between group and not grouped
@@ -1834,13 +1854,25 @@
1834 is_clang = misc.detect_if_cpp_compiler_is_clang(compiler)1854 is_clang = misc.detect_if_cpp_compiler_is_clang(compiler)
1835 is_lc = misc.detect_cpp_std_lib_dependence(compiler) == '-lc++'1855 is_lc = misc.detect_cpp_std_lib_dependence(compiler) == '-lc++'
18361856
1857
1837 # list of the variable to set in the make_opts file1858 # list of the variable to set in the make_opts file
1838 for_update= {'DEFAULT_CPP_COMPILER':compiler,1859 for_update= {'DEFAULT_CPP_COMPILER':compiler,
1839 'MACFLAG':'-mmacosx-version-min=10.7' if is_clang and is_lc else '',1860 'MACFLAG':'-mmacosx-version-min=10.7' if is_clang and is_lc else '',
1840 'STDLIB': '-lc++' if is_lc else '-lstdc++',1861 'STDLIB': '-lc++' if is_lc else '-lstdc++',
1841 'STDLIB_FLAG': '-stdlib=libc++' if is_lc and is_clang else ''1862 'STDLIB_FLAG': '-stdlib=libc++' if is_lc and is_clang else ''
1842 }1863 }
1843 1864
1865 # for MOJAVE remove the MACFLAG:
1866 if is_clang:
1867 import platform
1868 version, _, _ = platform.mac_ver()
1869 if not version:# not linux
1870 version = 14 # set version to remove MACFLAG
1871 else:
1872 version = int(version.split('.')[1])
1873 if version >= 14:
1874 for_update['MACFLAG'] = '-mmacosx-version-min=10.8' if is_lc else ''
1875
1844 if not root_dir:1876 if not root_dir:
1845 root_dir = self.dir_path1877 root_dir = self.dir_path
1846 make_opts = pjoin(root_dir, 'Source', 'make_opts')1878 make_opts = pjoin(root_dir, 'Source', 'make_opts')
@@ -3411,6 +3443,7 @@
3411 return s_and_t_channels3443 return s_and_t_channels
34123444
34133445
3446
3414#===============================================================================3447#===============================================================================
3415# ProcessExporterFortranME3448# ProcessExporterFortranME
3416#===============================================================================3449#===============================================================================
@@ -3419,7 +3452,16 @@
3419 MadEvent format."""3452 MadEvent format."""
34203453
3421 matrix_file = "matrix_madevent_v4.inc"3454 matrix_file = "matrix_madevent_v4.inc"
34223455
3456 # helper function for customise helas writter
3457 @staticmethod
3458 def custom_helas_call(call, arg):
3459 if arg['mass'] == '%(M)s,%(W)s,':
3460 arg['mass'] = '%(M)s, fk_%(W)s,'
3461 elif '%(W)s' in arg['mass']:
3462 raise Exception
3463 return call, arg
3464
3423 def copy_template(self, model):3465 def copy_template(self, model):
3424 """Additional actions needed for setup of Template3466 """Additional actions needed for setup of Template
3425 """3467 """
@@ -3443,6 +3485,7 @@
3443 3485
3444 3486
34453487
3488
34463489
34473490
3448 #===========================================================================3491 #===========================================================================
@@ -3933,10 +3976,29 @@
3933 # Extract helas calls3976 # Extract helas calls
3934 helas_calls = fortran_model.get_matrix_element_calls(\3977 helas_calls = fortran_model.get_matrix_element_calls(\
3935 matrix_element)3978 matrix_element)
3979
39363980
3937 replace_dict['helas_calls'] = "\n".join(helas_calls)3981 replace_dict['helas_calls'] = "\n".join(helas_calls)
39383982
39393983
3984 #adding the support for the fake width (forbidding too small width)
3985 mass_width = matrix_element.get_all_mass_widths()
3986 width_list = set([e[1] for e in mass_width])
3987
3988 replace_dict['fake_width_declaration'] = \
3989 (' double precision fk_%s \n' * len(width_list)) % tuple(width_list)
3990 replace_dict['fake_width_declaration'] += \
3991 (' save fk_%s \n' * len(width_list)) % tuple(width_list)
3992 fk_w_defs = []
3993 one_def = ' fk_%(w)s = SIGN(MAX(ABS(%(w)s), ABS(%(m)s*small_width_treatment)), %(w)s)'
3994 for m, w in mass_width:
3995 if w == 'zero':
3996 if ' fk_zero = 0d0' not in fk_w_defs:
3997 fk_w_defs.append(' fk_zero = 0d0')
3998 continue
3999 fk_w_defs.append(one_def %{'m':m, 'w':w})
4000 replace_dict['fake_width_definitions'] = '\n'.join(fk_w_defs)
4001
3940 # Extract version number and date from VERSION file4002 # Extract version number and date from VERSION file
3941 info_lines = self.get_mg5_info_lines()4003 info_lines = self.get_mg5_info_lines()
3942 replace_dict['info_lines'] = info_lines4004 replace_dict['info_lines'] = info_lines
39434005
=== modified file 'madgraph/iolibs/helas_call_writers.py'
--- madgraph/iolibs/helas_call_writers.py 2016-03-03 23:17:06 +0000
+++ madgraph/iolibs/helas_call_writers.py 2018-11-09 10:32:12 +0000
@@ -21,11 +21,13 @@
21import aloha.aloha_writers as aloha_writers21import aloha.aloha_writers as aloha_writers
22import aloha22import aloha
23from madgraph import MadGraph5Error23from madgraph import MadGraph5Error
24import madgraph.various.misc as misc
2425
25class HelasWriterError(Exception):26class HelasWriterError(Exception):
26 """Class for the error of this module """27 """Class for the error of this module """
27 pass28 pass
2829
30
29#===============================================================================31#===============================================================================
30# HelasCallWriter32# HelasCallWriter
31#===============================================================================33#===============================================================================
@@ -41,6 +43,19 @@
41 # Dictionaries from spin states to letters in Helas call43 # Dictionaries from spin states to letters in Helas call
42 mother_dict = {1: 'S', 2: 'O', -2: 'I', 3: 'V', 5: 'T', 4:'OR', -4:'IR',44 mother_dict = {1: 'S', 2: 'O', -2: 'I', 3: 'V', 5: 'T', 4:'OR', -4:'IR',
43 99:'P'}45 99:'P'}
46
47 @staticmethod
48 def customize_argument_for_all_other_helas_object(call, arg):
49 """ Place holder for PLUGIN/...
50 (used by madevent output for small width handling)
51 """
52 return call, arg
53
54 @staticmethod
55 def default_customize_argument_for_all_other_helas_object(call,arg):
56 return call,arg
57 #customize_argument_for_all_other_helas_object = fct_customize_argument_for_all_other_helas_object
58 #default_customize_argument_for_all_other_helas_object = fct_customize_argument_for_all_other_helas_object
4459
45 def default_setup(self):60 def default_setup(self):
4661
@@ -85,6 +100,11 @@
85100
86 return ['model', 'wavefunctions', 'amplitudes']101 return ['model', 'wavefunctions', 'amplitudes']
87102
103
104
105
106
107
88 def get_loop_amp_helas_calls(self, matrix_element):108 def get_loop_amp_helas_calls(self, matrix_element):
89 """Return a list of strings, corresponding to the Helas calls109 """Return a list of strings, corresponding to the Helas calls
90 for building loop amplitudes (AMPL) only."""110 for building loop amplitudes (AMPL) only."""
@@ -1197,6 +1217,7 @@
1197 arg['second_line'] = ampl+"="+ampl+"*(%(uvct)s)" 1217 arg['second_line'] = ampl+"="+ampl+"*(%(uvct)s)"
11981218
1199 # ALL ARGUMENT FORMATTED ###############################################1219 # ALL ARGUMENT FORMATTED ###############################################
1220 call, arg = HelasCallWriter.customize_argument_for_all_other_helas_object(call, arg)
1200 # Store the result.1221 # Store the result.
1201 call = call % arg1222 call = call % arg
1202 # Now we have a line correctly formatted1223 # Now we have a line correctly formatted
@@ -1209,6 +1230,8 @@
1209 else:1230 else:
1210 self.add_amplitude(argument.get_call_key(), call_function)1231 self.add_amplitude(argument.get_call_key(), call_function)
12111232
1233
1234
1212 def get_loop_amplitude_helas_calls(self, loop_matrix_element):1235 def get_loop_amplitude_helas_calls(self, loop_matrix_element):
1213 """ Returns a list of strings corresponding to the Helas calls for each 1236 """ Returns a list of strings corresponding to the Helas calls for each
1214 loop amplitude of this loop matrix element. This function is placed in1237 loop amplitude of this loop matrix element. This function is placed in
12151238
=== modified file 'madgraph/iolibs/template_files/matrix_loop_induced_madevent_group.inc'
--- madgraph/iolibs/template_files/matrix_loop_induced_madevent_group.inc 2017-11-01 10:08:52 +0000
+++ madgraph/iolibs/template_files/matrix_loop_induced_madevent_group.inc 2018-11-09 10:32:12 +0000
@@ -39,6 +39,7 @@
39C When set negative, the security above is removed39C When set negative, the security above is removed
40 DOUBLE PRECISION MULTICHANNEL_THRES40 DOUBLE PRECISION MULTICHANNEL_THRES
41 PARAMETER (MULTICHANNEL_THRES=1.0d-5)41 PARAMETER (MULTICHANNEL_THRES=1.0d-5)
42
42c43c
43c global (due to reading writting) 44c global (due to reading writting)
44c45c
4546
=== modified file 'madgraph/iolibs/template_files/matrix_madevent_group_v4.inc'
--- madgraph/iolibs/template_files/matrix_madevent_group_v4.inc 2017-11-01 10:08:52 +0000
+++ madgraph/iolibs/template_files/matrix_madevent_group_v4.inc 2018-11-09 10:32:12 +0000
@@ -68,6 +68,9 @@
68 68
69 REAL*8 POL(2)69 REAL*8 POL(2)
70 COMMON/TO_POLARIZATION/ POL70 COMMON/TO_POLARIZATION/ POL
71
72 double precision small_width_treatment
73 common/narrow_width/small_width_treatment
71 74
72 INTEGER ISUM_HEL75 INTEGER ISUM_HEL
73 LOGICAL MULTI_CHANNEL76 LOGICAL MULTI_CHANNEL
@@ -87,6 +90,7 @@
87C ----------90C ----------
88C BEGIN CODE91C BEGIN CODE
89C ----------92C ----------
93
90 NTRY(IMIRROR)=NTRY(IMIRROR)+194 NTRY(IMIRROR)=NTRY(IMIRROR)+1
91 THIS_NTRY(IMIRROR) = THIS_NTRY(IMIRROR)+1 95 THIS_NTRY(IMIRROR) = THIS_NTRY(IMIRROR)+1
92 DO I=1,NEXTERNAL96 DO I=1,NEXTERNAL
@@ -255,6 +259,11 @@
255C Needed for v4 models259C Needed for v4 models
256 COMPLEX*16 DUM0,DUM1260 COMPLEX*16 DUM0,DUM1
257 DATA DUM0, DUM1/(0d0, 0d0), (1d0, 0d0)/261 DATA DUM0, DUM1/(0d0, 0d0), (1d0, 0d0)/
262
263 %(fake_width_declaration)s
264 logical first
265 data first /.true./
266 save first
258C267C
259C FUNCTION268C FUNCTION
260C269C
@@ -265,6 +274,9 @@
265 Double Precision amp2(maxamps), jamp2(0:maxflow)274 Double Precision amp2(maxamps), jamp2(0:maxflow)
266 common/to_amps/ amp2, jamp2275 common/to_amps/ amp2, jamp2
267 include 'coupl.inc'276 include 'coupl.inc'
277
278 double precision small_width_treatment
279 common/narrow_width/small_width_treatment
268C 280C
269C COLOR DATA281C COLOR DATA
270C 282C
@@ -272,6 +284,12 @@
272C ----------284C ----------
273C BEGIN CODE285C BEGIN CODE
274C ----------286C ----------
287if (first) then
288 first=.false.
289 %(fake_width_definitions)s
290endif
291
292
275%(helas_calls)s293%(helas_calls)s
276%(jamp_lines)s294%(jamp_lines)s
277295
278296
=== modified file 'madgraph/iolibs/template_files/matrix_madevent_v4.inc'
--- madgraph/iolibs/template_files/matrix_madevent_v4.inc 2017-11-01 10:08:52 +0000
+++ madgraph/iolibs/template_files/matrix_madevent_v4.inc 2018-11-09 10:32:12 +0000
@@ -49,6 +49,7 @@
49 INTEGER IDUM, NGOOD, IGOOD(NCOMB), JHEL, J, JJ49 INTEGER IDUM, NGOOD, IGOOD(NCOMB), JHEL, J, JJ
50 REAL XRAN150 REAL XRAN1
51 EXTERNAL XRAN151 EXTERNAL XRAN1
52
52C 53C
53C GLOBAL VARIABLES54C GLOBAL VARIABLES
54C 55C
@@ -58,7 +59,7 @@
58 CHARACTER*101 HEL_BUFF59 CHARACTER*101 HEL_BUFF
59 COMMON/TO_HELICITY/ HEL_BUFF60 COMMON/TO_HELICITY/ HEL_BUFF
60 61
61 REAL*8 POL(2)62 REAL*8 POL(2)
62 COMMON/TO_POLARIZATION/ POL63 COMMON/TO_POLARIZATION/ POL
63 64
64 INTEGER ISUM_HEL65 INTEGER ISUM_HEL
@@ -68,6 +69,7 @@
68 DATA IDUM /-1/69 DATA IDUM /-1/
69 DATA XTRY, XREJ, NGOOD /0,0,0/70 DATA XTRY, XREJ, NGOOD /0,0,0/
70 SAVE YFRAC, IGOOD, JHEL71 SAVE YFRAC, IGOOD, JHEL
72
71%(helicity_lines)s73%(helicity_lines)s
72%(den_factor_line)s74%(den_factor_line)s
73C ----------75C ----------
@@ -235,6 +237,11 @@
235C Needed for v4 models237C Needed for v4 models
236 COMPLEX*16 DUM0,DUM1238 COMPLEX*16 DUM0,DUM1
237 DATA DUM0, DUM1/(0d0, 0d0), (1d0, 0d0)/239 DATA DUM0, DUM1/(0d0, 0d0), (1d0, 0d0)/
240
241 %(fake_width_declaration)s
242 logical first
243 data first /.true./
244 save first
238C245C
239C FUNCTION246C FUNCTION
240C247C
@@ -245,6 +252,9 @@
245 Double Precision amp2(maxamps), jamp2(0:maxflow)252 Double Precision amp2(maxamps), jamp2(0:maxflow)
246 common/to_amps/ amp2, jamp2253 common/to_amps/ amp2, jamp2
247 include 'coupl.inc'254 include 'coupl.inc'
255
256 double precision small_width_treatment
257 common/narrow_width/small_width_treatment
248C 258C
249C COLOR DATA259C COLOR DATA
250C 260C
@@ -252,6 +262,11 @@
252C ----------262C ----------
253C BEGIN CODE263C BEGIN CODE
254C ----------264C ----------
265if (first) then
266 first=.false.
267 %(fake_width_definitions)s
268endif
269
255%(helas_calls)s270%(helas_calls)s
256%(jamp_lines)s271%(jamp_lines)s
257272
258273
=== modified file 'madgraph/iolibs/template_files/matrix_madweight_group_v4.inc'
--- madgraph/iolibs/template_files/matrix_madweight_group_v4.inc 2012-02-29 12:53:46 +0000
+++ madgraph/iolibs/template_files/matrix_madweight_group_v4.inc 2018-11-09 10:32:12 +0000
@@ -22,6 +22,7 @@
22C ARGUMENTS 22C ARGUMENTS
23C 23C
24 REAL*8 P(0:3,NEXTERNAL),ANS24 REAL*8 P(0:3,NEXTERNAL),ANS
25
25C 26C
26C LOCAL VARIABLES 27C LOCAL VARIABLES
27C 28C
2829
=== modified file 'madgraph/madevent/sum_html.py'
--- madgraph/madevent/sum_html.py 2018-05-25 09:13:24 +0000
+++ madgraph/madevent/sum_html.py 2018-11-09 10:32:12 +0000
@@ -699,7 +699,7 @@
699 P_comb = Combine_results(Pdir)699 P_comb = Combine_results(Pdir)
700 700
701 if jobs:701 if jobs:
702 for job in filter(lambda j: j['p_dir'] == Pdir, jobs):702 for job in filter(lambda j: j['p_dir'] in Pdir, jobs):
703 P_comb.add_results(os.path.basename(job['dirname']),\703 P_comb.add_results(os.path.basename(job['dirname']),\
704 pjoin(job['dirname'],'results.dat'))704 pjoin(job['dirname'],'results.dat'))
705 elif folder_names:705 elif folder_names:
706706
=== modified file 'madgraph/various/banner.py'
--- madgraph/various/banner.py 2018-05-17 13:55:48 +0000
+++ madgraph/various/banner.py 2018-11-09 10:32:12 +0000
@@ -483,9 +483,9 @@
483 def charge_card(self, tag):483 def charge_card(self, tag):
484 """Build the python object associated to the card"""484 """Build the python object associated to the card"""
485 485
486 if tag == 'param_card':486 if tag in ['param_card', 'param']:
487 tag = 'slha'487 tag = 'slha'
488 elif tag == 'run_card':488 elif tag in ['run_card', 'run']:
489 tag = 'mgruncard' 489 tag = 'mgruncard'
490 elif tag == 'proc_card':490 elif tag == 'proc_card':
491 tag = 'mg5proccard' 491 tag = 'mg5proccard'
@@ -905,7 +905,10 @@
905 out = []905 out = []
906 for line in self:906 for line in self:
907 if line.startswith('define'):907 if line.startswith('define'):
908 name, content = line[7:].split('=',1)908 try:
909 name, content = line[7:].split('=',1)
910 except ValueError:
911 name, content = line[7:].split(None,1)
909 out.append((name, content))912 out.append((name, content))
910 return out 913 return out
911 else:914 else:
@@ -2372,9 +2375,9 @@
2372 else:2375 else:
2373 this_group = this_group[0]2376 this_group = this_group[0]
2374 if block_name in write_block:2377 if block_name in write_block:
2375 text += this_group.on_template % self2378 text += this_group.template_on % self
2376 for name in this_group.fields:2379 for name in this_group.fields:
2377 written.add(f)2380 written.add(name)
2378 if name in to_write:2381 if name in to_write:
2379 to_write.remove(name)2382 to_write.remove(name)
2380 else:2383 else:
@@ -2949,6 +2952,7 @@
2949 self.add_param('survey_splitting', -1, hidden=True, include=False, comment="for loop-induced control how many core are used at survey for the computation of a single iteration.")2952 self.add_param('survey_splitting', -1, hidden=True, include=False, comment="for loop-induced control how many core are used at survey for the computation of a single iteration.")
2950 self.add_param('survey_nchannel_per_job', 2, hidden=True, include=False, comment="control how many Channel are integrated inside a single job on cluster/multicore")2953 self.add_param('survey_nchannel_per_job', 2, hidden=True, include=False, comment="control how many Channel are integrated inside a single job on cluster/multicore")
2951 self.add_param('refine_evt_by_job', -1, hidden=True, include=False, comment="control the maximal number of events for the first iteration of the refine (larger means less jobs)")2954 self.add_param('refine_evt_by_job', -1, hidden=True, include=False, comment="control the maximal number of events for the first iteration of the refine (larger means less jobs)")
2955 self.add_param('small_width_treatment', 1e-6, hidden=True, comment="generation where the width is below VALUE times mass will be replace by VALUE times mass for the computation. The cross-section will be corrected assuming NWA. Not used for loop-induced process")
2952 2956
2953 # parameter allowing to define simple cut via the pdg2957 # parameter allowing to define simple cut via the pdg
2954 # Special syntax are related to those. (can not be edit directly)2958 # Special syntax are related to those. (can not be edit directly)
@@ -2970,7 +2974,7 @@
2970 self.add_param('etamax4pdg',[-1.], system=True) 2974 self.add_param('etamax4pdg',[-1.], system=True)
2971 self.add_param('mxxmin4pdg',[-1.], system=True)2975 self.add_param('mxxmin4pdg',[-1.], system=True)
2972 self.add_param('mxxpart_antipart', [False], system=True)2976 self.add_param('mxxpart_antipart', [False], system=True)
2973 2977
2974 2978
2975 def check_validity(self):2979 def check_validity(self):
2976 """ """2980 """ """
29772981
=== modified file 'madgraph/various/lhe_parser.py'
--- madgraph/various/lhe_parser.py 2018-06-21 08:57:21 +0000
+++ madgraph/various/lhe_parser.py 2018-11-09 10:32:12 +0000
@@ -167,6 +167,7 @@
167 """associate to this particle the decay in the associate event"""167 """associate to this particle the decay in the associate event"""
168 168
169 return self.event.add_decay_to_particle(self.event_id, decay_event)169 return self.event.add_decay_to_particle(self.event_id, decay_event)
170
170 171
171 def __repr__(self):172 def __repr__(self):
172 return 'Particle("%s", event=%s)' % (str(self), self.event)173 return 'Particle("%s", event=%s)' % (str(self), self.event)
@@ -1263,26 +1264,28 @@
1263 1264
1264 if 'part' == status:1265 if 'part' == status:
1265 part = Particle(line, event=self)1266 part = Particle(line, event=self)
1266 if part.E != 0:1267 if part.E != 0 or part.status==-1:
1267 self.append(part)1268 self.append(part)
1268 elif self.nexternal:1269 elif self.nexternal:
1269 self.nexternal-=11270 self.nexternal-=1
1270 else:1271 else:
1271 if '</event>' in line:1272 if '</event>' in line:
1272 line = line.replace('</event>','',1)1273 line = line.replace('</event>','',1)
1273 self.tag += '%s\n' % line1274 self.tag += '%s\n' % line
1274 1275
1275 self.assign_mother()1276 self.assign_mother()
1276 1277
1278
1277 def assign_mother(self):1279 def assign_mother(self):
1278 """convert the number in actual particle"""1280 """convert the number in actual particle"""
1279 #Security if not incoming particle. Define a fake particle and set all particle as 1281 #Security if not incoming particle. Define a fake particle
1280 # decaying from that fake particle
1281 if all(p.status != -1 for p in self):1282 if all(p.status != -1 for p in self):
1283 if not self.nexternal:
1284 return
1282 if self.warning_order:1285 if self.warning_order:
1283 Event.warning_order = False1286 Event.warning_order = False
1284 logger.warning("Weird format for lhe format: no incoming particle... adding a fake one")1287 logger.warning("Weird format for lhe format: no incoming particle... adding a fake one")
1285 1288 raise Exception
1286 mother = Particle(event=self)1289 mother = Particle(event=self)
1287 mother.status = -11290 mother.status = -1
1288 mother.pid = 01291 mother.pid = 0
@@ -1291,8 +1294,10 @@
1291 mother.event_id = 01294 mother.event_id = 0
1292 self.nexternal += 11295 self.nexternal += 1
1293 for p in self[1:]:1296 for p in self[1:]:
1294 p.mother1 = 11297 if isinstance(p.mother1, int) and p.mother1 > 1:
1295 p.mother2 = 11298 p.mother1 += 1
1299 if isinstance(p.mother2, int) and p.mother2 > 1:
1300 p.mother2 += 1
1296 p.event_id += 11301 p.event_id += 1
1297 1302
1298 1303
12991304
=== modified file 'madgraph/various/systematics.py'
--- madgraph/various/systematics.py 2018-03-25 08:30:56 +0000
+++ madgraph/various/systematics.py 2018-11-09 10:32:12 +0000
@@ -885,6 +885,7 @@
885 for onewgt in cevent.wgts:885 for onewgt in cevent.wgts:
886 if not __debug__ and (dyn== -1 and Dmur==1 and Dmuf==1 and pdf==self.orig_pdf):886 if not __debug__ and (dyn== -1 and Dmur==1 and Dmuf==1 and pdf==self.orig_pdf):
887 wgt += onewgt.ref_wgt 887 wgt += onewgt.ref_wgt
888 continue
888 889
889 if dyn == -1:890 if dyn == -1:
890 mur2 = onewgt.scales2[1]891 mur2 = onewgt.scales2[1]
891892
=== modified file 'models/check_param_card.py'
--- models/check_param_card.py 2018-04-26 14:24:04 +0000
+++ models/check_param_card.py 2018-11-09 10:32:12 +0000
@@ -178,7 +178,7 @@
178 178
179 def rename_keys(self, change_keys):179 def rename_keys(self, change_keys):
180 180
181 misc.sprint(self.param_dict, change_keys, [p.lhacode for p in self])181
182 for old_key, new_key in change_keys.items():182 for old_key, new_key in change_keys.items():
183 183
184 assert old_key in self.param_dict184 assert old_key in self.param_dict
@@ -509,15 +509,15 @@
509 if model_value.imag > 1e-5 * model_value.real:509 if model_value.imag > 1e-5 * model_value.real:
510 raise Exception, "Width should be real number: particle %s (%s) has mass: %s" 510 raise Exception, "Width should be real number: particle %s (%s) has mass: %s"
511 model_value = model_value.real511 model_value = model_value.real
512 if not misc.equal(model_value, param_value, 4):512 if not misc.equal(abs(model_value), param_value, 4):
513 modify = True513 modify = True
514 if loglevel == 20:514 if loglevel == 20:
515 logger.info('For consistency, the width of particle %s (%s) is changed to %s.' % (lhacode, particle.get('name'), model_value), '$MG:BOLD')515 logger.info('For consistency, the width of particle %s (%s) is changed to %s.' % (lhacode, particle.get('name'), model_value), '$MG:BOLD')
516 else:516 else:
517 logger.log(loglevel,'For consistency, the width of particle %s (%s) is changed to %s.' % (lhacode, particle.get('name'), model_value))517 logger.log(loglevel,'For consistency, the width of particle %s (%s) is changed to %s.' % (lhacode, particle.get('name'), model_value))
518 #logger.debug('was %s', param_value)518 #logger.debug('was %s', param_value)
519 if model_value != param_value: 519 if abs(model_value) != param_value:
520 self.get('decay').get(abs(particle.get_pdg_code())).value = model_value520 self.get('decay').get(abs(particle.get_pdg_code())).value = abs(model_value)
521521
522 return modify522 return modify
523523
@@ -640,6 +640,11 @@
640 logger.warning('information about \"%s %s" is missing (full block missing) using default value: %s.' %\640 logger.warning('information about \"%s %s" is missing (full block missing) using default value: %s.' %\
641 (block, lhaid, value))641 (block, lhaid, value))
642 value = str(value).lower()642 value = str(value).lower()
643 #special handling for negative mass -> set width negative
644 if block == 'decay':
645 if self['mass'].get(tuple(lhaid)).value < 0:
646 value = '-%s' % value
647
643 fout.writelines(' %s = %s' % (variable, ('%e'%float(value)).replace('e','d')))648 fout.writelines(' %s = %s' % (variable, ('%e'%float(value)).replace('e','d')))
644 if need_mp:649 if need_mp:
645 fout.writelines(' mp__%s = %s_16' % (variable, value))650 fout.writelines(' mp__%s = %s_16' % (variable, value))
@@ -1036,7 +1041,6 @@
1036 data.append(info[k])1041 data.append(info[k])
1037 else:1042 else:
1038 data.append(0.)1043 data.append(0.)
1039 misc.sprint(name, bench, data)
1040 ff.write(formatting % tuple([name] + bench + data))1044 ff.write(formatting % tuple([name] + bench + data))
1041 1045
1042 if not path:1046 if not path:
10431047
=== modified file 'models/import_ufo.py'
--- models/import_ufo.py 2018-06-15 13:28:53 +0000
+++ models/import_ufo.py 2018-11-09 10:32:12 +0000
@@ -112,6 +112,8 @@
112 data = urllib.urlopen(cluster_path)112 data = urllib.urlopen(cluster_path)
113 except Exception:113 except Exception:
114 continue114 continue
115 if data.getcode() != 200:
116 continue
115 break117 break
116 else:118 else:
117 raise MadGraph5Error, '''Model not found locally and Impossible to connect any of us servers.119 raise MadGraph5Error, '''Model not found locally and Impossible to connect any of us servers.
@@ -138,7 +140,8 @@
138 target = None 140 target = None
139 if 'PYTHONPATH' in os.environ and not local_dir:141 if 'PYTHONPATH' in os.environ and not local_dir:
140 for directory in os.environ['PYTHONPATH'].split(':'):142 for directory in os.environ['PYTHONPATH'].split(':'):
141 if 'UFO' in os.path.basename(directory) and os.path.exists(directory):143 if 'UFO' in os.path.basename(directory) and os.path.exists(directory) and\
144 misc.glob('*/couplings.py', path=directory):
142 target= directory 145 target= directory
143 if target is None:146 if target is None:
144 target = pjoin(MG5DIR, 'models') 147 target = pjoin(MG5DIR, 'models')
@@ -646,14 +649,15 @@
646 continue649 continue
647 names = [interaction['lorentz'][i] for i in to_lor[key]]650 names = [interaction['lorentz'][i] for i in to_lor[key]]
648 names.sort()651 names.sort()
649 652 if self.lorentz_info[names[0]].get('structure') == 'external':
653 continue
650 # get name of the new lorentz654 # get name of the new lorentz
651 if tuple(names) in self.lorentz_combine:655 if tuple(names) in self.lorentz_combine:
652 # already created new loretnz656 # already created new loretnz
653 new_name = self.lorentz_combine[tuple(names)]657 new_name = self.lorentz_combine[tuple(names)]
654 else:658 else:
655 new_name = self.add_merge_lorentz(names)659 new_name = self.add_merge_lorentz(names)
656 660
657 # remove the old couplings 661 # remove the old couplings
658 color, coup = key662 color, coup = key
659 to_remove = [(color, lor) for lor in to_lor[key]] 663 to_remove = [(color, lor) for lor in to_lor[key]]
@@ -696,7 +700,12 @@
696 # load the associate lorentz expression700 # load the associate lorentz expression
697 new_struct = ' + '.join([self.lorentz_info[n].get('structure') for n in names])701 new_struct = ' + '.join([self.lorentz_info[n].get('structure') for n in names])
698 spins = self.lorentz_info[names[0]].get('spins')702 spins = self.lorentz_info[names[0]].get('spins')
699 new_lor = self.add_lorentz(new_name, spins, new_struct)703 formfactors = sum([ self.lorentz_info[n].get('formfactors') for n in names \
704 if hasattr(self.lorentz_info[n], 'formfactors') \
705 and self.lorentz_info[n].get('formfactors') \
706 ],[])
707
708 new_lor = self.add_lorentz(new_name, spins, new_struct, formfactors)
700 self.lorentz_info[new_name] = new_lor709 self.lorentz_info[new_name] = new_lor
701 710
702 return new_name711 return new_name
@@ -1368,15 +1377,19 @@
1368 1377
1369 return '' if sign ==1 else '-'1378 return '' if sign ==1 else '-'
13701379
1371 def add_lorentz(self, name, spins , expr):1380 def add_lorentz(self, name, spins , expr, formfact=None):
1372 """ Add a Lorentz expression which is not present in the UFO """1381 """ Add a Lorentz expression which is not present in the UFO """
13731382
1383 logger.debug('MG5 converter defines %s to %s', name, expr)
1374 assert name not in [l.name for l in self.model['lorentz']]1384 assert name not in [l.name for l in self.model['lorentz']]
1375 with misc.TMP_variable(self.ufomodel.object_library, 'all_lorentz', 1385 with misc.TMP_variable(self.ufomodel.object_library, 'all_lorentz',
1376 self.model['lorentz']):1386 self.model['lorentz']):
1377 new = self.model['lorentz'][0].__class__(name = name,1387 new = self.model['lorentz'][0].__class__(name = name,
1378 spins = spins,1388 spins = spins,
1379 structure = expr)1389 structure = expr)
1390 if formfact:
1391 new.formfactors = formfact
1392
1380 assert name in [l.name for l in self.model['lorentz']]1393 assert name in [l.name for l in self.model['lorentz']]
1381 assert name not in [l.name for l in self.ufomodel.all_lorentz]1394 assert name not in [l.name for l in self.ufomodel.all_lorentz]
1382 #self.model['lorentz'].append(new) # already done by above command1395 #self.model['lorentz'].append(new) # already done by above command
@@ -1973,7 +1986,7 @@
1973 null_parameters.append(name)1986 null_parameters.append(name)
1974 elif value == 1:1987 elif value == 1:
1975 one_parameters.append(name)1988 one_parameters.append(name)
1976 1989
1977 return null_parameters, one_parameters1990 return null_parameters, one_parameters
1978 1991
1979 def apply_conditional_simplifications(self, modified_params,1992 def apply_conditional_simplifications(self, modified_params,
@@ -2301,7 +2314,9 @@
2301 particle['width'] = 'ZERO'2314 particle['width'] = 'ZERO'
2302 if particle['width'] in one_parameters:2315 if particle['width'] in one_parameters:
2303 one_parameters.remove(particle['width']) 2316 one_parameters.remove(particle['width'])
2304 2317 if particle['mass'] in one_parameters:
2318 one_parameters.remove(particle['mass'])
2319
2305 for pdg, particle in self['particle_dict'].items():2320 for pdg, particle in self['particle_dict'].items():
2306 if particle['mass'] in zero_parameters:2321 if particle['mass'] in zero_parameters:
2307 particle['mass'] = 'ZERO'2322 particle['mass'] = 'ZERO'
@@ -2340,6 +2355,13 @@
2340 for coupling in coupling_list:2355 for coupling in coupling_list:
2341 for use in re_pat.findall(coupling.expr):2356 for use in re_pat.findall(coupling.expr):
2342 used.add(use)2357 used.add(use)
2358
2359 # check in form-factor
2360 for lor in self['lorentz']:
2361 if hasattr(lor, 'formfactors') and lor.formfactors:
2362 for ff in lor.formfactors:
2363 for use in re_pat.findall(ff.value):
2364 used.add(use)
2343 else:2365 else:
2344 used = set([i for i in special_parameters if i])2366 used = set([i for i in special_parameters if i])
2345 2367
@@ -2484,16 +2506,26 @@
2484 if any( n.startswith('d') for n in names ):2506 if any( n.startswith('d') for n in names ):
2485 new_struct += '-' + ' - '.join(['1.*(%s)' %self.lorentz_info[n[1:]].get('structure') for n in names if n.startswith('d')])2507 new_struct += '-' + ' - '.join(['1.*(%s)' %self.lorentz_info[n[1:]].get('structure') for n in names if n.startswith('d')])
2486 spins = self.lorentz_info[names[0][1:]].get('spins')2508 spins = self.lorentz_info[names[0][1:]].get('spins')
2487 new_lor = self.add_lorentz(new_name, spins, new_struct)2509 formfact = sum([ self.lorentz_info[n[1:]].get('formfactors') for n in names \
2510 if hasattr(self.lorentz_info[n[1:]], 'formfactors') \
2511 and self.lorentz_info[n[1:]].get('formfactors') \
2512 ],[])
2513
2514
2515
2516
2517 new_lor = self.add_lorentz(new_name, spins, new_struct, formfact)
2488 self.lorentz_info[new_name] = new_lor2518 self.lorentz_info[new_name] = new_lor
2489 2519
2490 return new_name2520 return new_name
2491 2521
2492 def add_lorentz(self, name, spin, struct):2522 def add_lorentz(self, name, spin, struct, formfact=None):
2493 """adding lorentz structure to the current model"""2523 """adding lorentz structure to the current model"""
2494 new = self['lorentz'][0].__class__(name = name,2524 new = self['lorentz'][0].__class__(name = name,
2495 spins = spin,2525 spins = spin,
2496 structure = struct)2526 structure = struct)
2527 if formfact:
2528 new.formfactors = formfact
2497 self['lorentz'].append(new)2529 self['lorentz'].append(new)
2498 self.create_lorentz_dict()2530 self.create_lorentz_dict()
2499 2531
25002532
=== modified file 'tests/IOTests.py'
--- tests/IOTests.py 2016-09-07 13:39:06 +0000
+++ tests/IOTests.py 2018-11-09 10:32:12 +0000
@@ -296,20 +296,26 @@
296 """Clean up the file created. Called at the end of the test run."""296 """Clean up the file created. Called at the end of the test run."""
297 297
298 pathsToClean = [self.temporary_folder]298 pathsToClean = [self.temporary_folder]
299 299
300 if not self.clean_function is None:300 if not self.clean_function is None:
301 paths, prevent_cleanUp = self.clean_function(*args, **kwargs)301 paths, prevent_cleanUp = self.clean_function(*args, **kwargs)
302 pathsToClean.extend(paths)302 if isinstance(paths, str):
303 pathsToClean.append(paths)
304 else:
305 pathsToClean.extend(paths)
303 if prevent_cleanUp:306 if prevent_cleanUp:
304 print colored%(31,307 print colored%(31,
305 "Clean up of the following of temporary folders prevented:")308 "Clean up of the following of temporary folders prevented:")
306 for path in pathsToClean:309 for path in pathsToClean:
307 print colored%(31," > %s"%str(path))310 print colored%(31," > %s"%str(path))
308311
309 try:312 try:
310 for path in pathsToClean:313 for path in pathsToClean:
311 shutil.rmtree(path)314 if os.path.isdir(path):
312 except OSError:315 shutil.rmtree(path)
316 else:
317 os.remove(path)
318 except OSError,error:
313 pass 319 pass
314320
315#===============================================================================321#===============================================================================
316322
=== modified file 'tests/acceptance_tests/test_cmd_madloop.py'
--- tests/acceptance_tests/test_cmd_madloop.py 2018-04-24 22:39:50 +0000
+++ tests/acceptance_tests/test_cmd_madloop.py 2018-11-09 10:32:12 +0000
@@ -492,7 +492,7 @@
492 # Select the Tensor Integral to include in the test492 # Select the Tensor Integral to include in the test
493 misc.deactivate_dependence('pjfry', cmd = self.interface, log='stdout')493 misc.deactivate_dependence('pjfry', cmd = self.interface, log='stdout')
494 misc.deactivate_dependence('samurai', cmd = self.interface, log='stdout') 494 misc.deactivate_dependence('samurai', cmd = self.interface, log='stdout')
495 misc.activate_dependence('golem', cmd = self.interface, log='stdout')495 misc.deactivate_dependence('golem', cmd = self.interface, log='stdout')
496 misc.activate_dependence('ninja', cmd = self.interface, log='stdout',MG5dir=MG5DIR)496 misc.activate_dependence('ninja', cmd = self.interface, log='stdout',MG5dir=MG5DIR)
497 misc.activate_dependence('collier', cmd = self.interface, log='stdout',MG5dir=MG5DIR)497 misc.activate_dependence('collier', cmd = self.interface, log='stdout',MG5dir=MG5DIR)
498498
@@ -533,7 +533,7 @@
533 # Select the Tensor Integral to include in the test533 # Select the Tensor Integral to include in the test
534 misc.deactivate_dependence('pjfry', cmd = interface, log='stdout')534 misc.deactivate_dependence('pjfry', cmd = interface, log='stdout')
535 misc.deactivate_dependence('samurai', cmd = interface, log='stdout') 535 misc.deactivate_dependence('samurai', cmd = interface, log='stdout')
536 misc.activate_dependence('golem', cmd = interface, log='stdout')536 misc.deactivate_dependence('golem', cmd = interface, log='stdout')
537 misc.activate_dependence('ninja', cmd = interface, log='stdout',MG5dir=MG5DIR)537 misc.activate_dependence('ninja', cmd = interface, log='stdout',MG5dir=MG5DIR)
538538
539 run_cmd('generate g g > t t~ [virt=QCD]')539 run_cmd('generate g g > t t~ [virt=QCD]')
540540
=== modified file 'tests/acceptance_tests/test_export_fks.py'
--- tests/acceptance_tests/test_export_fks.py 2017-05-15 10:39:31 +0000
+++ tests/acceptance_tests/test_export_fks.py 2018-11-09 10:32:12 +0000
@@ -75,7 +75,7 @@
75 # Select the Tensor Integral to include in the test75 # Select the Tensor Integral to include in the test
76 misc.deactivate_dependence('pjfry', cmd = interface, log='stdout')76 misc.deactivate_dependence('pjfry', cmd = interface, log='stdout')
77 misc.deactivate_dependence('samurai', cmd = interface, log='stdout') 77 misc.deactivate_dependence('samurai', cmd = interface, log='stdout')
78 misc.activate_dependence('golem', cmd = interface, log='stdout')78 misc.deactivate_dependence('golem', cmd = interface, log='stdout')
79 misc.activate_dependence('ninja', cmd = interface, log='stdout',MG5dir=MG5DIR) 79 misc.activate_dependence('ninja', cmd = interface, log='stdout',MG5dir=MG5DIR)
80 80
81 run_cmd('import model %s' % model)81 run_cmd('import model %s' % model)
8282
=== modified file 'tests/acceptance_tests/test_madspin.py'
--- tests/acceptance_tests/test_madspin.py 2018-05-25 09:13:24 +0000
+++ tests/acceptance_tests/test_madspin.py 2018-11-09 10:32:12 +0000
@@ -93,7 +93,7 @@
93 subprocess.call([pjoin(MG5DIR, 'MadSpin', 'madspin'),93 subprocess.call([pjoin(MG5DIR, 'MadSpin', 'madspin'),
94 pjoin(self.path, 'test_hepmc')],94 pjoin(self.path, 'test_hepmc')],
95 cwd=pjoin(self.path),95 cwd=pjoin(self.path),
96 stdout=stdout,stderr=stdout)96 stdout=stdout,stderr=stderr)
97 self.assertTrue(os.path.exists(pjoin(self.path, 'test_decayed.lhe.gz')))97 self.assertTrue(os.path.exists(pjoin(self.path, 'test_decayed.lhe.gz')))
98 lhe = lhe_parser.EventFile(pjoin(self.path, 'test_decayed.lhe.gz'))98 lhe = lhe_parser.EventFile(pjoin(self.path, 'test_decayed.lhe.gz'))
99 self.assertEqual(10, len(lhe))99 self.assertEqual(10, len(lhe))
@@ -111,4 +111,56 @@
111 111
112 self.assertEqual(nb_dec, 116)112 self.assertEqual(nb_dec, 116)
113 self.assertEqual(nb_photon, 116)113 self.assertEqual(nb_photon, 116)
114
115\ No newline at end of file114\ No newline at end of file
115
116 def test_lhe_none_decay(self):
117 """ """
118
119 cwd = os.getcwd()
120
121 files.cp(pjoin(MG5DIR, 'tests', 'input_files', 'test_spinmode_none.lhe.gz'), self.path)
122
123
124 fsock = open(pjoin(self.path, 'test_hepmc'),'w')
125 text = """
126 set spinmode none
127 import ./test_spinmode_none.lhe.gz
128 decay z > mu+ mu-
129 launch
130 """
131
132 fsock.write(text)
133 fsock.close()
134
135 import subprocess
136 if logging.getLogger('madgraph').level <= 20:
137 stdout=None
138 stderr=None
139 else:
140 devnull =open(os.devnull,'w')
141 stdout=devnull
142 stderr=devnull
143
144 subprocess.call([pjoin(MG5DIR, 'MadSpin', 'madspin'),
145 pjoin(self.path, 'test_hepmc')],
146 cwd=pjoin(self.path),
147 stdout=stdout,stderr=stderr)
148
149 self.assertTrue(os.path.exists(pjoin(self.path, 'test_spinmode_none_decayed.lhe.gz')))
150 lhe = lhe_parser.EventFile(pjoin(self.path, 'test_spinmode_none_decayed.lhe.gz'))
151 self.assertEqual(100, len(lhe))
152
153 nb_dec = 0
154 nb_muon = 0
155 for event in lhe:
156 muon_in = 0
157 self.assertEqual(event.nexternal, len(event))
158 for particle in event:
159 if particle.pdg == 23:
160 self.assertEqual(particle.status,2)
161 nb_dec += 1
162 if particle.pdg == 13:
163 nb_muon += 1
164 muon_in +=1
165 self.assertEqual(muon_in, 1)
166 self.assertEqual(nb_dec, 189)
167 self.assertEqual(nb_muon, 100)
116\ No newline at end of file168\ No newline at end of file
117169
=== modified file 'tests/acceptance_tests/test_model_equivalence.py'
--- tests/acceptance_tests/test_model_equivalence.py 2018-06-15 08:18:31 +0000
+++ tests/acceptance_tests/test_model_equivalence.py 2018-11-09 10:32:12 +0000
@@ -146,15 +146,14 @@
146 ufo_model.pass_particles_name_in_mg_default()146 ufo_model.pass_particles_name_in_mg_default()
147 147
148 # import MG4 model148 # import MG4 model
149
149 model = base_objects.Model()150 model = base_objects.Model()
150 if not MG4DIR:151 if not MG4DIR:
151 raise MadGraph5Error, "Please provide a valid MG/ME path with -d"152 raise MadGraph5Error, "Please provide a valid MG/ME path with -d"
152 v4_path = os.path.join(MG4DIR, 'models', 'mssm_v4')153 v4_path = os.path.join(MG4DIR, 'models', 'mssm_v4')
153 if not os.path.isdir(v4_path):154 if not os.path.isdir(v4_path):
154 v4_path = os.path.join(MG4DIR, 'Models', 'mssm')155 import_ufo.import_model_from_db('mssm_v4', local_dir=True)
155 if not os.path.isdir(v4_path):156
156 raise MadGraph5Error, \
157 "Please provide a valid MG/ME path with -d"
158157
159 model.set('particles', files.read_from_file(158 model.set('particles', files.read_from_file(
160 os.path.join(v4_path,'particles.dat'),159 os.path.join(v4_path,'particles.dat'),
161160
=== removed file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%GOLEM_interface.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%GOLEM_interface.f 2016-08-07 07:16:07 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%GOLEM_interface.f 1970-01-01 00:00:00 +0000
@@ -1,748 +0,0 @@
1 SUBROUTINE GOLEMLOOP(NLOOPLINE,PL,M2L,RANK,RES,STABLE)
2C
3C Generated by MadGraph5_aMC@NLO v. %(version)s, %(date)s
4C By the MadGraph5_aMC@NLO Development Team
5C Visit launchpad.net/madgraph5 and amcatnlo.web.cern.ch
6C
7C Interface between MG5 and Golem95.
8C The Golem95 version should be higher than 1.3.0.
9C It supports RANK = NLOOPLINE + 1 tensor integrals when 1 <
10C NLOOPLINE < 6.
11C
12C Process: d~ u > w+ QED<=1 WEIGHTED<=2 [ all = QCD ]
13C Process: s~ c > w+ QED<=1 WEIGHTED<=2 [ all = QCD ]
14C
15C
16C MODULES
17C
18 USE MATRICE_S
19 USE FORM_FACTOR_TYPE, ONLY: FORM_FACTOR
20 USE PRECISION_GOLEM, ONLY: KI
21 USE TENS_COMB
22 USE TENS_REC
23 USE FORM_FACTOR_1P, ONLY: A10
24 USE FORM_FACTOR_2P, ONLY: A20
25 USE FORM_FACTOR_3P, ONLY: A30
26 USE FORM_FACTOR_4P, ONLY: A40
27 USE FORM_FACTOR_5P, ONLY: A50
28 USE FORM_FACTOR_6P, ONLY: A60
29 IMPLICIT NONE
30C
31C CONSTANTS
32C
33 INTEGER NEXTERNAL
34 PARAMETER (NEXTERNAL=3)
35 LOGICAL CHECKPCONSERVATION
36 PARAMETER (CHECKPCONSERVATION=.TRUE.)
37 REAL*8 NORMALIZATION
38 PARAMETER (NORMALIZATION = 1.D0/(16.D0*3.14159265358979323846D0*
39 $ *2))
40 REAL(KI),DIMENSION(0:3),PARAMETER::NULL_VEC = (/0.0_KI,0.0_KI
41 $ ,0.0_KI,0.0_KI/)
42C GOLEM_RUN_MODE = 1: Use directly MadLoop tensorial coefficients
43C GOLEM_RUN_MODE = 2: Reconstruct the tensorial coefficeints
44C directly from
45C numerator using golem internal reconstruction routine
46C GOLEM_RUN_MODE = 3: Cross-checked reconstructed coefficients
47C against
48C MadLoop internal ones.
49 INTEGER GOLEM_RUN_MODE
50 PARAMETER (GOLEM_RUN_MODE=1)
51C The following is the acceptance threshold used for
52C GOLEM_RUN_MODE = 3
53 REAL*8 COEF_CHECK_THRS
54 DATA COEF_CHECK_THRS/1.0D-13/
55 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
56
57 LOGICAL PASS_COEF_CHECK
58C
59C ARGUMENTS
60C
61 INTEGER NLOOPLINE, RANK
62 REAL*8 PL(0:3,NLOOPLINE)
63 REAL*8 PCT(0:3,0:NLOOPLINE-1), ABSPCT(0:3)
64 REAL*8 REF_P
65 REAL(KI) PGOLEM(NLOOPLINE,0:3)
66 COMPLEX*16 M2L(NLOOPLINE)
67 COMPLEX(KI) M2LGOLEM(NLOOPLINE)
68 COMPLEX*16 RES(3)
69 LOGICAL STABLE
70C
71C LOCAL VARIABLES
72C
73 INTEGER I, J, K
74 TYPE(FORM_FACTOR)::RES_GOLEM
75
76 COMPLEX(KI)::COEFFS0,COEFFS0_REC
77 TYPE(COEFF_TYPE_1)::COEFFS1,COEFFS1_REC
78 TYPE(COEFF_TYPE_2)::COEFFS2,COEFFS2_REC
79 TYPE(COEFF_TYPE_3)::COEFFS3,COEFFS3_REC
80 TYPE(COEFF_TYPE_4)::COEFFS4,COEFFS4_REC
81 TYPE(COEFF_TYPE_5)::COEFFS5,COEFFS5_REC
82 TYPE(COEFF_TYPE_6)::COEFFS6,COEFFS6_REC
83
84C The pinch propagator optimization is not used, so for now it is
85C always 0.
86 INTEGER PINCH
87C
88C EXTERNAL FUNCTIONS
89C
90 COMPLEX(KI) GOLEM_LOOPNUM
91 EXTERNAL GOLEM_LOOPNUM
92 LOGICAL COMPARE_COEFS_0
93 LOGICAL COMPARE_COEFS_1
94 LOGICAL COMPARE_COEFS_2
95 LOGICAL COMPARE_COEFS_3
96 LOGICAL COMPARE_COEFS_4
97 LOGICAL COMPARE_COEFS_5
98 LOGICAL COMPARE_COEFS_6
99C
100C GLOBAL VARIABLES
101C
102 INCLUDE 'coupl.inc'
103 INTEGER CTMODE
104 REAL*8 LSCALE
105 COMMON/CT/LSCALE,CTMODE
106
107 INTEGER ID,SQSOINDEX,R
108 COMMON/LOOP/ID,SQSOINDEX,R
109
110 LOGICAL CTINIT, TIRINIT, GOLEMINIT, SAMURAIINIT, NINJAINIT
111 $ ,COLLIERINIT
112 COMMON/REDUCTIONCODEINIT/CTINIT, TIRINIT,GOLEMINIT,SAMURAIINIT
113 $ ,NINJAINIT,COLLIERINIT
114
115 INTEGER NLOOPGROUPS
116 PARAMETER (NLOOPGROUPS=1)
117 INTEGER NSQUAREDSO
118 PARAMETER (NSQUAREDSO=1)
119 INCLUDE 'loop_max_coefs.inc'
120
121 COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
122 COMMON/LCOEFS/LOOPCOEFS
123C ----------
124C BEGIN CODE
125C ----------
126
127C The CT initialization is also performed here if not done already
128C because it calls MPINIT of OneLOop which is necessary on some
129C system
130 IF (CTINIT) THEN
131 CTINIT=.FALSE.
132 CALL INITCT()
133 ENDIF
134
135C INITIALIZE GOLEM IF NEEDED
136 IF (GOLEMINIT) THEN
137 GOLEMINIT=.FALSE.
138 CALL INITGOLEM()
139 ENDIF
140
141C No stability test intrisic to Golem95 now
142 STABLE=.TRUE.
143
144C This initialization must be done for each reduction because we
145C have not setup anyoptimization using pinched propagators yet.
146 CALL INITGOLEM95(NLOOPLINE)
147 PINCH = 0
148
149C YOU CAN FIND THE DETAILS ABOUT THE DIFFERENT CTMODE AT THE
150C BEGINNING OF THE FILE CTS_CUTS.F90 IN THE CUTTOOLS DISTRIBUTION
151
152C CONVERT THE MASSES TO BE COMPLEX
153 DO I=1,NLOOPLINE
154 M2LGOLEM(I)=M2L(I)
155 ENDDO
156
157C CONVERT THE MOMENTA FLOWING IN THE LOOP LINES TO CT CONVENTIONS
158 DO I=0,3
159 ABSPCT(I) = 0.D0
160 DO J=0,(NLOOPLINE-1)
161 PCT(I,J)=0.D0
162 ENDDO
163 ENDDO
164 DO I=0,3
165 DO J=1,NLOOPLINE
166 PCT(I,0)=PCT(I,0)+PL(I,J)
167 ABSPCT(I)=ABSPCT(I)+ABS(PL(I,J))
168 ENDDO
169 ENDDO
170 REF_P = MAX(ABSPCT(0), ABSPCT(1),ABSPCT(2),ABSPCT(3))
171 DO I=0,3
172 ABSPCT(I) = MAX(REF_P*1E-6, ABSPCT(I))
173 ENDDO
174 IF (CHECKPCONSERVATION.AND.REF_P.GT.1D-8) THEN
175 IF ((PCT(0,0)/ABSPCT(0)).GT.1.D-6) THEN
176 WRITE(*,*) 'energy is not conserved ',PCT(0,0)
177 STOP 'energy is not conserved'
178 ELSEIF ((PCT(1,0)/ABSPCT(1)).GT.1.D-6) THEN
179 WRITE(*,*) 'px is not conserved ',PCT(1,0)
180 STOP 'px is not conserved'
181 ELSEIF ((PCT(2,0)/ABSPCT(2)).GT.1.D-6) THEN
182 WRITE(*,*) 'py is not conserved ',PCT(2,0)
183 STOP 'py is not conserved'
184 ELSEIF ((PCT(3,0)/ABSPCT(3)).GT.1.D-6) THEN
185 WRITE(*,*) 'pz is not conserved ',PCT(3,0)
186 STOP 'pz is not conserved'
187 ENDIF
188 ENDIF
189 DO I=0,3
190 DO J=1,(NLOOPLINE-1)
191 DO K=1,J
192 PCT(I,J)=PCT(I,J)+PL(I,K)
193 ENDDO
194 ENDDO
195 ENDDO
196
197C Now convert the loop momenta to Golem95 conventions
198 DO I=0,3
199 PGOLEM(1,I)=0.0E0_KI
200 DO J=2,NLOOPLINE
201 PGOLEM(J,I)=PCT(I,J-1)
202 ENDDO
203 ENDDO
204
205C Fill in the kinematic s-matrix while taking care of on-shell
206C limits.
207 CALL SETUP_KIN_MATRIX(NLOOPLINE,PGOLEM,M2LGOLEM)
208C Construct the golem internal matrices derived from the kinetic
209C one.
210 CALL PREPARESMATRIX()
211
212C Fill in the golem coefficents and compute the loop
213 IF(GOLEM_RUN_MODE.EQ.2)THEN
214 RES_GOLEM = EVALUATE_B(GOLEM_LOOPNUM,PGOLEM,0,RANK)
215 ELSE
216 PASS_COEF_CHECK=.TRUE.
217 SELECT CASE(RANK)
218 CASE(0)
219 CALL FILL_GOLEM_COEFFS_0(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS0)
220 IF(GOLEM_RUN_MODE.EQ.3)THEN
221 COEFFS0_REC = GOLEM_LOOPNUM(NULL_VEC,0.0_KI)
222 PASS_COEF_CHECK=COMPARE_COEFS_0(COEFFS0,COEFFS0_REC)
223 ENDIF
224 CASE(1)
225 CALL FILL_GOLEM_COEFFS_1(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS1)
226 IF(GOLEM_RUN_MODE.EQ.3)THEN
227 CALL RECONSTRUCT1(GOLEM_LOOPNUM,COEFFS1_REC)
228 PASS_COEF_CHECK=COMPARE_COEFS_1(COEFFS1,COEFFS1_REC)
229 ENDIF
230 CASE(2)
231 CALL FILL_GOLEM_COEFFS_2(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS2)
232 IF(GOLEM_RUN_MODE.EQ.3)THEN
233 CALL RECONSTRUCT2(GOLEM_LOOPNUM,COEFFS2_REC)
234 PASS_COEF_CHECK=COMPARE_COEFS_2(COEFFS2,COEFFS2_REC)
235 ENDIF
236 CASE(3)
237 CALL FILL_GOLEM_COEFFS_3(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS3)
238 IF(GOLEM_RUN_MODE.EQ.3)THEN
239 CALL RECONSTRUCT3(GOLEM_LOOPNUM,COEFFS3_REC)
240 PASS_COEF_CHECK=COMPARE_COEFS_3(COEFFS3,COEFFS3_REC)
241 ENDIF
242 CASE(4)
243 CALL FILL_GOLEM_COEFFS_4(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS4)
244 IF(GOLEM_RUN_MODE.EQ.3)THEN
245 CALL RECONSTRUCT4(GOLEM_LOOPNUM,COEFFS4_REC)
246 PASS_COEF_CHECK=COMPARE_COEFS_4(COEFFS4,COEFFS4_REC)
247 ENDIF
248 CASE(5)
249 CALL FILL_GOLEM_COEFFS_5(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS5)
250 IF(GOLEM_RUN_MODE.EQ.3)THEN
251 CALL RECONSTRUCT5(GOLEM_LOOPNUM,COEFFS5_REC)
252 PASS_COEF_CHECK=COMPARE_COEFS_5(COEFFS5,COEFFS5_REC)
253 ENDIF
254 CASE(6)
255 CALL FILL_GOLEM_COEFFS_6(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS6)
256 IF(GOLEM_RUN_MODE.EQ.3)THEN
257 CALL RECONSTRUCT6(GOLEM_LOOPNUM,COEFFS6_REC)
258 PASS_COEF_CHECK=COMPARE_COEFS_6(COEFFS6,COEFFS6_REC)
259 ENDIF
260 CASE DEFAULT
261 WRITE(*,*)'Not yet implemented in Golem95 for rank= ',RANK
262 STOP
263 END SELECT
264
265 IF(.NOT.PASS_COEF_CHECK)THEN
266 WRITE(*,*)'Coefs mismatch for ID ',ID,' and rank ',RANK
267 WRITE(*,*)'Coefs form MadLoop5:'
268 SELECT CASE(RANK)
269 CASE(0)
270 WRITE(*,*)'Constant coef = ',COEFFS0
271 CASE(1)
272 CALL PRINT_COEFFS(COEFFS1)
273 CASE(2)
274 CALL PRINT_COEFFS(COEFFS2)
275 CASE(3)
276 CALL PRINT_COEFFS(COEFFS3)
277 CASE(4)
278 CALL PRINT_COEFFS(COEFFS4)
279 CASE(5)
280 CALL PRINT_COEFFS(COEFFS5)
281 CASE(6)
282 CALL PRINT_COEFFS(COEFFS6)
283 END SELECT
284 WRITE(*,*)'Coefs reconstructed by Golem95:'
285 SELECT CASE(RANK)
286 CASE(0)
287 WRITE(*,*)'Constant coef = ',COEFFS0_REC
288 CASE(1)
289 CALL PRINT_COEFFS(COEFFS1_REC)
290 CASE(2)
291 CALL PRINT_COEFFS(COEFFS2_REC)
292 CASE(3)
293 CALL PRINT_COEFFS(COEFFS3_REC)
294 CASE(4)
295 CALL PRINT_COEFFS(COEFFS4_REC)
296 CASE(5)
297 CALL PRINT_COEFFS(COEFFS5_REC)
298 CASE(6)
299 CALL PRINT_COEFFS(COEFFS6_REC)
300 END SELECT
301 STOP
302 ENDIF
303
304 SELECT CASE(NLOOPLINE)
305 CASE(1)
306 WRITE(*,*)'Golem95 cannot handle with tadpole yet'
307 STOP
308 CASE(2)
309 SELECT CASE(RANK)
310 CASE(0)
311 RES_GOLEM = COEFFS0*A20(PINCH)
312 CASE(1)
313 RES_GOLEM = CONTRACT2_1(COEFFS1,PGOLEM,PINCH)
314 CASE(2)
315 RES_GOLEM = CONTRACT2_2(COEFFS2,PGOLEM,PINCH)
316 CASE(3)
317 RES_GOLEM = CONTRACT2_3(COEFFS3,PGOLEM,PINCH)
318 CASE DEFAULT
319 WRITE(*,*)'Golem95 cannot handle with: N,r = ',2,RANK
320 STOP
321 END SELECT
322 CASE(3)
323 SELECT CASE(RANK)
324 CASE(0)
325 RES_GOLEM = COEFFS0*A30(PINCH)
326 CASE(1)
327 RES_GOLEM = CONTRACT3_1(COEFFS1,PGOLEM,PINCH)
328 CASE(2)
329 RES_GOLEM = CONTRACT3_2(COEFFS2,PGOLEM,PINCH)
330 CASE(3)
331 RES_GOLEM = CONTRACT3_3(COEFFS3,PGOLEM,PINCH)
332 CASE(4)
333 RES_GOLEM = CONTRACT3_4(COEFFS4,PGOLEM,PINCH)
334 CASE DEFAULT
335 WRITE(*,*)'Golem95 cannot handle with: N,r = ',3,RANK
336 STOP
337 END SELECT
338 CASE(4)
339 SELECT CASE(RANK)
340 CASE(0)
341 RES_GOLEM = COEFFS0*A40(PINCH)
342 CASE(1)
343 RES_GOLEM = CONTRACT4_1(COEFFS1,PGOLEM,PINCH)
344 CASE(2)
345 RES_GOLEM = CONTRACT4_2(COEFFS2,PGOLEM,PINCH)
346 CASE(3)
347 RES_GOLEM = CONTRACT4_3(COEFFS3,PGOLEM,PINCH)
348 CASE(4)
349 RES_GOLEM = CONTRACT4_4(COEFFS4,PGOLEM,PINCH)
350 CASE(5)
351 RES_GOLEM = CONTRACT4_5(COEFFS5,PGOLEM,PINCH)
352 CASE DEFAULT
353 WRITE(*,*)'Golem95 cannot handle with: N,r = ',4,RANK
354 STOP
355 END SELECT
356 CASE(5)
357 SELECT CASE(RANK)
358 CASE(0)
359 RES_GOLEM = COEFFS0*A50(PINCH)
360 CASE(1)
361 RES_GOLEM = CONTRACT5_1(COEFFS1,PGOLEM,PINCH)
362 CASE(2)
363 RES_GOLEM = CONTRACT5_2(COEFFS2,PGOLEM,PINCH)
364 CASE(3)
365 RES_GOLEM = CONTRACT5_3(COEFFS3,PGOLEM,PINCH)
366 CASE(4)
367 RES_GOLEM = CONTRACT5_4(COEFFS4,PGOLEM,PINCH)
368 CASE(5)
369 RES_GOLEM = CONTRACT5_5(COEFFS5,PGOLEM,PINCH)
370 CASE(6)
371 RES_GOLEM = CONTRACT5_6(COEFFS6,PGOLEM,PINCH)
372 CASE DEFAULT
373 WRITE(*,*)'Golem95 cannot handle with: N,r = ',5,RANK
374 STOP
375 END SELECT
376 CASE(6)
377 SELECT CASE(RANK)
378 CASE(0)
379 RES_GOLEM = COEFFS0*A60(PINCH)
380 CASE(1)
381 RES_GOLEM = CONTRACT6_1(COEFFS1,PGOLEM,PINCH)
382 CASE(2)
383 RES_GOLEM = CONTRACT6_2(COEFFS2,PGOLEM,PINCH)
384 CASE(3)
385 RES_GOLEM = CONTRACT6_3(COEFFS3,PGOLEM,PINCH)
386 CASE(4)
387 RES_GOLEM = CONTRACT6_4(COEFFS4,PGOLEM,PINCH)
388 CASE(5)
389 RES_GOLEM = CONTRACT6_5(COEFFS5,PGOLEM,PINCH)
390 CASE(6)
391 RES_GOLEM = CONTRACT6_6(COEFFS6,PGOLEM,PINCH)
392 CASE DEFAULT
393 WRITE(*,*)'Golem95 cannot handle with: N,r = ',6,RANK
394 STOP
395 END SELECT
396 CASE DEFAULT
397 WRITE(*,*)'Golem95 cannot handle with: N = ',NLOOPLINE
398 STOP
399 END SELECT
400 ENDIF
401
402 RES(1)=NORMALIZATION*2.0D0*DBLE(RES_GOLEM%%C+2.0*LOG(MU_R)
403 $ *RES_GOLEM%%B+2.0*LOG(MU_R)**2*RES_GOLEM%%A)
404 RES(2)=NORMALIZATION*2.0D0*DBLE(RES_GOLEM%%B+2.0*LOG(MU_R)
405 $ *RES_GOLEM%%A)
406 RES(3)=NORMALIZATION*2.0D0*DBLE(RES_GOLEM%%A)
407C WRITE(*,*) 'Loop ID',ID,' =',RES(1),RES(2),RES(3)
408
409C Finally free golem memory and cache
410 CALL EXITGOLEM95()
411
412 END
413
414 FUNCTION GOLEM_LOOPNUM(Q,MU2)
415 USE PRECISION_GOLEM, ONLY: KI
416 REAL(KI),DIMENSION(0:3),INTENT(IN)::Q
417 REAL(KI),INTENT(IN)::MU2
418 COMPLEX(KI)::GOLEM_LOOPNUM
419
420 COMPLEX*16 QQ(0:3),NUM
421 INTEGER I
422
423 DO I=0,3
424 QQ(I)=CMPLX(Q(I),0.0D0,KIND=16)
425 ENDDO
426
427 CALL LOOPNUM(QQ,NUM)
428 GOLEM_LOOPNUM=NUM
429 RETURN
430 END FUNCTION
431
432 SUBROUTINE INITGOLEM()
433C
434C INITIALISATION OF GOLEM
435C
436C
437C MODULE
438C
439 USE PARAMETRE
440C
441C LOCAL VARIABLES
442C
443 REAL*8 THRS
444 LOGICAL EXT_NUM_FOR_R1
445C
446C GLOBAL VARIABLES
447C
448 INCLUDE 'MadLoopParams.inc'
449C ----------
450C BEGIN CODE
451C ----------
452
453C DEFAULT PARAMETERS FOR GOLEM
454C -------------------------------
455C One can chose here to have either just the rational R1 piece
456C or everything but the R2
457 RAT_OR_TOT_PAR = TOT
458
459 END
460
461 SUBROUTINE SETUP_KIN_MATRIX(NLOOPLINE,PGOLEM,M2L)
462C
463C MODULE
464C
465 USE MATRICE_S
466 USE PRECISION_GOLEM, ONLY: KI
467C
468C ARGUMENTS
469C
470 INTEGER NLOOPLINE
471 REAL(KI) PGOLEM(NLOOPLINE,0:3)
472 COMPLEX(KI) M2L(NLOOPLINE)
473C
474C LOCAL VARIABLES
475C
476 INTEGER I,J
477 COMPLEX*16 S_MAT_FROM_MG(NLOOPLINE,NLOOPLINE)
478C ----------
479C BEGIN CODE
480C ----------
481
482 CALL BUILD_KINEMATIC_MATRIX(NLOOPLINE,PGOLEM,M2L,S_MAT_FROM_MG)
483
484 DO I=1,NLOOPLINE
485 DO J=1,NLOOPLINE
486 S_MAT(I,J)=S_MAT_FROM_MG(I,J)
487 ENDDO
488 ENDDO
489
490 END
491
492 FUNCTION COMPARE_COEFS_0(COEFS_A,COEFS_B)
493
494 USE PRECISION_GOLEM, ONLY: KI
495 COMPLEX(KI) COEFS_A,COEFS_B
496 REAL*8 COEF_CHECK_THRS
497 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
498 REAL*8 DENOM,NUM
499 LOGICAL COMPARE_COEFS_0
500
501 NUM = ABS(COEFS_A-COEFS_B)
502 DENOM = ABS(COEFS_A+COEFS_B)
503 IF(DENOM.GT.0D0)THEN
504 COMPARE_COEFS_0=((NUM/DENOM).LT.COEF_CHECK_THRS)
505 ELSE
506 COMPARE_COEFS_0=(NUM.LT.COEF_CHECK_THRS)
507 ENDIF
508
509 END
510
511 FUNCTION COMPARE_COEFS_1(COEFS_A,COEFS_B)
512
513 USE TENS_REC
514 TYPE(COEFF_TYPE_1)COEFS_A,COEFS_B
515 REAL*8 COEF_CHECK_THRS
516 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
517 REAL*8 DENOM,NUM
518 LOGICAL COMPARE_COEFS_1
519
520 NUM = ABS(COEFS_A%%C0-COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1-COEFS_B%%C1))
521 DENOM = ABS(COEFS_A%%C0+COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1+COEFS_B%%C1)
522 $ )
523
524 IF(DENOM.GT.0D0)THEN
525 COMPARE_COEFS_1=((NUM/DENOM).LT.COEF_CHECK_THRS)
526 ELSE
527 COMPARE_COEFS_1=(NUM.LT.COEF_CHECK_THRS)
528 ENDIF
529
530 END
531
532 FUNCTION COMPARE_COEFS_2(COEFS_A,COEFS_B)
533
534 USE TENS_REC
535 TYPE(COEFF_TYPE_2) COEFS_A,COEFS_B
536 REAL*8 COEF_CHECK_THRS
537 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
538 REAL*8 DENOM,NUM
539 LOGICAL COMPARE_COEFS_2
540
541 NUM = ABS(COEFS_A%%C0-COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1-COEFS_B%%C1))
542 $ +SUM(ABS(COEFS_A%%C2-COEFS_B%%C2))
543 DENOM = ABS(COEFS_A%%C0+COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1+COEFS_B%%C1)
544 $ )+SUM(ABS(COEFS_A%%C2+COEFS_B%%C2))
545 IF(DENOM.GT.0D0)THEN
546 COMPARE_COEFS_2=((NUM/DENOM).LT.COEF_CHECK_THRS)
547 ELSE
548 COMPARE_COEFS_2=(NUM.LT.COEF_CHECK_THRS)
549 ENDIF
550
551 END
552
553 FUNCTION COMPARE_COEFS_3(COEFS_A,COEFS_B)
554
555 USE TENS_REC
556 TYPE(COEFF_TYPE_3) COEFS_A, COEFS_B
557 REAL*8 COEF_CHECK_THRS
558 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
559 REAL*8 DENOM,NUM
560 LOGICAL COMPARE_COEFS_3
561
562 NUM = ABS(COEFS_A%%C0-COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1-COEFS_B%%C1))
563 $ +SUM(ABS(COEFS_A%%C2-COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3-COEFS_B%%C3))
564 DENOM = ABS(COEFS_A%%C0+COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1+COEFS_B%%C1)
565 $ )+SUM(ABS(COEFS_A%%C2+COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3+COEFS_B%%C3)
566 $ )
567 IF(DENOM.GT.0D0)THEN
568 COMPARE_COEFS_3=((NUM/DENOM).LT.COEF_CHECK_THRS)
569 ELSE
570 COMPARE_COEFS_3=(NUM.LT.COEF_CHECK_THRS)
571 ENDIF
572
573 END
574
575 FUNCTION COMPARE_COEFS_4(COEFS_A,COEFS_B)
576
577 USE TENS_REC
578 TYPE(COEFF_TYPE_4) COEFS_A, COEFS_B
579 REAL*8 COEF_CHECK_THRS
580 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
581 REAL*8 DENOM,NUM
582 LOGICAL COMPARE_COEFS_4
583
584 NUM = ABS(COEFS_A%%C0-COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1-COEFS_B%%C1))
585 $ +SUM(ABS(COEFS_A%%C2-COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3-COEFS_B%%C3))
586 $ +SUM(ABS(COEFS_A%%C4-COEFS_B%%C4))
587 DENOM = ABS(COEFS_A%%C0+COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1+COEFS_B%%C1)
588 $ )+SUM(ABS(COEFS_A%%C2+COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3+COEFS_B%%C3)
589 $ )+SUM(ABS(COEFS_A%%C4+COEFS_B%%C4))
590 IF(DENOM.GT.0D0)THEN
591 COMPARE_COEFS_4=((NUM/DENOM).LT.COEF_CHECK_THRS)
592 ELSE
593 COMPARE_COEFS_4=(NUM.LT.COEF_CHECK_THRS)
594 ENDIF
595
596 END
597
598 FUNCTION COMPARE_COEFS_5(COEFS_A,COEFS_B)
599
600 USE TENS_REC
601 TYPE(COEFF_TYPE_5) COEFS_A,COEFS_B
602 REAL*8 COEF_CHECK_THRS
603 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
604 REAL*8 DENOM,NUM
605 LOGICAL COMPARE_COEFS_5
606
607 NUM = ABS(COEFS_A%%C0-COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1-COEFS_B%%C1))
608 $ +SUM(ABS(COEFS_A%%C2-COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3-COEFS_B%%C3))
609 $ +SUM(ABS(COEFS_A%%C4-COEFS_B%%C4))
610 DENOM = ABS(COEFS_A%%C0+COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1+COEFS_B%%C1)
611 $ )+SUM(ABS(COEFS_A%%C2+COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3+COEFS_B%%C3)
612 $ )+SUM(ABS(COEFS_A%%C4+COEFS_B%%C4))
613 IF(DENOM.GT.0D0)THEN
614 COMPARE_COEFS_5=((NUM/DENOM).LT.COEF_CHECK_THRS)
615 ELSE
616 COMPARE_COEFS_5=(NUM.LT.COEF_CHECK_THRS)
617 ENDIF
618
619 END
620
621 FUNCTION COMPARE_COEFS_6(COEFS_A,COEFS_B)
622
623 USE TENS_REC
624 TYPE(COEFF_TYPE_6) COEFS_A,COEFS_B
625 REAL*8 COEF_CHECK_THRS
626 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
627 REAL*8 DENOM,NUM
628 LOGICAL COMPARE_COEFS_6
629
630 NUM = ABS(COEFS_A%%C0-COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1-COEFS_B%%C1))
631 $ +SUM(ABS(COEFS_A%%C2-COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3-COEFS_B%%C3))
632 $ +SUM(ABS(COEFS_A%%C4-COEFS_B%%C4))
633 DENOM = ABS(COEFS_A%%C0+COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1+COEFS_B%%C1)
634 $ )+SUM(ABS(COEFS_A%%C2+COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3+COEFS_B%%C3)
635 $ )+SUM(ABS(COEFS_A%%C4+COEFS_B%%C4))
636 IF(DENOM.GT.0D0)THEN
637 COMPARE_COEFS_6=((NUM/DENOM).LT.COEF_CHECK_THRS)
638 ELSE
639 COMPARE_COEFS_6=(NUM.LT.COEF_CHECK_THRS)
640 ENDIF
641
642 END
643
644
645 SUBROUTINE FILL_GOLEM_COEFFS_0(ML_COEFS,GOLEM_COEFS)
646 USE PRECISION_GOLEM, ONLY: KI
647 INCLUDE 'coef_specs.inc'
648 INCLUDE 'loop_max_coefs.inc'
649 COMPLEX*16 ML_COEFS(0:LOOPMAXCOEFS-1)
650 COMPLEX(KI) GOLEM_COEFS
651 GOLEM_COEFS=ML_COEFS(0)
652 END
653
654 SUBROUTINE FILL_GOLEM_COEFFS_1(ML_COEFS,GOLEM_COEFS)
655 USE TENS_REC, ONLY: COEFF_TYPE_1
656 INCLUDE 'coef_specs.inc'
657 INCLUDE 'loop_max_coefs.inc'
658 COMPLEX*16 ML_COEFS(0:LOOPMAXCOEFS-1)
659 TYPE(COEFF_TYPE_1) GOLEM_COEFS
660C Constant coefficient
661 GOLEM_COEFS%%C0=ML_COEFS(0)
662C Coefficient q(0)
663 GOLEM_COEFS%%C1(1,1)=-ML_COEFS(1)
664C Coefficient q(1)
665 GOLEM_COEFS%%C1(2,1)=-ML_COEFS(2)
666C Coefficient q(2)
667 GOLEM_COEFS%%C1(3,1)=-ML_COEFS(3)
668C Coefficient q(3)
669 GOLEM_COEFS%%C1(4,1)=-ML_COEFS(4)
670 END
671
672 SUBROUTINE FILL_GOLEM_COEFFS_2(ML_COEFS,GOLEM_COEFS)
673 USE TENS_REC, ONLY: COEFF_TYPE_2
674 INCLUDE 'coef_specs.inc'
675 INCLUDE 'loop_max_coefs.inc'
676 COMPLEX*16 ML_COEFS(0:LOOPMAXCOEFS-1)
677 TYPE(COEFF_TYPE_2) GOLEM_COEFS
678C Constant coefficient
679 GOLEM_COEFS%%C0=ML_COEFS(0)
680C Coefficient q(0)
681 GOLEM_COEFS%%C1(1,1)=-ML_COEFS(1)
682C Coefficient q(0)^2
683 GOLEM_COEFS%%C1(1,2)= ML_COEFS(5)
684C Coefficient q(1)
685 GOLEM_COEFS%%C1(2,1)=-ML_COEFS(2)
686C Coefficient q(1)^2
687 GOLEM_COEFS%%C1(2,2)= ML_COEFS(7)
688C Coefficient q(2)
689 GOLEM_COEFS%%C1(3,1)=-ML_COEFS(3)
690C Coefficient q(2)^2
691 GOLEM_COEFS%%C1(3,2)= ML_COEFS(10)
692C Coefficient q(3)
693 GOLEM_COEFS%%C1(4,1)=-ML_COEFS(4)
694C Coefficient q(3)^2
695 GOLEM_COEFS%%C1(4,2)= ML_COEFS(14)
696C Coefficient q(0)*q(1)
697 GOLEM_COEFS%%C2(1,1)= ML_COEFS(6)
698C Coefficient q(0)*q(2)
699 GOLEM_COEFS%%C2(2,1)= ML_COEFS(8)
700C Coefficient q(0)*q(3)
701 GOLEM_COEFS%%C2(3,1)= ML_COEFS(11)
702C Coefficient q(1)*q(2)
703 GOLEM_COEFS%%C2(4,1)= ML_COEFS(9)
704C Coefficient q(1)*q(3)
705 GOLEM_COEFS%%C2(5,1)= ML_COEFS(12)
706C Coefficient q(2)*q(3)
707 GOLEM_COEFS%%C2(6,1)= ML_COEFS(13)
708 END
709
710 SUBROUTINE FILL_GOLEM_COEFFS_3(ML_COEFS,GOLEM_COEFS)
711 USE TENS_REC, ONLY: COEFF_TYPE_3
712 INCLUDE 'coef_specs.inc'
713 INCLUDE 'loop_max_coefs.inc'
714 COMPLEX*16 ML_COEFS(0:LOOPMAXCOEFS-1)
715 TYPE(COEFF_TYPE_3) GOLEM_COEFS
716C Dummy routine for FILL_GOLEM_COEFS_3
717 STOP 'ERROR: 3 > 2'
718 END
719
720 SUBROUTINE FILL_GOLEM_COEFFS_4(ML_COEFS,GOLEM_COEFS)
721 USE TENS_REC, ONLY: COEFF_TYPE_4
722 INCLUDE 'coef_specs.inc'
723 INCLUDE 'loop_max_coefs.inc'
724 COMPLEX*16 ML_COEFS(0:LOOPMAXCOEFS-1)
725 TYPE(COEFF_TYPE_4) GOLEM_COEFS
726C Dummy routine for FILL_GOLEM_COEFS_4
727 STOP 'ERROR: 4 > 2'
728 END
729
730 SUBROUTINE FILL_GOLEM_COEFFS_5(ML_COEFS,GOLEM_COEFS)
731 USE TENS_REC, ONLY: COEFF_TYPE_5
732 INCLUDE 'coef_specs.inc'
733 INCLUDE 'loop_max_coefs.inc'
734 COMPLEX*16 ML_COEFS(0:LOOPMAXCOEFS-1)
735 TYPE(COEFF_TYPE_5) GOLEM_COEFS
736C Dummy routine for FILL_GOLEM_COEFS_5
737 STOP 'ERROR: 5 > 2'
738 END
739
740 SUBROUTINE FILL_GOLEM_COEFFS_6(ML_COEFS,GOLEM_COEFS)
741 USE TENS_REC, ONLY: COEFF_TYPE_6
742 INCLUDE 'coef_specs.inc'
743 INCLUDE 'loop_max_coefs.inc'
744 COMPLEX*16 ML_COEFS(0:LOOPMAXCOEFS-1)
745 TYPE(COEFF_TYPE_6) GOLEM_COEFS
746C Dummy routine for FILL_GOLEM_COEFS_6
747 STOP 'ERROR: 6 > 2'
748 END
7490
=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%TIR_interface.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%TIR_interface.f 2016-08-07 07:16:07 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%TIR_interface.f 2018-11-09 10:32:12 +0000
@@ -96,17 +96,9 @@
96 ENDIF96 ENDIF
9797
98 IF (MLREDUCTIONLIB(I_LIB).EQ.4) THEN98 IF (MLREDUCTIONLIB(I_LIB).EQ.4) THEN
99C Using Golem9599C Golem95 not available
100C PDEN is dummy for Golem95 so we just initialize it to zero100 WRITE(*,*) 'ERROR:: Golem95 is not interfaced.'
101C here so as to use it for the function SWITCHORDER101 STOP
102 DO I=0,3
103 DO J=1,NLOOPLINE-1
104 PDEN(I,J)=0.0D0
105 ENDDO
106 ENDDO
107 CALL SWITCH_ORDER(CTMODE,NLOOPLINE,PL,PDEN,M2L)
108 CALL GOLEMLOOP(NLOOPLINE,PL,M2L,RANK,RES,STABLE)
109 RETURN
110 ENDIF102 ENDIF
111103
112C INITIALIZE TIR IF NEEDED104C INITIALIZE TIR IF NEEDED
113105
=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_matrix.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_matrix.f 2017-08-14 07:04:15 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_matrix.f 2018-11-09 10:32:12 +0000
@@ -244,7 +244,7 @@
244C A FLAG TO DENOTE WHETHER THE CORRESPONDING LOOPLIBS ARE244C A FLAG TO DENOTE WHETHER THE CORRESPONDING LOOPLIBS ARE
245C AVAILABLE OR NOT245C AVAILABLE OR NOT
246 LOGICAL LOOPLIBS_AVAILABLE(NLOOPLIB)246 LOGICAL LOOPLIBS_AVAILABLE(NLOOPLIB)
247 DATA LOOPLIBS_AVAILABLE/.TRUE.,.FALSE.,.TRUE.,.TRUE.,.FALSE.247 DATA LOOPLIBS_AVAILABLE/.TRUE.,.FALSE.,.TRUE.,.FALSE.,.FALSE.
248 $ ,.TRUE.,.TRUE./248 $ ,.TRUE.,.TRUE./
249 COMMON/LOOPLIBS_AV/ LOOPLIBS_AVAILABLE249 COMMON/LOOPLIBS_AV/ LOOPLIBS_AVAILABLE
250C A FLAG TO DENOTE WHETHER THE CORRESPONDING DIRECTION TESTS250C A FLAG TO DENOTE WHETHER THE CORRESPONDING DIRECTION TESTS
251251
=== removed file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%GOLEM_interface.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%GOLEM_interface.f 2016-08-07 07:16:07 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%GOLEM_interface.f 1970-01-01 00:00:00 +0000
@@ -1,748 +0,0 @@
1 SUBROUTINE GOLEMLOOP(NLOOPLINE,PL,M2L,RANK,RES,STABLE)
2C
3C Generated by MadGraph5_aMC@NLO v. %(version)s, %(date)s
4C By the MadGraph5_aMC@NLO Development Team
5C Visit launchpad.net/madgraph5 and amcatnlo.web.cern.ch
6C
7C Interface between MG5 and Golem95.
8C The Golem95 version should be higher than 1.3.0.
9C It supports RANK = NLOOPLINE + 1 tensor integrals when 1 <
10C NLOOPLINE < 6.
11C
12C Process: u d~ > w+ QED<=1 WEIGHTED<=2 [ all = QCD ]
13C Process: c s~ > w+ QED<=1 WEIGHTED<=2 [ all = QCD ]
14C
15C
16C MODULES
17C
18 USE MATRICE_S
19 USE FORM_FACTOR_TYPE, ONLY: FORM_FACTOR
20 USE PRECISION_GOLEM, ONLY: KI
21 USE TENS_COMB
22 USE TENS_REC
23 USE FORM_FACTOR_1P, ONLY: A10
24 USE FORM_FACTOR_2P, ONLY: A20
25 USE FORM_FACTOR_3P, ONLY: A30
26 USE FORM_FACTOR_4P, ONLY: A40
27 USE FORM_FACTOR_5P, ONLY: A50
28 USE FORM_FACTOR_6P, ONLY: A60
29 IMPLICIT NONE
30C
31C CONSTANTS
32C
33 INTEGER NEXTERNAL
34 PARAMETER (NEXTERNAL=3)
35 LOGICAL CHECKPCONSERVATION
36 PARAMETER (CHECKPCONSERVATION=.TRUE.)
37 REAL*8 NORMALIZATION
38 PARAMETER (NORMALIZATION = 1.D0/(16.D0*3.14159265358979323846D0*
39 $ *2))
40 REAL(KI),DIMENSION(0:3),PARAMETER::NULL_VEC = (/0.0_KI,0.0_KI
41 $ ,0.0_KI,0.0_KI/)
42C GOLEM_RUN_MODE = 1: Use directly MadLoop tensorial coefficients
43C GOLEM_RUN_MODE = 2: Reconstruct the tensorial coefficeints
44C directly from
45C numerator using golem internal reconstruction routine
46C GOLEM_RUN_MODE = 3: Cross-checked reconstructed coefficients
47C against
48C MadLoop internal ones.
49 INTEGER GOLEM_RUN_MODE
50 PARAMETER (GOLEM_RUN_MODE=1)
51C The following is the acceptance threshold used for
52C GOLEM_RUN_MODE = 3
53 REAL*8 COEF_CHECK_THRS
54 DATA COEF_CHECK_THRS/1.0D-13/
55 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
56
57 LOGICAL PASS_COEF_CHECK
58C
59C ARGUMENTS
60C
61 INTEGER NLOOPLINE, RANK
62 REAL*8 PL(0:3,NLOOPLINE)
63 REAL*8 PCT(0:3,0:NLOOPLINE-1), ABSPCT(0:3)
64 REAL*8 REF_P
65 REAL(KI) PGOLEM(NLOOPLINE,0:3)
66 COMPLEX*16 M2L(NLOOPLINE)
67 COMPLEX(KI) M2LGOLEM(NLOOPLINE)
68 COMPLEX*16 RES(3)
69 LOGICAL STABLE
70C
71C LOCAL VARIABLES
72C
73 INTEGER I, J, K
74 TYPE(FORM_FACTOR)::RES_GOLEM
75
76 COMPLEX(KI)::COEFFS0,COEFFS0_REC
77 TYPE(COEFF_TYPE_1)::COEFFS1,COEFFS1_REC
78 TYPE(COEFF_TYPE_2)::COEFFS2,COEFFS2_REC
79 TYPE(COEFF_TYPE_3)::COEFFS3,COEFFS3_REC
80 TYPE(COEFF_TYPE_4)::COEFFS4,COEFFS4_REC
81 TYPE(COEFF_TYPE_5)::COEFFS5,COEFFS5_REC
82 TYPE(COEFF_TYPE_6)::COEFFS6,COEFFS6_REC
83
84C The pinch propagator optimization is not used, so for now it is
85C always 0.
86 INTEGER PINCH
87C
88C EXTERNAL FUNCTIONS
89C
90 COMPLEX(KI) GOLEM_LOOPNUM
91 EXTERNAL GOLEM_LOOPNUM
92 LOGICAL COMPARE_COEFS_0
93 LOGICAL COMPARE_COEFS_1
94 LOGICAL COMPARE_COEFS_2
95 LOGICAL COMPARE_COEFS_3
96 LOGICAL COMPARE_COEFS_4
97 LOGICAL COMPARE_COEFS_5
98 LOGICAL COMPARE_COEFS_6
99C
100C GLOBAL VARIABLES
101C
102 INCLUDE 'coupl.inc'
103 INTEGER CTMODE
104 REAL*8 LSCALE
105 COMMON/CT/LSCALE,CTMODE
106
107 INTEGER ID,SQSOINDEX,R
108 COMMON/LOOP/ID,SQSOINDEX,R
109
110 LOGICAL CTINIT, TIRINIT, GOLEMINIT, SAMURAIINIT, NINJAINIT
111 $ ,COLLIERINIT
112 COMMON/REDUCTIONCODEINIT/CTINIT, TIRINIT,GOLEMINIT,SAMURAIINIT
113 $ ,NINJAINIT,COLLIERINIT
114
115 INTEGER NLOOPGROUPS
116 PARAMETER (NLOOPGROUPS=1)
117 INTEGER NSQUAREDSO
118 PARAMETER (NSQUAREDSO=1)
119 INCLUDE 'loop_max_coefs.inc'
120
121 COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
122 COMMON/LCOEFS/LOOPCOEFS
123C ----------
124C BEGIN CODE
125C ----------
126
127C The CT initialization is also performed here if not done already
128C because it calls MPINIT of OneLOop which is necessary on some
129C system
130 IF (CTINIT) THEN
131 CTINIT=.FALSE.
132 CALL INITCT()
133 ENDIF
134
135C INITIALIZE GOLEM IF NEEDED
136 IF (GOLEMINIT) THEN
137 GOLEMINIT=.FALSE.
138 CALL INITGOLEM()
139 ENDIF
140
141C No stability test intrisic to Golem95 now
142 STABLE=.TRUE.
143
144C This initialization must be done for each reduction because we
145C have not setup anyoptimization using pinched propagators yet.
146 CALL INITGOLEM95(NLOOPLINE)
147 PINCH = 0
148
149C YOU CAN FIND THE DETAILS ABOUT THE DIFFERENT CTMODE AT THE
150C BEGINNING OF THE FILE CTS_CUTS.F90 IN THE CUTTOOLS DISTRIBUTION
151
152C CONVERT THE MASSES TO BE COMPLEX
153 DO I=1,NLOOPLINE
154 M2LGOLEM(I)=M2L(I)
155 ENDDO
156
157C CONVERT THE MOMENTA FLOWING IN THE LOOP LINES TO CT CONVENTIONS
158 DO I=0,3
159 ABSPCT(I) = 0.D0
160 DO J=0,(NLOOPLINE-1)
161 PCT(I,J)=0.D0
162 ENDDO
163 ENDDO
164 DO I=0,3
165 DO J=1,NLOOPLINE
166 PCT(I,0)=PCT(I,0)+PL(I,J)
167 ABSPCT(I)=ABSPCT(I)+ABS(PL(I,J))
168 ENDDO
169 ENDDO
170 REF_P = MAX(ABSPCT(0), ABSPCT(1),ABSPCT(2),ABSPCT(3))
171 DO I=0,3
172 ABSPCT(I) = MAX(REF_P*1E-6, ABSPCT(I))
173 ENDDO
174 IF (CHECKPCONSERVATION.AND.REF_P.GT.1D-8) THEN
175 IF ((PCT(0,0)/ABSPCT(0)).GT.1.D-6) THEN
176 WRITE(*,*) 'energy is not conserved ',PCT(0,0)
177 STOP 'energy is not conserved'
178 ELSEIF ((PCT(1,0)/ABSPCT(1)).GT.1.D-6) THEN
179 WRITE(*,*) 'px is not conserved ',PCT(1,0)
180 STOP 'px is not conserved'
181 ELSEIF ((PCT(2,0)/ABSPCT(2)).GT.1.D-6) THEN
182 WRITE(*,*) 'py is not conserved ',PCT(2,0)
183 STOP 'py is not conserved'
184 ELSEIF ((PCT(3,0)/ABSPCT(3)).GT.1.D-6) THEN
185 WRITE(*,*) 'pz is not conserved ',PCT(3,0)
186 STOP 'pz is not conserved'
187 ENDIF
188 ENDIF
189 DO I=0,3
190 DO J=1,(NLOOPLINE-1)
191 DO K=1,J
192 PCT(I,J)=PCT(I,J)+PL(I,K)
193 ENDDO
194 ENDDO
195 ENDDO
196
197C Now convert the loop momenta to Golem95 conventions
198 DO I=0,3
199 PGOLEM(1,I)=0.0E0_KI
200 DO J=2,NLOOPLINE
201 PGOLEM(J,I)=PCT(I,J-1)
202 ENDDO
203 ENDDO
204
205C Fill in the kinematic s-matrix while taking care of on-shell
206C limits.
207 CALL SETUP_KIN_MATRIX(NLOOPLINE,PGOLEM,M2LGOLEM)
208C Construct the golem internal matrices derived from the kinetic
209C one.
210 CALL PREPARESMATRIX()
211
212C Fill in the golem coefficents and compute the loop
213 IF(GOLEM_RUN_MODE.EQ.2)THEN
214 RES_GOLEM = EVALUATE_B(GOLEM_LOOPNUM,PGOLEM,0,RANK)
215 ELSE
216 PASS_COEF_CHECK=.TRUE.
217 SELECT CASE(RANK)
218 CASE(0)
219 CALL FILL_GOLEM_COEFFS_0(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS0)
220 IF(GOLEM_RUN_MODE.EQ.3)THEN
221 COEFFS0_REC = GOLEM_LOOPNUM(NULL_VEC,0.0_KI)
222 PASS_COEF_CHECK=COMPARE_COEFS_0(COEFFS0,COEFFS0_REC)
223 ENDIF
224 CASE(1)
225 CALL FILL_GOLEM_COEFFS_1(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS1)
226 IF(GOLEM_RUN_MODE.EQ.3)THEN
227 CALL RECONSTRUCT1(GOLEM_LOOPNUM,COEFFS1_REC)
228 PASS_COEF_CHECK=COMPARE_COEFS_1(COEFFS1,COEFFS1_REC)
229 ENDIF
230 CASE(2)
231 CALL FILL_GOLEM_COEFFS_2(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS2)
232 IF(GOLEM_RUN_MODE.EQ.3)THEN
233 CALL RECONSTRUCT2(GOLEM_LOOPNUM,COEFFS2_REC)
234 PASS_COEF_CHECK=COMPARE_COEFS_2(COEFFS2,COEFFS2_REC)
235 ENDIF
236 CASE(3)
237 CALL FILL_GOLEM_COEFFS_3(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS3)
238 IF(GOLEM_RUN_MODE.EQ.3)THEN
239 CALL RECONSTRUCT3(GOLEM_LOOPNUM,COEFFS3_REC)
240 PASS_COEF_CHECK=COMPARE_COEFS_3(COEFFS3,COEFFS3_REC)
241 ENDIF
242 CASE(4)
243 CALL FILL_GOLEM_COEFFS_4(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS4)
244 IF(GOLEM_RUN_MODE.EQ.3)THEN
245 CALL RECONSTRUCT4(GOLEM_LOOPNUM,COEFFS4_REC)
246 PASS_COEF_CHECK=COMPARE_COEFS_4(COEFFS4,COEFFS4_REC)
247 ENDIF
248 CASE(5)
249 CALL FILL_GOLEM_COEFFS_5(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS5)
250 IF(GOLEM_RUN_MODE.EQ.3)THEN
251 CALL RECONSTRUCT5(GOLEM_LOOPNUM,COEFFS5_REC)
252 PASS_COEF_CHECK=COMPARE_COEFS_5(COEFFS5,COEFFS5_REC)
253 ENDIF
254 CASE(6)
255 CALL FILL_GOLEM_COEFFS_6(LOOPCOEFS(0,SQSOINDEX,ID),COEFFS6)
256 IF(GOLEM_RUN_MODE.EQ.3)THEN
257 CALL RECONSTRUCT6(GOLEM_LOOPNUM,COEFFS6_REC)
258 PASS_COEF_CHECK=COMPARE_COEFS_6(COEFFS6,COEFFS6_REC)
259 ENDIF
260 CASE DEFAULT
261 WRITE(*,*)'Not yet implemented in Golem95 for rank= ',RANK
262 STOP
263 END SELECT
264
265 IF(.NOT.PASS_COEF_CHECK)THEN
266 WRITE(*,*)'Coefs mismatch for ID ',ID,' and rank ',RANK
267 WRITE(*,*)'Coefs form MadLoop5:'
268 SELECT CASE(RANK)
269 CASE(0)
270 WRITE(*,*)'Constant coef = ',COEFFS0
271 CASE(1)
272 CALL PRINT_COEFFS(COEFFS1)
273 CASE(2)
274 CALL PRINT_COEFFS(COEFFS2)
275 CASE(3)
276 CALL PRINT_COEFFS(COEFFS3)
277 CASE(4)
278 CALL PRINT_COEFFS(COEFFS4)
279 CASE(5)
280 CALL PRINT_COEFFS(COEFFS5)
281 CASE(6)
282 CALL PRINT_COEFFS(COEFFS6)
283 END SELECT
284 WRITE(*,*)'Coefs reconstructed by Golem95:'
285 SELECT CASE(RANK)
286 CASE(0)
287 WRITE(*,*)'Constant coef = ',COEFFS0_REC
288 CASE(1)
289 CALL PRINT_COEFFS(COEFFS1_REC)
290 CASE(2)
291 CALL PRINT_COEFFS(COEFFS2_REC)
292 CASE(3)
293 CALL PRINT_COEFFS(COEFFS3_REC)
294 CASE(4)
295 CALL PRINT_COEFFS(COEFFS4_REC)
296 CASE(5)
297 CALL PRINT_COEFFS(COEFFS5_REC)
298 CASE(6)
299 CALL PRINT_COEFFS(COEFFS6_REC)
300 END SELECT
301 STOP
302 ENDIF
303
304 SELECT CASE(NLOOPLINE)
305 CASE(1)
306 WRITE(*,*)'Golem95 cannot handle with tadpole yet'
307 STOP
308 CASE(2)
309 SELECT CASE(RANK)
310 CASE(0)
311 RES_GOLEM = COEFFS0*A20(PINCH)
312 CASE(1)
313 RES_GOLEM = CONTRACT2_1(COEFFS1,PGOLEM,PINCH)
314 CASE(2)
315 RES_GOLEM = CONTRACT2_2(COEFFS2,PGOLEM,PINCH)
316 CASE(3)
317 RES_GOLEM = CONTRACT2_3(COEFFS3,PGOLEM,PINCH)
318 CASE DEFAULT
319 WRITE(*,*)'Golem95 cannot handle with: N,r = ',2,RANK
320 STOP
321 END SELECT
322 CASE(3)
323 SELECT CASE(RANK)
324 CASE(0)
325 RES_GOLEM = COEFFS0*A30(PINCH)
326 CASE(1)
327 RES_GOLEM = CONTRACT3_1(COEFFS1,PGOLEM,PINCH)
328 CASE(2)
329 RES_GOLEM = CONTRACT3_2(COEFFS2,PGOLEM,PINCH)
330 CASE(3)
331 RES_GOLEM = CONTRACT3_3(COEFFS3,PGOLEM,PINCH)
332 CASE(4)
333 RES_GOLEM = CONTRACT3_4(COEFFS4,PGOLEM,PINCH)
334 CASE DEFAULT
335 WRITE(*,*)'Golem95 cannot handle with: N,r = ',3,RANK
336 STOP
337 END SELECT
338 CASE(4)
339 SELECT CASE(RANK)
340 CASE(0)
341 RES_GOLEM = COEFFS0*A40(PINCH)
342 CASE(1)
343 RES_GOLEM = CONTRACT4_1(COEFFS1,PGOLEM,PINCH)
344 CASE(2)
345 RES_GOLEM = CONTRACT4_2(COEFFS2,PGOLEM,PINCH)
346 CASE(3)
347 RES_GOLEM = CONTRACT4_3(COEFFS3,PGOLEM,PINCH)
348 CASE(4)
349 RES_GOLEM = CONTRACT4_4(COEFFS4,PGOLEM,PINCH)
350 CASE(5)
351 RES_GOLEM = CONTRACT4_5(COEFFS5,PGOLEM,PINCH)
352 CASE DEFAULT
353 WRITE(*,*)'Golem95 cannot handle with: N,r = ',4,RANK
354 STOP
355 END SELECT
356 CASE(5)
357 SELECT CASE(RANK)
358 CASE(0)
359 RES_GOLEM = COEFFS0*A50(PINCH)
360 CASE(1)
361 RES_GOLEM = CONTRACT5_1(COEFFS1,PGOLEM,PINCH)
362 CASE(2)
363 RES_GOLEM = CONTRACT5_2(COEFFS2,PGOLEM,PINCH)
364 CASE(3)
365 RES_GOLEM = CONTRACT5_3(COEFFS3,PGOLEM,PINCH)
366 CASE(4)
367 RES_GOLEM = CONTRACT5_4(COEFFS4,PGOLEM,PINCH)
368 CASE(5)
369 RES_GOLEM = CONTRACT5_5(COEFFS5,PGOLEM,PINCH)
370 CASE(6)
371 RES_GOLEM = CONTRACT5_6(COEFFS6,PGOLEM,PINCH)
372 CASE DEFAULT
373 WRITE(*,*)'Golem95 cannot handle with: N,r = ',5,RANK
374 STOP
375 END SELECT
376 CASE(6)
377 SELECT CASE(RANK)
378 CASE(0)
379 RES_GOLEM = COEFFS0*A60(PINCH)
380 CASE(1)
381 RES_GOLEM = CONTRACT6_1(COEFFS1,PGOLEM,PINCH)
382 CASE(2)
383 RES_GOLEM = CONTRACT6_2(COEFFS2,PGOLEM,PINCH)
384 CASE(3)
385 RES_GOLEM = CONTRACT6_3(COEFFS3,PGOLEM,PINCH)
386 CASE(4)
387 RES_GOLEM = CONTRACT6_4(COEFFS4,PGOLEM,PINCH)
388 CASE(5)
389 RES_GOLEM = CONTRACT6_5(COEFFS5,PGOLEM,PINCH)
390 CASE(6)
391 RES_GOLEM = CONTRACT6_6(COEFFS6,PGOLEM,PINCH)
392 CASE DEFAULT
393 WRITE(*,*)'Golem95 cannot handle with: N,r = ',6,RANK
394 STOP
395 END SELECT
396 CASE DEFAULT
397 WRITE(*,*)'Golem95 cannot handle with: N = ',NLOOPLINE
398 STOP
399 END SELECT
400 ENDIF
401
402 RES(1)=NORMALIZATION*2.0D0*DBLE(RES_GOLEM%%C+2.0*LOG(MU_R)
403 $ *RES_GOLEM%%B+2.0*LOG(MU_R)**2*RES_GOLEM%%A)
404 RES(2)=NORMALIZATION*2.0D0*DBLE(RES_GOLEM%%B+2.0*LOG(MU_R)
405 $ *RES_GOLEM%%A)
406 RES(3)=NORMALIZATION*2.0D0*DBLE(RES_GOLEM%%A)
407C WRITE(*,*) 'Loop ID',ID,' =',RES(1),RES(2),RES(3)
408
409C Finally free golem memory and cache
410 CALL EXITGOLEM95()
411
412 END
413
414 FUNCTION GOLEM_LOOPNUM(Q,MU2)
415 USE PRECISION_GOLEM, ONLY: KI
416 REAL(KI),DIMENSION(0:3),INTENT(IN)::Q
417 REAL(KI),INTENT(IN)::MU2
418 COMPLEX(KI)::GOLEM_LOOPNUM
419
420 COMPLEX*16 QQ(0:3),NUM
421 INTEGER I
422
423 DO I=0,3
424 QQ(I)=CMPLX(Q(I),0.0D0,KIND=16)
425 ENDDO
426
427 CALL LOOPNUM(QQ,NUM)
428 GOLEM_LOOPNUM=NUM
429 RETURN
430 END FUNCTION
431
432 SUBROUTINE INITGOLEM()
433C
434C INITIALISATION OF GOLEM
435C
436C
437C MODULE
438C
439 USE PARAMETRE
440C
441C LOCAL VARIABLES
442C
443 REAL*8 THRS
444 LOGICAL EXT_NUM_FOR_R1
445C
446C GLOBAL VARIABLES
447C
448 INCLUDE 'MadLoopParams.inc'
449C ----------
450C BEGIN CODE
451C ----------
452
453C DEFAULT PARAMETERS FOR GOLEM
454C -------------------------------
455C One can chose here to have either just the rational R1 piece
456C or everything but the R2
457 RAT_OR_TOT_PAR = TOT
458
459 END
460
461 SUBROUTINE SETUP_KIN_MATRIX(NLOOPLINE,PGOLEM,M2L)
462C
463C MODULE
464C
465 USE MATRICE_S
466 USE PRECISION_GOLEM, ONLY: KI
467C
468C ARGUMENTS
469C
470 INTEGER NLOOPLINE
471 REAL(KI) PGOLEM(NLOOPLINE,0:3)
472 COMPLEX(KI) M2L(NLOOPLINE)
473C
474C LOCAL VARIABLES
475C
476 INTEGER I,J
477 COMPLEX*16 S_MAT_FROM_MG(NLOOPLINE,NLOOPLINE)
478C ----------
479C BEGIN CODE
480C ----------
481
482 CALL BUILD_KINEMATIC_MATRIX(NLOOPLINE,PGOLEM,M2L,S_MAT_FROM_MG)
483
484 DO I=1,NLOOPLINE
485 DO J=1,NLOOPLINE
486 S_MAT(I,J)=S_MAT_FROM_MG(I,J)
487 ENDDO
488 ENDDO
489
490 END
491
492 FUNCTION COMPARE_COEFS_0(COEFS_A,COEFS_B)
493
494 USE PRECISION_GOLEM, ONLY: KI
495 COMPLEX(KI) COEFS_A,COEFS_B
496 REAL*8 COEF_CHECK_THRS
497 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
498 REAL*8 DENOM,NUM
499 LOGICAL COMPARE_COEFS_0
500
501 NUM = ABS(COEFS_A-COEFS_B)
502 DENOM = ABS(COEFS_A+COEFS_B)
503 IF(DENOM.GT.0D0)THEN
504 COMPARE_COEFS_0=((NUM/DENOM).LT.COEF_CHECK_THRS)
505 ELSE
506 COMPARE_COEFS_0=(NUM.LT.COEF_CHECK_THRS)
507 ENDIF
508
509 END
510
511 FUNCTION COMPARE_COEFS_1(COEFS_A,COEFS_B)
512
513 USE TENS_REC
514 TYPE(COEFF_TYPE_1)COEFS_A,COEFS_B
515 REAL*8 COEF_CHECK_THRS
516 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
517 REAL*8 DENOM,NUM
518 LOGICAL COMPARE_COEFS_1
519
520 NUM = ABS(COEFS_A%%C0-COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1-COEFS_B%%C1))
521 DENOM = ABS(COEFS_A%%C0+COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1+COEFS_B%%C1)
522 $ )
523
524 IF(DENOM.GT.0D0)THEN
525 COMPARE_COEFS_1=((NUM/DENOM).LT.COEF_CHECK_THRS)
526 ELSE
527 COMPARE_COEFS_1=(NUM.LT.COEF_CHECK_THRS)
528 ENDIF
529
530 END
531
532 FUNCTION COMPARE_COEFS_2(COEFS_A,COEFS_B)
533
534 USE TENS_REC
535 TYPE(COEFF_TYPE_2) COEFS_A,COEFS_B
536 REAL*8 COEF_CHECK_THRS
537 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
538 REAL*8 DENOM,NUM
539 LOGICAL COMPARE_COEFS_2
540
541 NUM = ABS(COEFS_A%%C0-COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1-COEFS_B%%C1))
542 $ +SUM(ABS(COEFS_A%%C2-COEFS_B%%C2))
543 DENOM = ABS(COEFS_A%%C0+COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1+COEFS_B%%C1)
544 $ )+SUM(ABS(COEFS_A%%C2+COEFS_B%%C2))
545 IF(DENOM.GT.0D0)THEN
546 COMPARE_COEFS_2=((NUM/DENOM).LT.COEF_CHECK_THRS)
547 ELSE
548 COMPARE_COEFS_2=(NUM.LT.COEF_CHECK_THRS)
549 ENDIF
550
551 END
552
553 FUNCTION COMPARE_COEFS_3(COEFS_A,COEFS_B)
554
555 USE TENS_REC
556 TYPE(COEFF_TYPE_3) COEFS_A, COEFS_B
557 REAL*8 COEF_CHECK_THRS
558 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
559 REAL*8 DENOM,NUM
560 LOGICAL COMPARE_COEFS_3
561
562 NUM = ABS(COEFS_A%%C0-COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1-COEFS_B%%C1))
563 $ +SUM(ABS(COEFS_A%%C2-COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3-COEFS_B%%C3))
564 DENOM = ABS(COEFS_A%%C0+COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1+COEFS_B%%C1)
565 $ )+SUM(ABS(COEFS_A%%C2+COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3+COEFS_B%%C3)
566 $ )
567 IF(DENOM.GT.0D0)THEN
568 COMPARE_COEFS_3=((NUM/DENOM).LT.COEF_CHECK_THRS)
569 ELSE
570 COMPARE_COEFS_3=(NUM.LT.COEF_CHECK_THRS)
571 ENDIF
572
573 END
574
575 FUNCTION COMPARE_COEFS_4(COEFS_A,COEFS_B)
576
577 USE TENS_REC
578 TYPE(COEFF_TYPE_4) COEFS_A, COEFS_B
579 REAL*8 COEF_CHECK_THRS
580 COMMON/COEF_CHECK_THRS/COEF_CHECK_THRS
581 REAL*8 DENOM,NUM
582 LOGICAL COMPARE_COEFS_4
583
584 NUM = ABS(COEFS_A%%C0-COEFS_B%%C0)+SUM(ABS(COEFS_A%%C1-COEFS_B%%C1))
585 $ +SUM(ABS(COEFS_A%%C2-COEFS_B%%C2))+SUM(ABS(COEFS_A%%C3-COEFS_B%%C3))