← Back to branch summary

~madteam/mg5amcnlo/series2.0

« back to all changes in this revision

Viewing changes to madgraph/iolibs/template_files/madevent_symmetry.f

merged with 2.3 rev286

Show diffs side-by-side

added added

removed removed

Lines of Context:
madgraph/iolibs/template_files/madevent_symmetry.f
60
60
c      write(*,*) 'Enter compression (0=none, 1=sym, 2=BW, 3=full)'
61
61
c      read(*,*) icomp
62
62
c      if (icomp .gt. 3 .or. icomp .lt. 0) icomp=0
63
 
      if (icomp .eq. 0) then
64
 
         write(*,*) 'No compression, summing every diagram and ',
65
 
     $        'every B.W.'
66
 
      elseif (icomp .eq. 1) then
67
 
         write(*,*) 'Using symmetry but summing every B.W. '
68
 
      elseif (icomp .eq. 2) then
69
 
         write(*,*) 'Assuming B.W. but summing every diagram. '
70
 
      elseif (icomp .eq. 3) then
71
 
         write(*,*) 'Full compression. Using symmetry and assuming B.W.'
72
 
      else
73
 
         write(*,*) 'Unknown compression',icomp
74
 
         stop
75
 
      endif
 
63
c      if (icomp .eq. 0) then
 
64
c         write(*,*) 'No compression, summing every diagram and ',
 
65
c     $        'every B.W.'
 
66
c      elseif (icomp .eq. 1) then
 
67
c         write(*,*) 'Using symmetry but summing every B.W. '
 
68
c      elseif (icomp .eq. 2) then
 
69
c         write(*,*) 'Assuming B.W. but summing every diagram. '
 
70
c      elseif (icomp .eq. 3) then
 
71
c         write(*,*) 'Full compression. Using symmetry and assuming B.W.'
 
72
c      else
 
73
c         write(*,*) 'Unknown compression',icomp
 
74
c         stop
 
75
c      endif
76
76
 
77
77
      call load_para(npara,param,value)
78
 
      call get_integer(npara,param,value," nhel ",nhel_survey,0)
79
 
c     If different card options set for nhel_refine and nhel_survey:
80
 
      call get_integer(npara,param,value," nhel_survey ",nhel_survey,
81
 
     $     1*nhel_survey)
82
78
      call get_logical(npara,param,value," gridpack ",gridpack,.false.)
83
79
      call get_real(npara,param,value," bwcutoff ",bwcutoff,5d0)
84
80
      call get_real(npara,param,value," ebeam1 ",ebeam(1),0d0)
106
102
         stot=m1**2+m2**2+2*(pi1(0)*pi2(0)-pi1(3)*pi2(3))
107
103
      endif
108
104
 
109
 
      write(*,*) 'Read parameters:'
110
 
      write(*,*) 'nhel_survey = ',nhel_survey
111
 
      write(*,*) 'gridpack    = ',gridpack
112
 
      write(*,*) 'bwcutoff    = ',bwcutoff
113
 
      write(*,*) 'sqrt(stot)  = ',sqrt(stot)
114
 
 
115
 
      call printout
 
105
c      call printout
116
106
      include 'props.inc'
117
107
 
118
108
c
128
118
      enddo
129
119
      close(25)
130
120
 
131
 
      call write_input(j, nhel_survey)
132
121
      call write_bash(mapconfig,use_config,prwidth,icomp,iforest,sprop)
133
122
      end
134
123
 
135
 
      subroutine write_input(nconfigs, nhel_survey)
136
 
c***************************************************************************
137
 
c     Writes out input file for approximate calculation based on the
138
 
c     number of active configurations
139
 
c***************************************************************************
140
 
      implicit none
141
 
c
142
 
c     Constants
143
 
c
144
 
      include 'genps.inc'
145
 
      include 'nexternal.inc'
146
 
      include 'run_config.inc'
147
 
      integer    nhel_survey
148
 
c      integer   npoint_tot,         npoint_min
149
 
c      parameter (npoint_tot=50000, npoint_min=1000)
150
 
c
151
 
c     Arguments
152
 
c
153
 
      integer nconfigs
154
 
c
155
 
c     local
156
 
c
157
 
      integer npoints,itmax
158
 
      double precision acc
159
 
c-----
160
 
c  Begin Code
161
 
c-----
162
 
      write(*,*) 'Give npoints, max iter, and accuracy'
163
 
      read(*,*)  npoints,itmax,acc
164
 
 
165
 
      open (unit=26, file = 'input_app.txt', status='unknown',
166
 
     $     err=99)
167
 
      write(26,*) npoints,itmax,3,
168
 
     &     '     !Number of events and max and min iterations'      
169
 
      write(26,'(f8.4,a)') acc, '    !Accuracy'
170
 
      write(26,*) ' 2       !Grid Adjustment 0=none, 2=adjust'
171
 
      write(26,*) ' 1       !Suppress Amplitude 1=yes'
172
 
      write(26,*) nhel_survey,'       !Helicity Sum/event 0=exact'
173
 
      close(26)
174
 
      return
175
 
 99   close(26)
176
 
      write(*,*) 'Error opening input_app.txt'
177
 
      end
178
 
 
179
 
      subroutine open_bash_file(lun)
180
 
c***********************************************************************
181
 
c     Opens bash file for looping including standard header info
182
 
c     which can be used with pbs, or stand alone
183
 
c***********************************************************************
184
 
      implicit none
185
 
c
186
 
c     Constants
187
 
c
188
 
      include 'maxparticles.inc'
189
 
      include 'run_config.inc'
190
 
c
191
 
c     Arguments
192
 
c
193
 
      integer lun
194
 
c
195
 
c     Local
196
 
c
197
 
      character*30 fname
198
 
      integer ic, npos
199
 
      character*10 formstr
200
 
 
201
 
      data ic/0/
202
 
c-----
203
 
c  Begin Code
204
 
c-----
205
 
      ic=ic+1
206
 
      fname='ajob'
207
 
c     Write ic with correct number of digits
208
 
      npos=int(dlog10(dble(ic)))+1
209
 
      write(formstr,'(a,i1,a)') '(I',npos,')'
210
 
      write(fname(5:(5+npos-1)),formstr) ic
211
 
      open (unit=lun, file = fname, status='unknown')
212
 
      write(lun,15) '#!/bin/bash'
213
 
c      write(lun,15) '#PBS -q ' // PBS_QUE
214
 
c      write(lun,15) '#PBS -o PBS.log'
215
 
c      write(lun,15) '#PBS -e PBS.err'
216
 
c      write(lun,15) 'if [[ "$PBS_O_WORKDIR" != "" ]]; then' 
217
 
c      write(lun,15) '    cd $PBS_O_WORKDIR'
218
 
c      write(lun,15) 'fi'
219
 
      write(lun,15) 'k=run1_app.log'
220
 
      write(lun,15) 'script=' // fname
221
 
c      write(lun,15) 'rm -f wait.$script >& /dev/null'
222
 
c      write(lun,15) 'touch run.$script'
223
 
      write(lun,'(a$)') 'for i in '
224
 
 15   format(a)
225
 
      end
226
124
 
227
125
      subroutine write_bash(mapconfig,use_config, prwidth, jcomp,iforest,
228
126
     $     sprop)
273
171
      open (unit=27, file = 'symfact.dat', status='unknown')
274
172
      nsym=int(dlog10(dble(mapconfig(0))))+3
275
173
 
276
 
      call open_bash_file(26)
277
174
      ic = 0      
278
175
c     ncode is number of digits needed for the code
279
176
      ncode=int(dlog10(3d0)*(max_particles-3))+1
280
177
      do i=1,mapconfig(0)
281
 
         print *,'Writing config ',mapconfig(i)
282
178
         if (use_config(i) .gt. 0) then
283
179
            call bw_conflict(i,iforest(1,-max_branch,i),lconflict,
284
180
     $           sprop(1,-max_branch,i), gForceBW(-max_branch,i))
296
192
                     if(lconflict(-j)) then
297
193
c                        write(*,*) 'Got conflict ',-nbw,j
298
194
                        iarray(nbw)=1 !Cuts on BW
299
 
                        if (nbw .gt. imax) then
300
 
                           write(*,*) 'Too many BW w conflicts',nbw,imax
301
 
                        endif
 
195
c                        if (nbw .gt. imax) then
 
196
c                           write(*,*) 'Too many BW w conflicts',nbw,imax
 
197
c                        endif
302
198
                     endif
303
199
                  endif
304
200
               enddo
309
205
               call enCode(icode,iarray,ibase,imax)
310
206
               if(failConfig(i,iarray,iforest(1,-max_branch,i),
311
207
     $           sprop(1,-max_branch,i),gForceBW(-max_branch,i))) then 
312
 
                  print *,'Skipping impossible config ',mapconfig(i),'.',icode
313
208
                  goto 100
314
209
               endif
315
210
               ic=ic+1
316
 
               if (ic .gt. ChanPerJob) then
317
 
                  call close_bash_file(26)
318
 
                  call open_bash_file(26)
319
 
                  ic = 1
320
 
               endif
321
211
               nconf=int(dlog10(dble(mapconfig(i))))+1
322
212
c               write(*,*) 'mapping',ic,mapconfig(i),icode               
323
213
               if (icode .eq. 0) then
324
214
c                 Create format string based on number of digits
325
215
                  write(formstr,'(a,i1,a)') '(I',nconf,'$)'
326
 
                  write(26,formstr) mapconfig(i)
 
216
                  write(*,formstr) mapconfig(i)
327
217
c                 Write symmetry factors
328
218
                  write(formstr2,'(a,i1,a)') '(2i',nsym,')'
329
219
                  write(27,formstr2) mapconfig(i),use_config(i)
337
227
                     write(formstr,'(a,i2,a,i1,a)') '(F',nconf+ncode+1,
338
228
     $                    '.',ncode,'$)'
339
229
                  endif
340
 
                  write(26,formstr) dconfig
 
230
                  write(*,formstr) dconfig
341
231
c                 Write symmetry factors
342
232
                  nconf=int(dlog10(dble(mapconfig(i))))+1
343
233
                  if(nconf+ncode+1.lt.10) then
350
240
                  dconfig=mapconfig(i)+icode*1d0/10**ncode
351
241
                  write(27,formstr2) dconfig,use_config(i)
352
242
               endif
353
 
               write(26,'(a$)') ' '
 
243
               write(*,'(a$)') ' '
354
244
 100           call bw_increment_array(iarray,imax,ibase,done)
355
245
            enddo
356
246
         else
358
248
            write(27,formstr2) mapconfig(i),use_config(i)
359
249
         endif
360
250
      enddo
361
 
      call close_bash_file(26)
362
251
      close(27)
363
252
      if(ic.eq.0) then
364
253
c        Stop generation with error message
367
256
         if(.not.file_exists) filename = '../' // filename
368
257
         open(unit=26,file=filename,status='unknown')
369
258
         write(26,*)'No Phase Space. Please check particle masses.'
370
 
         write(*,*)'Error: No valid channels found. ',
371
 
     $        'Please check particle masses.'
372
 
         close(26)
 
259
c         write(*,*)'Error: No valid channels found. ',
 
260
c     $        'Please check particle masses.'
373
261
      endif
374
262
      end
375
263
 
424
312
c-----
425
313
      include 'props.inc'   !Propagator mass and width information prmass,prwidth
426
314
      include 'pmass.inc'   !External particle masses
427
 
      write(*,*) 'Checking for BW in config number ',iconfig
 
315
c      write(*,*) 'Checking for BW in config number ',iconfig
428
316
c
429
317
c     Reset variables
430
318
c      
458
346
            if (xmass(-i) .gt. prmass(-i,iconfig) .and.
459
347
     $           iden_part(-i).eq.0) then !Can't be on shell, and not radiation
460
348
               lconflict(-i)=.true.
461
 
               write(*,*) "Found Conflict", iconfig,i,
462
 
     $              prmass(-i,iconfig),xmass(-i)
 
349
c               write(*,*) "Found Conflict", iconfig,i,
 
350
c     $              prmass(-i,iconfig),xmass(-i)
463
351
            endif
464
352
         endif
465
353
         if (iden_part(-i).eq.0) then
475
363
         if (lconflict(-j)) then 
476
364
            lconflict(itree(1,-j)) = .true.
477
365
            lconflict(itree(2,-j)) = .true.
478
 
            write(*,*) 'Adding conflict ',itree(1,-j),itree(2,-j)
 
366
c            write(*,*) 'Adding conflict ',itree(1,-j),itree(2,-j)
479
367
         endif
480
368
      enddo
481
369
c
482
370
c     If not enough energy, mark all BWs as conflicting
483
371
c
484
372
      if(stot.lt.mtot**2)then
485
 
         write(*,*) 'Not enough energy, set all BWs as conflicting'
 
373
c         write(*,*) 'Not enough energy, set all BWs as conflicting'
486
374
         do j=i,1,-1
487
375
            lconflict(-j) = .true.
488
376
         enddo
496
384
               lconflict(-j) = .false.
497
385
c               write(*,*) 'No conflict BW',iconfig,j
498
386
               continue
499
 
            else
500
 
               write(*,*) 'Conflicting BW',iconfig,j
 
387
c            else
 
388
c               write(*,*) 'Conflicting BW',iconfig,j
501
389
            endif
502
390
         endif
503
391
      enddo                  
588
476
         xwidth(-i)=prwidth(-i,iconfig)
589
477
         if (xwidth(-i) .gt. 0d0) then
590
478
            nbw=nbw+1
591
 
            if (iarray(nbw) .eq. 1) then
592
 
               if(xmass(-i).gt.prmass(-i,iconfig)+5d0*xwidth(-i)) then
593
 
                  failConfig=.true.
594
 
                  return
595
 
               else
596
 
                  xmass(-i)=max(xmass(-i),prmass(-i,iconfig)-5d0*xwidth(-i))
597
 
               endif
598
 
            else if(forcebw(-i) .eq. 1) then
 
479
            if(forcebw(-i) .eq. 1) then
 
480
c               if (iarray(nbw) .ne. 1) then
 
481
c                  write(*,*) "fail due to iarray", iarray(nbw)
 
482
c                  failConfig=.true.
 
483
c                  return
 
484
c               endif
599
485
               if(xmass(-i).gt.prmass(-i,iconfig)+bwcutoff*xwidth(-i)) then
600
486
                  failConfig=.true.
601
487
                  return
603
489
                  xmass(-i)=max(xmass(-i),prmass(-i,iconfig)-
604
490
     $                 bwcutoff*xwidth(-i))
605
491
               endif
 
492
            else if (iarray(nbw) .eq. 1) then
 
493
               if(xmass(-i).gt.prmass(-i,iconfig)+5d0*xwidth(-i)) then
 
494
                  failConfig=.true.
 
495
                  return
 
496
               else
 
497
                  xmass(-i)=max(xmass(-i),prmass(-i,iconfig)-5d0*xwidth(-i))
 
498
               endif
 
499
 
606
500
            endif
607
501
         endif
608
502
         mtot=mtot+xmass(-i)
617
511
      return
618
512
      end
619
513
 
620
 
      subroutine close_bash_file(lun)
621
 
c***********************************************************************
622
 
c     Opens bash file for looping including standard header info
623
 
c     which can be used with pbs, or stand alone
624
 
c***********************************************************************
625
 
      implicit none
626
 
c
627
 
c     Constants
628
 
c
629
 
c
630
 
c     Constants
631
 
c
632
 
c     Arguments
633
 
c
634
 
      integer lun
635
 
c
636
 
c     global
637
 
c
638
 
      logical gridpack
639
 
      common/to_gridpack/gridpack
640
 
c
641
 
c     local
642
 
c
643
 
      character*30 fname
644
 
      integer ic
645
 
 
646
 
      data ic/0/
647
 
c-----
648
 
c  Begin Code
649
 
c-----
650
 
 
651
 
      write(lun,'(a)') '; do'
652
 
c
653
 
c     Now write the commands
654
 
c      
655
 
c      write(lun,20) 'echo $i >& run.$script'
656
 
      write(lun,20) 'j=G$i'
657
 
      write(lun,20) 'if [[ ! -e $j ]]; then'
658
 
      write(lun,25) 'mkdir $j'
659
 
      write(lun,20) 'fi'
660
 
      write(lun,20) 'cd $j'
661
 
      write(lun,20) 'rm -f ftn25 ftn26 ftn99'
662
 
      write(lun,20) 'rm -f $k'
663
 
      write(lun,20) 'cat ../input_app.txt >& input_app.txt'
664
 
      write(lun,20) 'echo $i >> input_app.txt'
665
 
      write(lun,20) 'for((try=1;try<=10;try+=1)); '
666
 
      write(lun,20) 'do'
667
 
      write(lun,20) '../madevent > $k <input_app.txt'
668
 
      write(lun,20) 'if [ -s $k ]'
669
 
      write(lun,20) 'then'
670
 
      write(lun,20) '    break'
671
 
      write(lun,20) 'else'
672
 
      write(lun,20) 'sleep 1'
673
 
c      write(lun,20) 'rm -rf $k; ../madevent > $k <input_app.txt'
674
 
      write(lun,20) 'fi'
675
 
      write(lun,20) 'done'
676
 
      write(lun,20) 'rm -f ftn25 ftn99'
677
 
      if(.not.gridpack) write(lun,20) 'rm -f ftn26'
678
 
      write(lun,20) 'echo "" >> $k; echo "ls status:" >> $k; ls >> $k'
679
 
      write(lun,20) 'cp $k log.txt'
680
 
      write(lun,20) 'cd ../'
681
 
      write(lun,15) 'done'
682
 
 15   format(a)
683
 
 20   format(5x,a)
684
 
 25   format(10x,a)
685
 
      close(lun)
686
 
      end
687
 
 
688
 
 
689
 
 
690
514
      subroutine bw_increment_array(iarray,imax,ibase,done)
691
515
c************************************************************************
692
516
c     Increments iarray     
730
554
      done = .not. found
731
555
      end
732
556
 
733
 
      subroutine store_events()
 
557
      subroutine store_events(grid)
734
558
c**********************************************************************
735
559
c     Dummy routine
736
560
c**********************************************************************
 
561
      integer grid
737
562
      end
738
563
 
739
564
      double precision function dsig(pp,wgt,imode)