2
c*****************************************************************************
3
c Given identical particles, and the configurations. This program identifies
4
c identical configurations and specifies which ones can be skipped
5
c*****************************************************************************
11
include 'nexternal.inc'
12
include '../../Source/run_config.inc'
13
include 'nFKSconfigs.inc'
14
include 'fks_info.inc'
19
double precision ZERO,one
20
parameter (ZERO = 0d0)
23
parameter(maxswitch=99)
27
integer itree(2,-max_branch:-1)
29
integer i,j, k, n, nsym,l,ii,jj
30
double precision diff,xi_i_fks
31
double precision pmass(nexternal)
33
integer biforest(2,-max_branch:-1,lmaxconfigs)
34
integer fksmother,fksgrandmother,fksaunt,compare
35
integer fksconfiguration,mapbconf(0:lmaxconfigs)
36
integer r2b(lmaxconfigs),b2r(lmaxconfigs)
37
logical searchforgranny,is_beta_cms,is_granny_sch,topdown,non_prop
38
integer nbranch,ns_channel,nt_channel
40
integer fks_j_from_i(nexternal,0:nexternal)
41
& ,particle_type(nexternal),pdg_type(nexternal)
42
common /c_fks_inc/fks_j_from_i,particle_type,pdg_type
43
double precision fxl(15),wfxl(15),limit(15),wlimit(15)
44
double precision lxp(15,0:3,nexternal+1),xp(15,0:3,nexternal+1)
45
double precision fks_Sij
46
double precision check,tolerance,zh,h_damp
47
parameter (tolerance=1.d-4)
48
integer kk,ll,bs,bs_min,bs_max,iconfig_in
50
integer nsofttests,ncolltests,nerr,imax,iflag,iret
52
c Local for generating amps
54
double precision p(0:3,99), wgt, x(99), fx
55
double complex wgt1(2)
56
double precision p1(0:3,99),xx(maxinvar)
57
integer ninvar, ndim, minconfig, maxconfig
59
integer ncall,itmax,nconfigs,ntry, ngraphs
60
integer ic(nexternal,maxswitch), jc(12),nswitch
61
double precision saveamp(maxamps)
65
double precision totmass
67
double precision xi_i_fks_fix_save,y_ij_fks_fix_save
68
double precision xi_i_fks_fix,y_ij_fks_fix
69
common/cxiyfix/xi_i_fks_fix,y_ij_fks_fix
73
Double Precision amp2(maxamps), jamp2(0:maxamps)
74
common/to_amps/ amp2, jamp2
77
logical calculatedBorn
78
common/ccalculatedBorn/calculatedBorn
81
common/fks_indices/i_fks,j_fks
83
double precision p1_cnt(0:3,nexternal,-2:2)
84
double precision wgt_cnt(-2:2)
85
double precision pswgt_cnt(-2:2)
86
double precision jac_cnt(-2:2)
87
common/counterevnts/p1_cnt,wgt_cnt,pswgt_cnt,jac_cnt
89
double precision p_born(0:3,nexternal-1)
92
double precision xi_i_fks_ev,y_ij_fks_ev
93
double precision p_i_fks_ev(0:3),p_i_fks_cnt(0:3,-2:2)
94
common/fksvariables/xi_i_fks_ev,y_ij_fks_ev,p_i_fks_ev,p_i_fks_cnt
96
double precision xi_i_fks_cnt(-2:2)
97
common /cxiifkscnt/xi_i_fks_cnt
100
common/crotategranny/rotategranny
102
logical softtest,colltest
103
common/sctests/softtest,colltest
106
common /toxexternal/ xexternal
108
c Particle types (=color) of i_fks, j_fks and fks_mother
109
integer i_type,j_type,m_type
110
double precision ch_i,ch_j,ch_m
111
common/cparticle_types/i_type,j_type,m_type,ch_i,ch_j,ch_m
112
double precision particle_charge(nexternal), particle_charge_born(nexternal-1)
113
common /c_charges/particle_charge
114
common /c_charges_born/particle_charge_born
120
double precision dsig,ran2
122
external check_swap, fks_Sij
124
c define here the maximum fraction of failures to consider the test
126
double precision max_fail, fail_frac
127
parameter (max_fail=0.3d0)
131
logical multi_channel
132
common/to_matrix/isum_hel, multi_channel
134
integer fks_conf_number,fks_loop_min,fks_loop_max,fks_loop
136
COMMON/C_NFKSPROCESS/NFKSPROCESS
141
integer orders(nsplitorders)
142
double precision fxl_split(15,amp_split_size),wfxl_split(15,amp_split_size)
143
double precision limit_split(15,amp_split_size), wlimit_split(15,amp_split_size)
145
c born configuration stuff
146
include 'born_ngraphs.inc'
147
include 'born_conf.inc'
148
LOGICAL IS_A_J(NEXTERNAL),IS_A_LP(NEXTERNAL),IS_A_LM(NEXTERNAL)
149
LOGICAL IS_A_PH(NEXTERNAL)
150
COMMON /TO_SPECISA/IS_A_J,IS_A_LP,IS_A_LM,IS_A_PH
153
common /c_new_point/new_point
158
if (fks_configs.eq.1) then
159
if (pdg_type_d(1,fks_i_d(1)).eq.-21) then
160
write (*,*) 'Process generated with [LOonly=QCD]. '/
166
write(*,*)'Enter xi_i, y_ij to be used in coll/soft tests'
167
write(*,*)' Enter -2 to generate them randomly'
168
read(*,*)xi_i_fks_fix_save,y_ij_fks_fix_save
170
write(*,*)'Enter number of tests for soft and collinear limits'
171
read(*,*)nsofttests,ncolltests
173
write(*,*)'Sum over helicity (0), or random helicity (1)'
176
call setrun !Sets up run parameters
177
call setpara('param_card.dat') !Sets up couplings and masses
178
call setcuts !Sets up cuts
180
c When doing hadron-hadron collision reduce the effect collision energy.
181
c Note that tests are always performed at fixed energy with Bjorken x=1.
183
include 'pmass.inc' ! make sure to set the masses after the model has been included
184
do i=nincoming+1,nexternal
185
if (is_a_j(i) .and. i.ne.nexternal) then
186
totmass = totmass + max(ptj,pmass(i))
187
elseif ((is_a_lp(i).or.is_a_lm(i)) .and. i.ne.nexternal) then
188
totmass = totmass + max(mll/2d0,mll_sf/2d0,ptl,pmass(i))
190
totmass = totmass + pmass(i)
193
if (lpp(1).ne.0) ebeam(1)=max(ebeam(1)/20d0,totmass)
194
if (lpp(2).ne.0) ebeam(2)=max(ebeam(2)/20d0,totmass)
197
write (*,*) 'Give FKS configuration number ("0" loops over all)'
198
read (*,*) fks_conf_number
200
if (fks_conf_number.eq.0) then
202
fks_loop_max=fks_configs
204
fks_loop_min=fks_conf_number
205
fks_loop_max=fks_conf_number
208
do fks_loop=fks_loop_min,fks_loop_max
211
write (*,*) '================================================='
213
write (*,*) 'NEW FKS CONFIGURATION:'
215
call fks_inc_chooser()
216
call leshouche_inc_chooser()
217
write (*,*) 'FKS configuration number is ',nFKSprocess
218
write (*,*) 'FKS partons are: i=',i_fks,' j=',j_fks
219
write (*,*) 'with PDGs: i=',PDG_type(i_fks),' j='
229
c Set color types of i_fks, j_fks and fks_mother.
230
i_type=particle_type(i_fks)
231
j_type=particle_type(j_fks)
232
ch_i=particle_charge(i_fks)
233
ch_j=particle_charge(j_fks)
234
call get_mother_col_charge(i_type,ch_i,j_type,ch_j,m_type,ch_m)
238
c Get momentum configuration
241
c Set xexternal to true to use the x's from external vegas in the
242
c x_to_f_arg subroutine
247
write(*,*)'Enter graph number (iconfig), '
248
& //"'0' loops over all graphs"
251
if (iconfig_in.eq.0) then
254
elseif (iconfig_in.eq.-1) then
262
do iconfig=bs_min,bs_max ! Born configurations
266
call setfksfactor(.false.)
277
call generate_momenta(ndim,iconfig,wgt,x,p)
278
calculatedBorn=.false.
279
do while (( wgt.lt.0 .or. p(0,1).le.0d0 .or. p_born(0,1).le.0d0
280
& ) .and. ntry .lt. 1000)
286
call generate_momenta(ndim,iconfig,wgt,x,p)
287
calculatedBorn=.false.
291
if (ntry.ge.1000) then
292
write (*,*) 'No points passed cuts...'
293
write (12,*) 'ERROR: no points passed cuts...'
294
& //' Cannot perform ME tests properly for config',iconfig
298
call sborn(p_born,wgt1)
309
call get_helicity(i_fks,j_fks)
310
do iamp=1,amp_split_size
312
fxl_split(i,iamp) = 0d0
313
wfxl_split(i,iamp) = 0d0
314
limit_split(i,iamp) = 0d0
315
wlimit_split(i,iamp) = 0d0
319
if(nsofttests.le.10)then
324
y_ij_fks_fix=y_ij_fks_fix_save
332
call generate_momenta(ndim,iconfig,wgt,x,p)
333
do while (( wgt.lt.0 .or. p(0,1).le.0d0) .and. ntry.lt.1000)
339
call generate_momenta(ndim,iconfig,wgt,x,p)
342
if(nsofttests.le.10)write (*,*) 'ntry',ntry
345
call generate_momenta(ndim,iconfig,wgt,x,p)
346
calculatedBorn=.false.
347
call set_cms_stuff(0)
348
call sreal(p1_cnt(0,1,0),zero,y_ij_fks_ev,fx)
351
! keep track of the separate splitorders
352
do iamp=1,amp_split_size
353
fxl_split(i,iamp) = amp_split(iamp)*jac_cnt(0)
354
wfxl_split(i,iamp)=jac_cnt(0)
356
calculatedBorn=.false.
357
call set_cms_stuff(-100)
358
call sreal(p,xi_i_fks_ev,y_ij_fks_ev,fx)
361
! keep track of the separate splitorders
362
do iamp=1,amp_split_size
363
limit_split(i,iamp) = amp_split(iamp)*wgt
364
wlimit_split(i,iamp) = wgt
369
lxp(i,l,k)=p1_cnt(l,k,0)
373
xp(i,l,nexternal+1)=p_i_fks_ev(l)
374
lxp(i,l,nexternal+1)=p_i_fks_cnt(l,0)
376
xi_i_fks_fix=xi_i_fks_fix/10d0
379
if(nsofttests.le.10)then
380
write (*,*) 'Soft limit:'
381
write (*,*) ' Sum of all contributions:'
383
call xprintout(6,limit(i),fxl(i))
385
! check the contributions coming from each splitorders
386
! only look at the non vanishing ones
387
do iamp=1, amp_split_size
388
if (limit_split(1,iamp).ne.0d0.or.fxl_split(1,iamp).ne.0d0) then
389
write(*,*) ' Split-order', iamp
390
call amp_split_pos_to_orders(iamp,orders)
391
do i = 1, nsplitorders
392
write(*,*) ' ',ordernames(i), ':', orders(i)
395
call xprintout(6,limit_split(i,iamp),fxl_split(i,iamp))
397
call checkres2(limit_split(1,iamp),fxl_split(1,iamp),
398
& wlimit_split(1,iamp),wfxl_split(1,iamp),xp,lxp,
399
& iflag,imax,j,nexternal,i_fks,j_fks,iret)
400
write(*,*) 'RETURN CODE', iret
406
write(80,*)'****************************'
414
call xprintout(80,xp(i,l,k),lxp(i,l,k))
420
call checkres2(limit,fxl,wlimit,wfxl,xp,lxp,
421
& iflag,imax,j,nexternal,i_fks,j_fks,iret)
423
! check the contributions coming from each splitorders
424
! only look at the non vanishing ones
425
do iamp=1, amp_split_size
426
if (limit_split(1,iamp).ne.0d0.or.fxl_split(1,iamp).ne.0d0) then
427
call checkres2(limit_split(1,iamp),fxl_split(1,iamp),
428
& wlimit_split(1,iamp),wfxl_split(1,iamp),xp,lxp,
429
& iflag,imax,j,nexternal,i_fks,j_fks,iret)
436
if(nsofttests.gt.10)then
437
write(*,*)'Soft tests done for (Born) config',iconfig
438
write(*,*)'Failures:',nerr
439
fail_frac= nerr/dble(nsofttests)
440
if (fail_frac.lt.max_fail) then
441
write(*,401) nFKSprocess, fail_frac
443
write(*,402) nFKSprocess, fail_frac
453
if (pmass(j_fks).ne.0d0) then
454
write (*,*) 'No collinear test for massive j_fks'
461
c Set rotategranny=.true. to align grandmother along the z axis, when
462
c grandmother is not the c.m. system (if granny=cms, this rotation coincides
463
c with the identity, and the following is harmless).
464
c WARNING: the setting of rotategranny changes the definition of xij_aor
465
c in genps_fks_test.f
471
call get_helicity(i_fks,j_fks)
472
do iamp=1,amp_split_size
474
fxl_split(i,iamp) = 0d0
475
wfxl_split(i,iamp) = 0d0
476
limit_split(i,iamp) = 0d0
477
wlimit_split(i,iamp) = 0d0
481
if(ncolltests.le.10)then
487
xi_i_fks_fix=xi_i_fks_fix_save
494
call generate_momenta(ndim,iconfig,wgt,x,p)
495
do while (( wgt.lt.0 .or. p(0,1).le.0d0) .and. ntry.lt.1000)
501
call generate_momenta(ndim,iconfig,wgt,x,p)
504
if(ncolltests.le.10)write (*,*) 'ntry',ntry
506
y_ij_fks_fix=1-0.1d0**i
508
call generate_momenta(ndim,iconfig,wgt,x,p)
509
calculatedBorn=.false.
510
call set_cms_stuff(1)
511
call sreal(p1_cnt(0,1,1),xi_i_fks_cnt(1),one,fx)
513
! keep track of the separate splitorders
514
do iamp=1,amp_split_size
515
fxl_split(i,iamp) = amp_split(iamp)*jac_cnt(1)
516
wfxl_split(i,iamp) = jac_cnt(1)
519
calculatedBorn=.false.
520
call set_cms_stuff(-100)
521
call sreal(p,xi_i_fks_ev,y_ij_fks_ev,fx)
524
! keep track of the separate splitorders
525
do iamp=1,amp_split_size
526
limit_split(i,iamp) = amp_split(iamp)*wgt
527
wlimit_split(i,iamp) = wgt
531
lxp(i,l,k)=p1_cnt(l,k,1)
536
lxp(i,l,nexternal+1)=p_i_fks_cnt(l,1)
537
xp(i,l,nexternal+1)=p_i_fks_ev(l)
540
if(ncolltests.le.10)then
541
write (*,*) 'Collinear limit:'
542
write (*,*) ' Sum of all contributions:'
544
call xprintout(6,limit(i),fxl(i))
546
! check the contributions coming from each splitorders
547
! only look at the non vanishing ones
548
do iamp=1, amp_split_size
549
if (limit_split(1,iamp).ne.0d0.or.fxl_split(1,iamp).ne.0d0) then
550
write(*,*) ' Split-order', iamp
551
call amp_split_pos_to_orders(iamp,orders)
552
do i = 1, nsplitorders
553
write(*,*) ' ',ordernames(i), ':', orders(i)
556
call xprintout(6,limit_split(i,iamp),fxl_split(i,iamp))
558
call checkres2(limit_split(1,iamp),fxl_split(1,iamp),
559
& wlimit_split(1,iamp),wfxl_split(1,iamp),xp,lxp,
560
& iflag,imax,j,nexternal,i_fks,j_fks,iret)
561
write(*,*) 'RETURN CODE', iret
566
write(80,*)'****************************'
574
call xprintout(80,xp(i,l,k),lxp(i,l,k))
580
call checkres2(limit,fxl,wlimit,wfxl,xp,lxp,
581
& iflag,imax,j,nexternal,i_fks,j_fks,iret)
583
! check the contributions coming from each splitorders
584
! only look at the non vanishing ones
585
do iamp=1, amp_split_size
586
if (limit_split(1,iamp).ne.0d0.or.fxl_split(1,iamp).ne.0d0) then
587
call checkres2(limit_split(1,iamp),fxl_split(1,iamp),
588
& wlimit_split(1,iamp),wfxl_split(1,iamp),xp,lxp,
589
& iflag,imax,j,nexternal,i_fks,j_fks,iret)
595
if(ncolltests.gt.10)then
596
write(*,*)'Collinear tests done for (Born) config', iconfig
597
write(*,*)'Failures:',nerr
598
fail_frac= nerr/dble(ncolltests)
599
if (fail_frac.lt.max_fail) then
600
write(*,501) nFKSprocess, fail_frac
602
write(*,502) nFKSprocess, fail_frac
608
enddo ! Loop over Born configurations
609
enddo ! Loop over nFKSprocess
613
401 format(' Soft test ',i2,' PASSED. Fraction of failures: ',
615
402 format(' Soft test ',I2,' FAILED. Fraction of failures: ',
617
501 format('Collinear test ',i2,' PASSED. Fraction of failures: ',
619
502 format('Collinear test ',I2,' FAILED. Fraction of failures: ',
628
subroutine clear_events()
632
subroutine store_events()
634
integer function n_unwgted()
638
subroutine outfun(pp,www)
640
include 'nexternal.inc'
641
real*8 pp(0:3,nexternal),www
643
write(*,*)'This routine should not be called here'