~madteam/mg5amcnlo/series2.0

« back to all changes in this revision

Viewing changes to vendor/IREGI/src/cti_reduce.f90

  • Committer: olivier Mattelaer
  • Date: 2015-03-05 00:14:16 UTC
  • mfrom: (258.1.9 2.3)
  • mto: (258.8.1 2.3)
  • mto: This revision was merged to the branch mainline in revision 259.
  • Revision ID: olivier.mattelaer@uclouvain.be-20150305001416-y9mzeykfzwnl9t0j
partial merge

Show diffs side-by-side

added added

removed removed

Lines of Context:
26
26
    LOGICAL,INTENT(OUT)::TSTABLE
27
27
    INTEGER,DIMENSION(1,1)::sol11
28
28
    REAL(KIND(1d0)),DIMENSION(1)::factor1
29
 
    INTEGER::first=0,i,j,numzerp,n,r,nr,init,ntot
 
29
    INTEGER::first=0,i,j,jk,numzerp,n,r,nr,init,ntot
30
30
    INTEGER,DIMENSION(NLOOPLINE)::zerp
31
 
    REAL(KIND(1d0)),DIMENSION(xiarraymax2)::syfactor
32
 
    INTEGER,DIMENSION(xiarraymax2,-1:NLOOPLINE)::sy
 
31
    ! f[i,n]=(i+n-1)!/(n-1)!/i!
 
32
    ! nmax=5,rmax=6,NLOOPLINE=nmax+1
 
33
    ! see syntensor in ti_reduce.f90
 
34
    ! xiarraymax2=(f[0,nmax])+(f[1,nmax])+(f[2,nmax]+f[0,nmax])
 
35
    ! +(f[3,nmax]+f[1,nmax])+(f[4,nmax]+f[2,nmax]+f[0,nmax])
 
36
    ! +(f[5,nmax]+f[3,nmax]+f[1,nmax])+(f[6,nmax]+f[4,nmax]+f[2,nmax]+f[0,nmax])
 
37
    ! when rmax=5,nmax=5 -> xiarraymax2=314
 
38
    ! when rmax=6,nmax=5 -> xiarraymax2=610
 
39
    INTEGER::xiarraymax2,xiarraymax,xiarraymax3
 
40
    !REAL(KIND(1d0)),DIMENSION(xiarraymax2)::syfactor
 
41
    REAL(KIND(1d0)),DIMENSION(:),ALLOCATABLE::syfactor
 
42
    !INTEGER,DIMENSION(xiarraymax2,-1:NLOOPLINE)::sy
 
43
    INTEGER,DIMENSION(:,:),ALLOCATABLE::sy
33
44
!    REAL(KIND(1d0)),PARAMETER::EPS=1d-10
34
45
    REAL(KIND(1d0))::temp
35
46
    INTEGER::ctmode_local
36
47
    LOGICAL::do_shift
37
 
    SAVE first
 
48
    SAVE first,xiarraymax,xiarraymax2,xiarraymax3,syfactor,sy
38
49
    IF(PRESENT(CTMODE))THEN
39
50
       ctmode_local=CTMODE
40
51
    ELSE
43
54
    MU_R_IREGI=MU
44
55
    IF(first.EQ.0)THEN
45
56
       IF(.NOT.print_banner)THEN
46
 
       WRITE(*,*)"#####################################################################################"
47
 
       WRITE(*,*)"#                                                                                   #"
48
 
       WRITE(*,*)"#                           IREGI-alpha-1.0.0                                       #"
49
 
       WRITE(*,*)"#    package for one-loop tensor Integral REduction with General propagator Indices #"
50
 
       WRITE(*,*)"#    Author: Hua-Sheng Shao (huasheng.shao@cern.ch/erdissshaw@gmail.com)            #"
51
 
       WRITE(*,*)"#    a) Physics School, Peking University, Beijing, China                           #"
52
 
       WRITE(*,*)"#    b) PH Department, TH Unit, CERN, Geneva, Switzerland                           #"
53
 
       WRITE(*,*)"#                                                                                   #"
54
 
       WRITE(*,*)"#####################################################################################"
55
 
       print_banner=.TRUE.
 
57
          INCLUDE "banner.inc"
 
58
          print_banner=.TRUE.
 
59
       ENDIF
 
60
       xiarraymax2=0
 
61
       DO i=0,MAXRANK_IREGI
 
62
          j=i/2+1
 
63
          DO jk=1,j
 
64
             xiarraymax2=xiarraymax2+&
 
65
                  xiarray_arg1(2*jk-2+MOD(i,2),MAXNLOOP_IREGI-1)
 
66
          ENDDO
 
67
       ENDDO
 
68
       xiarraymax=xiarray_arg1(MAXRANK_IREGI,MAXNLOOP_IREGI-1)
 
69
       xiarraymax3=xiarray_arg1(MAXRANK_IREGI,4)
 
70
       IF(.NOT.ALLOCATED(syfactor))THEN
 
71
          ALLOCATE(syfactor(xiarraymax2))
 
72
       ENDIF
 
73
       IF(.NOT.ALLOCATED(sy))THEN
 
74
          ALLOCATE(sy(xiarraymax2,-1:MAXNLOOP_IREGI))
56
75
       ENDIF
57
76
       ! initialization xiarray and metric
58
77
       CALL all_Integers(1,1,1,sol11,factor1)
99
118
       ENDIF
100
119
    ENDDO
101
120
    n=NLOOPLINE-numzerp
102
 
    IF(n.GE.6.OR.MAXRANK.GT.6)THEN
103
 
       WRITE(*,*)"ERROR: out of range of comp_tensor_integral_reduce (N<7,R<7)"
 
121
    IF(n.GE.MAXNLOOP_IREGI.OR.MAXRANK.GT.MAXRANK_IREGI)THEN
 
122
       WRITE(*,100)"ERROR: out of range of comp_tensor_integral_reduce (N<=",MAXNLOOP_IREGI,",R<=",MAXRANK_IREGI,")"
104
123
       STOP
105
124
    ENDIF
106
 
    CALL sytensor(n,MAXRANK,ntot,sy(1:xiarraymax2,-1:n),syfactor(1:xiarraymax2))
107
 
    CALL comp_symmetry(IMODE,NLOOPLINE,ntot,sy(1:xiarraymax2,-1:n),syfactor(1:xiarraymax2),&
 
125
    CALL sytensor(n,MAXRANK,xiarraymax,xiarraymax2,ntot,&
 
126
         sy(1:xiarraymax2,-1:n),syfactor(1:xiarraymax2))
 
127
    CALL comp_symmetry(IMODE,NLOOPLINE,ntot,xiarraymax,xiarraymax2,xiarraymax3,&
 
128
         sy(1:xiarraymax2,-1:n),syfactor(1:xiarraymax2),&
108
129
         numzerp,zerp,PDEN,M2L,NCOEFS,TICOEFS)
109
130
    TSTABLE=STABLE_IREGI
110
131
    IF(TSTABLE.AND.do_shift)THEN
114
135
       CALL SHIFT_MOM(PDEN_N,NCOEFS,TICOEFS_SAVE,TICOEFS)
115
136
    ENDIF
116
137
    RETURN
 
138
100 FORMAT(2X,A55,I2,A4,I2,A1)
117
139
  END SUBROUTINE comp_tensor_integral_reduce
118
140
 
119
 
  SUBROUTINE comp_symmetry(IMODE,NLOOPLINE,ntot,sy,syfactor,numzerp,zerp,PDEN,M2L,NCOEFS,coefs)
 
141
  SUBROUTINE comp_symmetry(IMODE,NLOOPLINE,ntot,xiarraymax,xiarraymax2,xiarraymax3,&
 
142
       sy,syfactor,numzerp,zerp,PDEN,M2L,NCOEFS,coefs)
120
143
    ! IMODE=0, IBP reduction  
121
144
    ! IMODE=1, PaVe reduction
122
145
    ! IMODE=2, PaVe reduction with stablility improved by IBP reduction 
123
146
    IMPLICIT NONE
124
 
    INTEGER,INTENT(IN)::ntot,numzerp,NLOOPLINE,IMODE,NCOEFS
 
147
    INTEGER,INTENT(IN)::ntot,numzerp,NLOOPLINE,IMODE,NCOEFS,xiarraymax,xiarraymax2,xiarraymax3
125
148
    INTEGER,DIMENSION(xiarraymax2,-1:NLOOPLINE-numzerp),INTENT(IN)::sy
 
149
    !INTEGER,DIMENSION(*,*),INTENT(IN)::sy
126
150
    REAL(KIND(1d0)),DIMENSION(xiarraymax2),INTENT(IN)::syfactor
127
151
    INTEGER,DIMENSION(NLOOPLINE),INTENT(IN)::zerp
128
152
    REAL(KIND(1d0)),DIMENSION(NLOOPLINE,0:3),INTENT(IN)::PDEN
137
161
    INTEGER::i,j,k,init,r,nntot,nt
138
162
    LOGICAL::lzero
139
163
    INTEGER,DIMENSION(0:NLOOPLINE)::nj
140
 
    INTEGER,DIMENSION(0:NLOOPLINE,1:5)::jlist
 
164
    INTEGER,DIMENSION(0:NLOOPLINE,1:MAXRANK_IREGI)::jlist
141
165
    REAL(KIND(1d0))::res
142
166
    COMPLEX(KIND(1d0)),DIMENSION(1:4)::scalar
 
167
    !REAL(KIND(1d0)),DIMENSION(xiarraymax)::coco
143
168
    REAL(KIND(1d0)),DIMENSION(xiarraymax3)::coco
144
 
    INTEGER,DIMENSION(5)::lor
 
169
    INTEGER,DIMENSION(MAXRANK_IREGI)::lor
145
170
    coefs(0:NCOEFS-1,1:4)=DCMPLX(0d0)
146
171
    j=1
147
172
    DO i=1,NLOOPLINE
177
202
             init=init+sol(j,4)
178
203
          ENDIF
179
204
          nj(0:nt)=0
180
 
          jlist(0:nt,1:5)=0
 
205
          jlist(0:nt,1:MAXRANK_IREGI)=0
181
206
          CALL recursive_symmetry(nt,1,r,PCLL(1:nt,0:3),sy(i,-1:nt),lor(1:r),&
182
 
               nj(0:nt),jlist(0:nt,1:5),res)
 
207
               nj(0:nt),jlist(0:nt,1:MAXRANK_IREGI),res)
183
208
          lzero=lzero.AND.(ABS(res).LT.EPS)
184
209
          IF(.NOT.ML5_CONVENTION)THEN
185
210
             coco(j)=res*factor(j)