~madteam/mg5amcnlo/series2.0

« back to all changes in this revision

Viewing changes to vendor/IREGI/src/funlib.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:
6
6
    IMPLICIT NONE
7
7
    INTEGER::fac
8
8
    INTEGER,INTENT(IN)::i
 
9
    INTEGER::init=0,j
 
10
    LOGICAL,DIMENSION(12)::lfacsave
 
11
    INTEGER,DIMENSION(12)::facsave
 
12
    SAVE facsave,init,lfacsave
 
13
    IF(init.EQ.0)THEN
 
14
       lfacsave(1:12)=.FALSE.
 
15
       init=1
 
16
    ENDIF
9
17
    IF(i.LT.0)THEN
10
18
       WRITE(*,*)"ERROR: i<0 in factorial with i=",i
11
19
       STOP
12
20
    ENDIF
 
21
    IF(i.GT.12)THEN
 
22
       WRITE(*,*)"ERROR: i > 12, please take long type (KIND=LINT) integer"
 
23
       STOP
 
24
    ENDIF
13
25
    IF(i.EQ.0)THEN
14
26
       fac=1
15
27
    ELSE
16
 
       fac=i*factorial(i-1)
 
28
       IF(lfacsave(i))THEN
 
29
          fac=facsave(i)
 
30
       ELSE
 
31
          fac=i*factorial(i-1)
 
32
          facsave(i)=fac
 
33
          lfacsave(i)=.TRUE.
 
34
       ENDIF
17
35
    ENDIF
18
36
  END FUNCTION factorial
19
37
 
51
69
    ENDIF
52
70
    k0=0
53
71
    DO j=0,i
54
 
       k=factorial(i-j+n-2)/factorial(n-2)/factorial(i-j)
 
72
       k=xiarray_arg1(i-j,n-1)
55
73
       sol(k0+1:k0+k,1)=j
56
74
       CALL calc_all_integers(n-1,k,i-j,sol(k0+1:k0+k,2:n))
57
75
       k0=k0+k
68
86
    INTEGER,INTENT(IN)::n,ntot,i
69
87
    INTEGER,DIMENSION(ntot,n),INTENT(OUT)::sol
70
88
    REAL(KIND(1d0)),DIMENSION(ntot),INTENT(OUT)::factor
71
 
    INTEGER::ifirst=0,j,k,ntemptot
 
89
    INTEGER::ifirst=0,j,jk,k,ntemptot
 
90
    INTEGER::maxxiarray,inum,nnum,maxxiarray1,maxxiarray2
 
91
    INTEGER::maxfactor_xiarray
72
92
    SAVE ifirst
73
93
    ! calculate xiarray_i_n first
74
94
    IF(ifirst.EQ.0)THEN
75
 
       DO j=0,10
76
 
          x1array(j,1)=j
77
 
       ENDDO
78
 
       ! n=2
79
 
       CALL calc_all_integers(2,1,0,xiarray_0_2(1:1,1:2))
80
 
       factor_xiarray_0_2(1)=1d0
81
 
       CALL calc_all_integers(2,2,1,xiarray_1_2(1:2,1:2))
82
 
       DO j=1,2
83
 
          factor_xiarray_1_2(j)=DBLE(factorial(1))
84
 
          DO k=1,2
85
 
             factor_xiarray_1_2(j)=factor_xiarray_1_2(j)/&
86
 
                  DBLE(factorial(xiarray_1_2(j,k)))
87
 
          ENDDO
88
 
       ENDDO
89
 
       CALL calc_all_integers(2,3,2,xiarray_2_2(1:3,1:2))
90
 
       DO j=1,3
91
 
          factor_xiarray_2_2(j)=DBLE(factorial(2))
92
 
          DO k=1,2
93
 
             factor_xiarray_2_2(j)=factor_xiarray_2_2(j)/&
94
 
                  DBLE(factorial(xiarray_2_2(j,k)))
95
 
          ENDDO
96
 
       ENDDO
97
 
       CALL calc_all_integers(2,4,3,xiarray_3_2(1:4,1:2))
98
 
       DO j=1,4
99
 
          factor_xiarray_3_2(j)=DBLE(factorial(3))
100
 
          DO k=1,2
101
 
             factor_xiarray_3_2(j)=factor_xiarray_3_2(j)/&
102
 
                  DBLE(factorial(xiarray_3_2(j,k)))
103
 
          ENDDO
104
 
       ENDDO
105
 
       CALL calc_all_integers(2,5,4,xiarray_4_2(1:5,1:2))
106
 
       DO j=1,5
107
 
          factor_xiarray_4_2(j)=DBLE(factorial(4))
108
 
          DO k=1,2
109
 
             factor_xiarray_4_2(j)=factor_xiarray_4_2(j)/&
110
 
                  DBLE(factorial(xiarray_4_2(j,k)))
111
 
          ENDDO
112
 
       ENDDO
113
 
       CALL calc_all_integers(2,6,5,xiarray_5_2(1:6,1:2))
114
 
       DO j=1,6
115
 
          factor_xiarray_5_2(j)=DBLE(factorial(5))
116
 
          DO k=1,2
117
 
             factor_xiarray_5_2(j)=factor_xiarray_5_2(j)/&
118
 
                  DBLE(factorial(xiarray_5_2(j,k)))
119
 
          ENDDO
120
 
       ENDDO
121
 
       CALL calc_all_integers(2,7,6,xiarray_6_2(1:7,1:2))
122
 
       DO j=1,7
123
 
          factor_xiarray_6_2(j)=DBLE(factorial(6))
124
 
          DO k=1,2
125
 
             factor_xiarray_6_2(j)=factor_xiarray_6_2(j)/&
126
 
                  DBLE(factorial(xiarray_6_2(j,k)))
127
 
          ENDDO
128
 
       ENDDO
129
 
       ntot_xiarray(0,2)=1
130
 
       ntot_xiarray(1,2)=2
131
 
       ntot_xiarray(2,2)=3
132
 
       ntot_xiarray(3,2)=4
133
 
       ntot_xiarray(4,2)=5
134
 
       ntot_xiarray(5,2)=6
135
 
       ntot_xiarray(6,2)=7
136
 
       ! n=3
137
 
       CALL calc_all_integers(3,1,0,xiarray_0_3(1:1,1:3))
138
 
       factor_xiarray_0_3(1)=1d0
139
 
       CALL calc_all_integers(3,3,1,xiarray_1_3(1:3,1:3))
140
 
       DO j=1,3
141
 
          factor_xiarray_1_3(j)=DBLE(factorial(1))
142
 
          DO k=1,3
143
 
             factor_xiarray_1_3(j)=factor_xiarray_1_3(j)/&
144
 
                  DBLE(factorial(xiarray_1_3(j,k)))
145
 
          ENDDO
146
 
       ENDDO
147
 
       CALL calc_all_integers(3,6,2,xiarray_2_3(1:6,1:3))
148
 
       DO j=1,6
149
 
          factor_xiarray_2_3(j)=DBLE(factorial(2))
150
 
          DO k=1,3
151
 
             factor_xiarray_2_3(j)=factor_xiarray_2_3(j)/&
152
 
                  DBLE(factorial(xiarray_2_3(j,k)))
153
 
          ENDDO
154
 
       ENDDO
155
 
       CALL calc_all_integers(3,10,3,xiarray_3_3(1:10,1:3))
156
 
       DO j=1,10
157
 
          factor_xiarray_3_3(j)=DBLE(factorial(3))
158
 
          DO k=1,3
159
 
             factor_xiarray_3_3(j)=factor_xiarray_3_3(j)/&
160
 
                  DBLE(factorial(xiarray_3_3(j,k)))
161
 
          ENDDO
162
 
       ENDDO
163
 
       CALL calc_all_integers(3,15,4,xiarray_4_3(1:15,1:3))
164
 
       DO j=1,15
165
 
          factor_xiarray_4_3(j)=DBLE(factorial(4))
166
 
          DO k=1,3
167
 
             factor_xiarray_4_3(j)=factor_xiarray_4_3(j)/&
168
 
                  DBLE(factorial(xiarray_4_3(j,k)))
169
 
          ENDDO
170
 
       ENDDO
171
 
       CALL calc_all_integers(3,21,5,xiarray_5_3(1:21,1:3))
172
 
       DO j=1,21
173
 
          factor_xiarray_5_3(j)=DBLE(factorial(5))
174
 
          DO k=1,3
175
 
             factor_xiarray_5_3(j)=factor_xiarray_5_3(j)/&
176
 
                  DBLE(factorial(xiarray_5_3(j,k)))
177
 
          ENDDO
178
 
       ENDDO
179
 
       CALL calc_all_integers(3,28,6,xiarray_6_3(1:28,1:3))
180
 
       DO j=1,28
181
 
          factor_xiarray_6_3(j)=DBLE(factorial(6))
182
 
          DO k=1,3
183
 
             factor_xiarray_6_3(j)=factor_xiarray_6_3(j)/&
184
 
                  DBLE(factorial(xiarray_6_3(j,k)))
185
 
          ENDDO
186
 
       ENDDO
187
 
       ntot_xiarray(0,3)=1
188
 
       ntot_xiarray(1,3)=3
189
 
       ntot_xiarray(2,3)=6
190
 
       ntot_xiarray(3,3)=10
191
 
       ntot_xiarray(4,3)=15
192
 
       ntot_xiarray(5,3)=21
193
 
       ntot_xiarray(6,3)=28
194
 
       ! n=4
195
 
       CALL calc_all_integers(4,1,0,xiarray_0_4(1:1,1:4))
196
 
       factor_xiarray_0_4(1)=1d0
197
 
       CALL calc_all_integers(4,4,1,xiarray_1_4(1:4,1:4))
198
 
       DO j=1,4
199
 
          factor_xiarray_1_4(j)=DBLE(factorial(1))
200
 
          DO k=1,4
201
 
             factor_xiarray_1_4(j)=factor_xiarray_1_4(j)/&
202
 
                  DBLE(factorial(xiarray_1_4(j,k)))
203
 
          ENDDO
204
 
       ENDDO
205
 
       CALL calc_all_integers(4,10,2,xiarray_2_4(1:10,1:4))
206
 
       DO j=1,10
207
 
          factor_xiarray_2_4(j)=DBLE(factorial(2))
208
 
          DO k=1,4
209
 
             factor_xiarray_2_4(j)=factor_xiarray_2_4(j)/&
210
 
                  DBLE(factorial(xiarray_2_4(j,k)))
211
 
          ENDDO
212
 
       ENDDO
213
 
       CALL calc_all_integers(4,20,3,xiarray_3_4(1:20,1:4))
214
 
       DO j=1,20
215
 
          factor_xiarray_3_4(j)=DBLE(factorial(3))
216
 
          DO k=1,4
217
 
             factor_xiarray_3_4(j)=factor_xiarray_3_4(j)/&
218
 
                  DBLE(factorial(xiarray_3_4(j,k)))
219
 
          ENDDO
220
 
       ENDDO
221
 
       CALL calc_all_integers(4,35,4,xiarray_4_4(1:35,1:4))
222
 
       DO j=1,35
223
 
          factor_xiarray_4_4(j)=DBLE(factorial(4))
224
 
          DO k=1,4
225
 
             factor_xiarray_4_4(j)=factor_xiarray_4_4(j)/&
226
 
                  DBLE(factorial(xiarray_4_4(j,k)))
227
 
          ENDDO
228
 
       ENDDO
229
 
       CALL calc_all_integers(4,56,5,xiarray_5_4(1:56,1:4))
230
 
       DO j=1,56
231
 
          factor_xiarray_5_4(j)=DBLE(factorial(5))
232
 
          DO k=1,4
233
 
             factor_xiarray_5_4(j)=factor_xiarray_5_4(j)/&
234
 
                  DBLE(factorial(xiarray_5_4(j,k)))
235
 
          ENDDO
236
 
       ENDDO
237
 
       CALL calc_all_integers(4,84,6,xiarray_6_4(1:84,1:4))
238
 
       DO j=1,84
239
 
          factor_xiarray_6_4(j)=DBLE(factorial(6))
240
 
          DO k=1,4
241
 
             factor_xiarray_6_4(j)=factor_xiarray_6_4(j)/&
242
 
                  DBLE(factorial(xiarray_6_4(j,k)))
243
 
          ENDDO
244
 
       ENDDO
245
 
       ntot_xiarray(0,4)=1
246
 
       ntot_xiarray(1,4)=4
247
 
       ntot_xiarray(2,4)=10
248
 
       ntot_xiarray(3,4)=20
249
 
       ntot_xiarray(4,4)=35
250
 
       ntot_xiarray(5,4)=56
251
 
       ntot_xiarray(6,4)=84
252
 
       ! n=5
253
 
       CALL calc_all_integers(5,1,0,xiarray_0_5(1:1,1:5))
254
 
       factor_xiarray_0_5(1)=1d0
255
 
       CALL calc_all_integers(5,5,1,xiarray_1_5(1:5,1:5))
256
 
       DO j=1,5
257
 
          factor_xiarray_1_5(j)=DBLE(factorial(1))
258
 
          DO k=1,5
259
 
             factor_xiarray_1_5(j)=factor_xiarray_1_5(j)/&
260
 
                  DBLE(factorial(xiarray_1_5(j,k)))
261
 
          ENDDO
262
 
       ENDDO
263
 
       CALL calc_all_integers(5,15,2,xiarray_2_5(1:15,1:5))
264
 
       DO j=1,15
265
 
          factor_xiarray_2_5(j)=DBLE(factorial(2))
266
 
          DO k=1,5
267
 
             factor_xiarray_2_5(j)=factor_xiarray_2_5(j)/&
268
 
                  DBLE(factorial(xiarray_2_5(j,k)))
269
 
          ENDDO
270
 
       ENDDO
271
 
       CALL calc_all_integers(5,35,3,xiarray_3_5(1:35,1:5))
272
 
       DO j=1,35
273
 
          factor_xiarray_3_5(j)=DBLE(factorial(3))
274
 
          DO k=1,5
275
 
             factor_xiarray_3_5(j)=factor_xiarray_3_5(j)/&
276
 
                  DBLE(factorial(xiarray_3_5(j,k)))
277
 
          ENDDO
278
 
       ENDDO
279
 
       CALL calc_all_integers(5,70,4,xiarray_4_5(1:70,1:5))
280
 
       DO j=1,70
281
 
          factor_xiarray_4_5(j)=DBLE(factorial(4))
282
 
          DO k=1,5
283
 
             factor_xiarray_4_5(j)=factor_xiarray_4_5(j)/&
284
 
                  DBLE(factorial(xiarray_4_5(j,k)))
285
 
          ENDDO
286
 
       ENDDO
287
 
       CALL calc_all_integers(5,126,5,xiarray_5_5(1:126,1:5))
288
 
       DO j=1,126
289
 
          factor_xiarray_5_5(j)=DBLE(factorial(5))
290
 
          DO k=1,5
291
 
             factor_xiarray_5_5(j)=factor_xiarray_5_5(j)/&
292
 
                  DBLE(factorial(xiarray_5_5(j,k)))
293
 
          ENDDO
294
 
       ENDDO
295
 
       CALL calc_all_integers(5,210,6,xiarray_6_5(1:210,1:5))
296
 
       DO j=1,210
297
 
          factor_xiarray_6_5(j)=DBLE(factorial(6))
298
 
          DO k=1,5
299
 
             factor_xiarray_6_5(j)=factor_xiarray_6_5(j)/&
300
 
                  DBLE(factorial(xiarray_6_5(j,k)))
301
 
          ENDDO
302
 
       ENDDO
303
 
       ntot_xiarray(0,5)=1
304
 
       ntot_xiarray(1,5)=5
305
 
       ntot_xiarray(2,5)=15
306
 
       ntot_xiarray(3,5)=35
307
 
       ntot_xiarray(4,5)=70
308
 
       ntot_xiarray(5,5)=126
309
 
       ntot_xiarray(6,5)=210
 
95
       maxxiarray=(MAXRANK_IREGI+1)*(MAXNLOOP_IREGI-1)
 
96
       IF(.NOT.ALLOCATED(xiarray))THEN
 
97
          ALLOCATE(xiarray(maxxiarray))
 
98
       ENDIF
 
99
       maxfactor_xiarray=(MAXRANK_IREGI+1)*(MAXNLOOP_IREGI-2)
 
100
       IF(.NOT.ALLOCATED(factor_xiarray))THEN
 
101
          ALLOCATE(factor_xiarray(MAXRANK_IREGI+2:maxfactor_xiarray+MAXRANK_IREGI+1))
 
102
       ENDIF
 
103
       ! x1+x2+...+xn==i
 
104
       ! C(i+n-1)^(n-1)=(i+n-1)!/(n-1)!/i!
 
105
       DO j=1,maxxiarray
 
106
          inum=MOD(j-1,MAXRANK_IREGI+1) ! i
 
107
          nnum=(j-1)/(MAXRANK_IREGI+1)+1 ! n
 
108
          maxxiarray1=xiarray_arg1(inum,nnum)
 
109
          maxxiarray2=nnum
 
110
          IF(.NOT.ALLOCATED(xiarray(j)%xiarray_i_n))THEN
 
111
             ALLOCATE(xiarray(j)%xiarray_i_n(maxxiarray1,maxxiarray2))
 
112
          ENDIF
 
113
          IF(nnum.EQ.1)THEN
 
114
             xiarray(j)%xiarray_i_n(1,1)=j-1
 
115
             CYCLE
 
116
          ENDIF
 
117
          CALL calc_all_integers(nnum,maxxiarray1,inum,&
 
118
               xiarray(j)%xiarray_i_n(1:maxxiarray1,1:maxxiarray2))
 
119
          IF(.NOT.ALLOCATED(factor_xiarray(j)%factor_xiarray_i_n))THEN
 
120
             ALLOCATE(factor_xiarray(j)%factor_xiarray_i_n(maxxiarray1))
 
121
          ENDIF
 
122
          DO jk=1,maxxiarray1
 
123
             factor_xiarray(j)%factor_xiarray_i_n(jk)=DBLE(factorial(inum))
 
124
             DO k=1,nnum
 
125
                factor_xiarray(j)%factor_xiarray_i_n(jk)=factor_xiarray(j)%factor_xiarray_i_n(jk)/&
 
126
                     DBLE(factorial(xiarray(j)%xiarray_i_n(jk,k)))
 
127
             ENDDO
 
128
          ENDDO
 
129
          ntot_xiarray(inum,nnum)=maxxiarray1
 
130
       ENDDO
310
131
       ifirst=1
311
132
    ENDIF
312
 
    IF(n.EQ.1.AND.i.GT.10)THEN
313
 
       WRITE(*,*)"ERROR:i is out of the range 10 for n=1 in all_integers"
 
133
    IF(n.GT.MAXNLOOP_IREGI-1.OR.n.LT.1)THEN
 
134
       WRITE(*,100)"ERROR: n is out of range 1<=n<=",MAXNLOOP_IREGI-1," in all_integers"
 
135
       STOP
 
136
    ENDIF
 
137
    IF(i.GT.MAXRANK_IREGI.OR.i.LT.0)THEN
 
138
       WRITE(*,100)"ERROR: r is out of range 0<=r<=",MAXRANK_IREGI," in all_integers"
314
139
       STOP
315
140
    ENDIF
316
141
    IF(n.EQ.1.AND.ntot.NE.1)THEN
317
142
       WRITE(*,*)"ERROR:ntot should be 1 when n=1 in all_integers"
318
143
       STOP
319
144
    ENDIF
320
 
    IF(n.GE.2.AND.n.LE.5.AND.i.GT.6)THEN
321
 
       WRITE(*,*)"ERROR: i is out of the range 6 for 2<=n<=5 in all_integers"
322
 
       STOP
323
 
    ENDIF
324
 
    IF(n.GE.2.AND.n.LE.5)THEN
 
145
    IF(n.GE.2.AND.n.LE.MAXNLOOP_IREGI-1)THEN
325
146
       ! Make it work in MadLoop, otherwise it is wrong
326
147
       IF(ntot.NE.ntot_xiarray(i,n))THEN
327
148
          WRITE(*,*)"ERROR: ntot is not correct in all_integers"
328
149
          STOP
329
150
       ENDIF
330
151
    ENDIF
331
 
    IF(n.GT.5.OR.n.LT.1)THEN
332
 
       WRITE(*,*)"ERROR: n is out of range 1<=n<=5 in all_integers"
333
 
       STOP
334
 
    ENDIF
 
152
    j=(n-1)*(MAXRANK_IREGI+1)+i+1
335
153
    SELECT CASE(n)
336
154
       CASE(1)
337
 
          sol(1:ntot,1:n)=x1array(i:i,1:n)
 
155
          sol(1:ntot,1:n)=xiarray(j)%xiarray_i_n(1:1,1:n)
338
156
          factor(1)=1d0
339
 
       CASE(2)
340
 
          SELECT CASE(i)
341
 
             CASE(0)
342
 
                sol(1:ntot,1:n)=xiarray_0_2(1:ntot,1:n)
343
 
                factor(1:ntot)=factor_xiarray_0_2(1:ntot)
344
 
             CASE(1)
345
 
                sol(1:ntot,1:n)=xiarray_1_2(1:ntot,1:n)
346
 
                factor(1:ntot)=factor_xiarray_1_2(1:ntot)
347
 
             CASE(2)
348
 
                sol(1:ntot,1:n)=xiarray_2_2(1:ntot,1:n)
349
 
                factor(1:ntot)=factor_xiarray_2_2(1:ntot)
350
 
             CASE(3)
351
 
                sol(1:ntot,1:n)=xiarray_3_2(1:ntot,1:n)
352
 
                factor(1:ntot)=factor_xiarray_3_2(1:ntot)
353
 
             CASE(4)
354
 
                sol(1:ntot,1:n)=xiarray_4_2(1:ntot,1:n)
355
 
                factor(1:ntot)=factor_xiarray_4_2(1:ntot)
356
 
             CASE(5)
357
 
                sol(1:ntot,1:n)=xiarray_5_2(1:ntot,1:n)
358
 
                factor(1:ntot)=factor_xiarray_5_2(1:ntot)
359
 
             CASE(6)
360
 
                sol(1:ntot,1:n)=xiarray_6_2(1:ntot,1:n)
361
 
                factor(1:ntot)=factor_xiarray_6_2(1:ntot)
362
 
          END SELECT
363
 
       CASE(3)
364
 
          SELECT CASE(i)
365
 
             CASE(0)
366
 
                sol(1:ntot,1:n)=xiarray_0_3(1:ntot,1:n)
367
 
                factor(1:ntot)=factor_xiarray_0_3(1:ntot)
368
 
             CASE(1)
369
 
                sol(1:ntot,1:n)=xiarray_1_3(1:ntot,1:n)
370
 
                factor(1:ntot)=factor_xiarray_1_3(1:ntot)
371
 
             CASE(2)
372
 
                sol(1:ntot,1:n)=xiarray_2_3(1:ntot,1:n)
373
 
                factor(1:ntot)=factor_xiarray_2_3(1:ntot)
374
 
             CASE(3)
375
 
                sol(1:ntot,1:n)=xiarray_3_3(1:ntot,1:n)
376
 
                factor(1:ntot)=factor_xiarray_3_3(1:ntot)
377
 
             CASE(4)
378
 
                sol(1:ntot,1:n)=xiarray_4_3(1:ntot,1:n)
379
 
                factor(1:ntot)=factor_xiarray_4_3(1:ntot)
380
 
             CASE(5)
381
 
                sol(1:ntot,1:n)=xiarray_5_3(1:ntot,1:n)
382
 
                factor(1:ntot)=factor_xiarray_5_3(1:ntot)
383
 
             CASE(6)
384
 
                sol(1:ntot,1:n)=xiarray_6_3(1:ntot,1:n)
385
 
                factor(1:ntot)=factor_xiarray_6_3(1:ntot)
386
 
          END SELECT
387
 
       CASE(4)
388
 
          SELECT CASE(i)
389
 
             CASE(0)
390
 
                sol(1:ntot,1:n)=xiarray_0_4(1:ntot,1:n)
391
 
                factor(1:ntot)=factor_xiarray_0_4(1:ntot)
392
 
             CASE(1)
393
 
                sol(1:ntot,1:n)=xiarray_1_4(1:ntot,1:n)
394
 
                factor(1:ntot)=factor_xiarray_1_4(1:ntot)
395
 
             CASE(2)
396
 
                sol(1:ntot,1:n)=xiarray_2_4(1:ntot,1:n)
397
 
                factor(1:ntot)=factor_xiarray_2_4(1:ntot)
398
 
             CASE(3)
399
 
                sol(1:ntot,1:n)=xiarray_3_4(1:ntot,1:n)
400
 
                factor(1:ntot)=factor_xiarray_3_4(1:ntot)
401
 
             CASE(4)
402
 
                sol(1:ntot,1:n)=xiarray_4_4(1:ntot,1:n)
403
 
                factor(1:ntot)=factor_xiarray_4_4(1:ntot)
404
 
             CASE(5)
405
 
                sol(1:ntot,1:n)=xiarray_5_4(1:ntot,1:n)
406
 
                factor(1:ntot)=factor_xiarray_5_4(1:ntot)
407
 
             CASE(6)
408
 
                sol(1:ntot,1:n)=xiarray_6_4(1:ntot,1:n)
409
 
                factor(1:ntot)=factor_xiarray_6_4(1:ntot)
410
 
          END SELECT
411
 
       CASE(5)
412
 
          SELECT CASE(i)
413
 
             CASE(0)
414
 
                sol(1:ntot,1:n)=xiarray_0_5(1:ntot,1:n)
415
 
                factor(1:ntot)=factor_xiarray_0_5(1:ntot)
416
 
             CASE(1)
417
 
                sol(1:ntot,1:n)=xiarray_1_5(1:ntot,1:n)
418
 
                factor(1:ntot)=factor_xiarray_1_5(1:ntot)
419
 
             CASE(2)
420
 
                sol(1:ntot,1:n)=xiarray_2_5(1:ntot,1:n)
421
 
                factor(1:ntot)=factor_xiarray_2_5(1:ntot)
422
 
             CASE(3)
423
 
                sol(1:ntot,1:n)=xiarray_3_5(1:ntot,1:n)
424
 
                factor(1:ntot)=factor_xiarray_3_5(1:ntot)
425
 
             CASE(4)
426
 
                sol(1:ntot,1:n)=xiarray_4_5(1:ntot,1:n)
427
 
                factor(1:ntot)=factor_xiarray_4_5(1:ntot)
428
 
             CASE(5)
429
 
                sol(1:ntot,1:n)=xiarray_5_5(1:ntot,1:n)
430
 
                factor(1:ntot)=factor_xiarray_5_5(1:ntot)
431
 
             CASE(6)
432
 
                sol(1:ntot,1:n)=xiarray_6_5(1:ntot,1:n)
433
 
                factor(1:ntot)=factor_xiarray_6_5(1:ntot)
434
 
          END SELECT
 
157
       CASE DEFAULT
 
158
          sol(1:ntot,1:n)=xiarray(j)%xiarray_i_n(1:ntot,1:n)
 
159
          factor(1:ntot)=factor_xiarray(j)%factor_xiarray_i_n(1:ntot)
435
160
    END SELECT
436
161
    RETURN
 
162
100 FORMAT(2X,A31,I2,A16)
437
163
  END SUBROUTINE all_Integers
438
164
 
439
165
  SUBROUTINE calc_factorial_pair
440
166
    IMPLICIT NONE
441
167
    INTEGER::i,j
442
 
    factorial_pair(1:10,0)=1d0
443
 
    DO i=1,10
444
 
       DO j=1,10
 
168
    factorial_pair(1:MAXINDICES_IREGI,0)=1d0
 
169
    DO i=1,MAXINDICES_IREGI
 
170
       DO j=1,MAXRANK_IREGI
445
171
          factorial_pair(i,j)=DBLE(factorial(i+j-1))&
446
172
               /DBLE(factorial(i-1))/DBLE(factorial(j))
447
173
       ENDDO
970
696
    ENDDO
971
697
    RETURN
972
698
  END SUBROUTINE SHIFT_MOM
 
699
 
 
700
  FUNCTION xiarray_arg1(i,n)
 
701
    IMPLICIT NONE
 
702
    INTEGER,INTENT(IN)::i,n
 
703
    INTEGER::xiarray_arg1
 
704
    INTEGER::imax,imin,j
 
705
    !INTEGER(KIND=LINT)::itemp
 
706
    IF(i+n-1.GT.12)THEN
 
707
       imax=MAX(n-1,i)
 
708
       imin=MIN(n-1,i)
 
709
       xiarray_arg1=1
 
710
       DO j=imax+1,i+n-1
 
711
          xiarray_arg1=xiarray_arg1*j
 
712
       ENDDO
 
713
       IF(imin.LE.12)THEN
 
714
          xiarray_arg1=xiarray_arg1/factorial(imin)
 
715
       ELSE
 
716
          DO j=1,imin
 
717
             xiarray_arg1=xiarray_arg1/j
 
718
          ENDDO
 
719
       ENDIF
 
720
    ELSE
 
721
       xiarray_arg1=factorial(i+n-1)/factorial(n-1)/factorial(i)
 
722
    ENDIF
 
723
    RETURN
 
724
  END FUNCTION xiarray_arg1
973
725
END MODULE FUNLIB