~maddevelopers/mg5amcnlo/2.7.1.3

« back to all changes in this revision

Viewing changes to tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/gg_wmtbx/CT_interface.f

  • Committer: olivier Mattelaer
  • Date: 2016-05-12 11:00:18 UTC
  • mfrom: (262.1.150 2.3.4)
  • Revision ID: olivier.mattelaer@uclouvain.be-20160512110018-sevb79f0wm4g8mpp
pass to 2.4.0

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
C     ===========================================
 
2
C     ===== Beginning of CutTools Interface =====
 
3
C     ===========================================
1
4
      SUBROUTINE ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,RES,STABLE)
2
5
C     
3
6
C     Generated by MadGraph5_aMC@NLO v. %(version)s, %(date)s
6
9
C     
7
10
C     Interface between MG5 and CutTools.
8
11
C     
9
 
C     Process: g g > w- t b~ QED=1 QCD=2 [ virt = QCD ]
 
12
C     Process: g g > w- t b~ QED<=1 QCD<=2 [ virt = QCD ]
10
13
C     
11
14
C     
12
15
C     CONSTANTS 
23
26
C     
24
27
      INTEGER NLOOPLINE, RANK
25
28
      REAL*8 PL(0:3,NLOOPLINE)
26
 
      REAL*8 PCT(0:3,0:NLOOPLINE-1)
 
29
      REAL*8 PCT(0:3,0:NLOOPLINE-1),ABSPCT(0:3)
 
30
      REAL*8 REF_P
27
31
      COMPLEX*16 M2L(NLOOPLINE)
28
32
      COMPLEX*16 M2LCT(0:NLOOPLINE-1)
29
33
      COMPLEX*16 RES(3)
33
37
C     
34
38
      COMPLEX*16 R1, ACC
35
39
      INTEGER I, J, K
36
 
      LOGICAL CTINIT, TIRINIT, GOLEMINIT
37
 
      COMMON/REDUCTIONCODEINIT/CTINIT,TIRINIT,GOLEMINIT
 
40
      LOGICAL CTINIT, TIRINIT, GOLEMINIT, SAMURAIINIT, NINJAINIT
 
41
      COMMON/REDUCTIONCODEINIT/CTINIT,TIRINIT,GOLEMINIT,SAMURAIINIT
 
42
     $ ,NINJAINIT
38
43
C     
39
44
C     EXTERNAL FUNCTIONS
40
45
C     
71
76
 
72
77
C     CONVERT THE MOMENTA FLOWING IN THE LOOP LINES TO CT CONVENTIONS
73
78
      DO I=0,3
 
79
        ABSPCT(I)=0.D0
74
80
        DO J=0,(NLOOPLINE-1)
75
81
          PCT(I,J)=0.D0
76
82
        ENDDO
78
84
      DO I=0,3
79
85
        DO J=1,NLOOPLINE
80
86
          PCT(I,0)=PCT(I,0)+PL(I,J)
 
87
          ABSPCT(I)=ABSPCT(I)+ABS(PL(I,J))
81
88
        ENDDO
82
89
      ENDDO
83
 
      IF (CHECKPCONSERVATION) THEN
84
 
        IF (PCT(0,0).GT.1.D-6) THEN
85
 
          WRITE(*,*) 'energy is not conserved ',PCT(0,0)
86
 
          STOP 'energy is not conserved'
87
 
        ELSEIF (PCT(1,0).GT.1.D-6) THEN
88
 
          WRITE(*,*) 'px is not conserved ',PCT(1,0)
89
 
          STOP 'px is not conserved'
90
 
        ELSEIF (PCT(2,0).GT.1.D-6) THEN
91
 
          WRITE(*,*) 'py is not conserved ',PCT(2,0)
92
 
          STOP 'py is not conserved'
93
 
        ELSEIF (PCT(3,0).GT.1.D-6) THEN
94
 
          WRITE(*,*) 'pz is not conserved ',PCT(3,0)
95
 
          STOP 'pz is not conserved'
 
90
      REF_P = MAX(ABSPCT(0), ABSPCT(1),ABSPCT(2),ABSPCT(3))
 
91
      DO I=0,3
 
92
        ABSPCT(I) = MAX(REF_P*1E-6, ABSPCT(I))
 
93
      ENDDO
 
94
      IF (CHECKPCONSERVATION.AND.REF_P.GT.1D-8) THEN
 
95
        IF ((PCT(0,0)/ABSPCT(0)).GT.1.D-6) THEN
 
96
          WRITE(*,*) 'energy is not conserved (flag:CT95)',PCT(0,0)
 
97
          STOP 'energy is not conserved (flag:CT95)'
 
98
        ELSEIF ((PCT(1,0)/ABSPCT(1)).GT.1.D-6) THEN
 
99
          WRITE(*,*) 'px is not conserved (flag:CT95)',PCT(1,0)
 
100
          STOP 'px is not conserved (flag:CT95)'
 
101
        ELSEIF ((PCT(2,0)/ABSPCT(2)).GT.1.D-6) THEN
 
102
          WRITE(*,*) 'py is not conserved (flag:CT95)',PCT(2,0)
 
103
          STOP 'py is not conserved (flag:CT95)'
 
104
        ELSEIF ((PCT(3,0)/ABSPCT(3)).GT.1.D-6) THEN
 
105
          WRITE(*,*) 'pz is not conserved (flag:CT95)',PCT(3,0)
 
106
          STOP 'pz is not conserved (flag:CT95)'
96
107
        ENDIF
97
108
      ENDIF
98
109
      DO I=0,3
108
119
      RES(1)=NORMALIZATION*2.0D0*DBLE(RES(1))
109
120
      RES(2)=NORMALIZATION*2.0D0*DBLE(RES(2))
110
121
      RES(3)=NORMALIZATION*2.0D0*DBLE(RES(3))
111
 
C     WRITE(*,*) 'Loop ID',ID,' =',RES(1),RES(2),RES(3)
 
122
C     WRITE(*,*) 'CutTools: Loop ID',ID,' =',RES(1),RES(2),RES(3)
112
123
      END
113
124
 
114
125
      SUBROUTINE ML5_0_INITCT()
151
162
 
152
163
      END
153
164
 
154
 
      SUBROUTINE ML5_0_LOOP_5(W1, W2, W3, W4, W5, M1, M2, M3, M4, M5
155
 
     $ ,  RANK, SQUAREDSOINDEX, LOOPNUM)
 
165
      SUBROUTINE ML5_0_BUILD_KINEMATIC_MATRIX(NLOOPLINE,P_LOOP,M2L
 
166
     $ ,S_MAT)
 
167
C     
 
168
C     Helper function that compute the loop kinematic matrix with
 
169
C      proper thresholds
 
170
C     NLOOPLINE : Number of loop lines
 
171
C     P_LOOP    : List of external momenta running in the loop, i.e.
 
172
C      q_i in the denominator (l_i+q_i)**2-m_i**2
 
173
C     M2L       : List of complex-valued masses running in the loop.
 
174
C     S_MAT(N,N): Kinematic matrix output.
 
175
C     
 
176
C     ARGUMENTS
 
177
C     
 
178
      INTEGER NLOOPLINE
 
179
      REAL*8 P_LOOP(NLOOPLINE,0:3)
 
180
      COMPLEX*16 M2L(NLOOPLINE)
 
181
      COMPLEX*16 S_MAT(NLOOPLINE,NLOOPLINE)
 
182
C     
 
183
C     GLOBAL VARIABLES
 
184
C     
 
185
      INCLUDE 'MadLoopParams.inc'
 
186
C     
 
187
C     LOCAL VARIABLES
 
188
C     
 
189
      INTEGER I,J,K
 
190
      COMPLEX*16 DIFFSQ
 
191
      REAL*8 REF_NORMALIZATION
 
192
 
 
193
C     ----------
 
194
C     BEGIN CODE
 
195
C     ----------
 
196
 
 
197
      DO I=1,NLOOPLINE
 
198
        DO J=1,NLOOPLINE
 
199
 
 
200
          IF(I.EQ.J)THEN
 
201
            S_MAT(I,J)=-(M2L(I)+M2L(J))
 
202
          ELSE
 
203
            DIFFSQ = (DCMPLX(P_LOOP(I,0),0.0D0)-DCMPLX(P_LOOP(J,0)
 
204
     $       ,0.0D0))**2
 
205
            DO K=1,3
 
206
              DIFFSQ = DIFFSQ - (DCMPLX(P_LOOP(I,K),0.0D0)
 
207
     $         -DCMPLX(P_LOOP(J,K),0.0D0))**2
 
208
            ENDDO
 
209
C           Default value of the kinematic matrix
 
210
            S_MAT(I,J)=DIFFSQ-M2L(I)-M2L(J)
 
211
C           And we now test various thresholds. Normaly, at most one
 
212
C            applies.
 
213
            IF(ABS(M2L(I)).NE.0.0D0)THEN
 
214
              IF(ABS((DIFFSQ-M2L(I))/M2L(I)).LT.OSTHRES)THEN
 
215
                S_MAT(I,J)=-M2L(J)
 
216
              ENDIF
 
217
            ENDIF
 
218
            IF(ABS(M2L(J)).NE.0.0D0)THEN
 
219
              IF(ABS((DIFFSQ-M2L(J))/M2L(J)).LT.OSTHRES)THEN
 
220
                S_MAT(I,J)=-M2L(I)
 
221
              ENDIF
 
222
            ENDIF
 
223
C           Choose what seems the most appropriate way to compare
 
224
C           massless onshellness.
 
225
            REF_NORMALIZATION=0.0D0
 
226
C           Here, we chose to base the threshold only on the energy
 
227
C            component
 
228
            DO K=0,0
 
229
              REF_NORMALIZATION = REF_NORMALIZATION + ABS(P_LOOP(I,K))
 
230
     $          + ABS(P_LOOP(J,K))
 
231
            ENDDO
 
232
            REF_NORMALIZATION = (REF_NORMALIZATION/2.0D0)**2
 
233
            IF(REF_NORMALIZATION.NE.0.0D0)THEN
 
234
              IF(ABS(DIFFSQ/REF_NORMALIZATION).LT.OSTHRES)THEN
 
235
                S_MAT(I,J)=-(M2L(I)+M2L(J))
 
236
              ENDIF
 
237
            ENDIF
 
238
          ENDIF
 
239
 
 
240
        ENDDO
 
241
      ENDDO
 
242
 
 
243
      END
 
244
 
 
245
      SUBROUTINE ML5_0_MP_BUILD_KINEMATIC_MATRIX(NLOOPLINE,P_LOOP,M2L
 
246
     $ ,S_MAT)
 
247
C     
 
248
C     Helper function that compute the loop kinematic matrix with
 
249
C      proper thresholds
 
250
C     NLOOPLINE : Number of loop lines
 
251
C     P_LOOP    : List of external momenta running in the loop, i.e.
 
252
C      q_i in the denominator (l_i+q_i)**2-m_i**2
 
253
C     M2L       : List of complex-valued masses running in the loop.
 
254
C     S_MAT(N,N): Kinematic matrix output.
 
255
C     
 
256
C     ARGUMENTS
 
257
C     
 
258
      INTEGER NLOOPLINE
 
259
      REAL*16 P_LOOP(NLOOPLINE,0:3)
 
260
      COMPLEX*32 M2L(NLOOPLINE)
 
261
      COMPLEX*32 S_MAT(NLOOPLINE,NLOOPLINE)
 
262
C     
 
263
C     GLOBAL VARIABLES
 
264
C     
 
265
      INCLUDE 'MadLoopParams.inc'
 
266
C     
 
267
C     LOCAL VARIABLES
 
268
C     
 
269
      INTEGER I,J,K
 
270
      COMPLEX*32 DIFFSQ
 
271
      REAL*16 REF_NORMALIZATION
 
272
 
 
273
C     ----------
 
274
C     BEGIN CODE
 
275
C     ----------
 
276
 
 
277
      DO I=1,NLOOPLINE
 
278
        DO J=1,NLOOPLINE
 
279
 
 
280
          IF(I.EQ.J)THEN
 
281
            S_MAT(I,J)=-(M2L(I)+M2L(J))
 
282
          ELSE
 
283
            DIFFSQ = (CMPLX(P_LOOP(I,0),0.0E0_16,KIND=16)
 
284
     $       -CMPLX(P_LOOP(J,0),0.0E0_16,KIND=16))**2
 
285
            DO K=1,3
 
286
              DIFFSQ = DIFFSQ - (CMPLX(P_LOOP(I,K),0.0E0_16,KIND=16)
 
287
     $         -CMPLX(P_LOOP(J,K),0.0E0_16,KIND=16))**2
 
288
            ENDDO
 
289
C           Default value of the kinematic matrix
 
290
            S_MAT(I,J)=DIFFSQ-M2L(I)-M2L(J)
 
291
C           And we now test various thresholds. Normaly, at most one
 
292
C            applies.
 
293
            IF(ABS(M2L(I)).NE.0.0E0_16)THEN
 
294
              IF(ABS((DIFFSQ-M2L(I))/M2L(I)).LT.OSTHRES)THEN
 
295
                S_MAT(I,J)=-M2L(J)
 
296
              ENDIF
 
297
            ENDIF
 
298
            IF(ABS(M2L(J)).NE.0.0E0_16)THEN
 
299
              IF(ABS((DIFFSQ-M2L(J))/M2L(J)).LT.OSTHRES)THEN
 
300
                S_MAT(I,J)=-M2L(I)
 
301
              ENDIF
 
302
            ENDIF
 
303
C           Choose what seems the most appropriate way to compare
 
304
C           massless onshellness.
 
305
            REF_NORMALIZATION=0.0E0_16
 
306
C           Here, we chose to base the threshold only on the energy
 
307
C            component
 
308
            DO K=0,0
 
309
              REF_NORMALIZATION = REF_NORMALIZATION + ABS(P_LOOP(I,K))
 
310
     $          + ABS(P_LOOP(J,K))
 
311
            ENDDO
 
312
            REF_NORMALIZATION = (REF_NORMALIZATION/2.0E0_16)**2
 
313
            IF(REF_NORMALIZATION.NE.0.0E0_16)THEN
 
314
              IF(ABS(DIFFSQ/REF_NORMALIZATION).LT.OSTHRES)THEN
 
315
                S_MAT(I,J)=-(M2L(I)+M2L(J))
 
316
              ENDIF
 
317
            ENDIF
 
318
          ENDIF
 
319
 
 
320
        ENDDO
 
321
      ENDDO
 
322
 
 
323
      END
 
324
 
 
325
 
 
326
 
 
327
      SUBROUTINE ML5_0_LOOP_5(W1, W2, W3, W4, W5, M1, M2, M3, M4, M5, 
 
328
     $  RANK, SQUAREDSOINDEX, LOOPNUM)
156
329
      INTEGER    NEXTERNAL
157
330
      PARAMETER (NEXTERNAL=5)
158
331
      INTEGER    NLOOPLINE
178
351
C     LOCAL VARIABLES 
179
352
C     
180
353
      REAL*8 PL(0:3,NLOOPLINE)
 
354
      REAL*16 MP_PL(0:3,NLOOPLINE)
181
355
      COMPLEX*16 M2L(NLOOPLINE)
182
356
      INTEGER PAIRING(NLOOPLINE),WE(5)
183
357
      INTEGER I, J, K, TEMP,I_LIB
204
378
 
205
379
      COMPLEX*16 W(20,NWAVEFUNCS)
206
380
      COMMON/ML5_0_W/W
 
381
      COMPLEX*32 MP_W(20,NWAVEFUNCS)
 
382
      COMMON/ML5_0_MP_W/MP_W
 
383
 
207
384
      REAL*8 LSCALE
208
385
      INTEGER CTMODE
209
386
      COMMON/ML5_0_CT/LSCALE,CTMODE
214
391
C     BEGIN CODE
215
392
C     ----------
216
393
 
 
394
C     Determine it uses qp or not
 
395
      DOING_QP = (CTMODE.GE.4)
 
396
 
217
397
      IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED).OR.GOODAMP(SQUAREDSOIND
218
 
     $ EX,LOOPNUM)) THEN
 
398
     $EX,LOOPNUM)) THEN
219
399
        WE(1)=W1
220
400
        WE(2)=W2
221
401
        WE(3)=W3
237
417
          TEMP=1
238
418
          DO J=1,NLOOPLINE
239
419
            PL(I,J)=0.D0
 
420
            IF (DOING_QP) THEN
 
421
              MP_PL(I,J)=0.0E+0_16
 
422
            ENDIF
240
423
            DO K=TEMP,(TEMP+PAIRING(J)-1)
241
424
              PL(I,J)=PL(I,J)-DBLE(W(1+I,WE(K)))
 
425
              IF (DOING_QP) THEN
 
426
                MP_PL(I,J)=MP_PL(I,J)-REAL(MP_W(1+I,WE(K)),KIND=16)
 
427
              ENDIF
242
428
            ENDDO
243
429
            TEMP=TEMP+PAIRING(J)
244
430
          ENDDO
255
441
            EXIT
256
442
          ENDIF
257
443
        ENDDO
258
 
C       Determine it uses qp or not
259
 
        DOING_QP=.FALSE.
260
 
        IF(CTMODE.GE.4)DOING_QP=.TRUE.
261
444
C       Choose the correct loop library
262
445
        CALL ML5_0_CHOOSE_LOOPLIB(LIBINDEX,NLOOPLINE,RANK,COMPLEX_MASS
263
 
     $   ,DOING_QP,I_LIB)
 
446
     $   ,ID,DOING_QP,I_LIB)
264
447
        IF(MLREDUCTIONLIB(I_LIB).EQ.1)THEN
265
 
C         USE CUTTOOLS
266
 
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1,SQUAREDSOIN
267
 
     $     DEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
 
448
C         CutTools is used
 
449
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1
 
450
     $     ,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
268
451
        ELSE
269
 
C         USE TIR
 
452
C         Tensor Integral Reduction is used 
270
453
          CALL ML5_0_TIRLOOP(SQUAREDSOINDEX,LOOPNUM,I_LIB,NLOOPLINE,PL
271
 
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX
272
 
     $     ,LOOPNUM))
 
454
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)
 
455
     $     ,S(SQUAREDSOINDEX,LOOPNUM))
273
456
        ENDIF
274
457
      ELSE
275
458
        LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
306
489
C     LOCAL VARIABLES 
307
490
C     
308
491
      REAL*8 PL(0:3,NLOOPLINE)
 
492
      REAL*16 MP_PL(0:3,NLOOPLINE)
309
493
      COMPLEX*16 M2L(NLOOPLINE)
310
494
      INTEGER PAIRING(NLOOPLINE),WE(5)
311
495
      INTEGER I, J, K, TEMP,I_LIB
332
516
 
333
517
      COMPLEX*16 W(20,NWAVEFUNCS)
334
518
      COMMON/ML5_0_W/W
 
519
      COMPLEX*32 MP_W(20,NWAVEFUNCS)
 
520
      COMMON/ML5_0_MP_W/MP_W
 
521
 
335
522
      REAL*8 LSCALE
336
523
      INTEGER CTMODE
337
524
      COMMON/ML5_0_CT/LSCALE,CTMODE
342
529
C     BEGIN CODE
343
530
C     ----------
344
531
 
 
532
C     Determine it uses qp or not
 
533
      DOING_QP = (CTMODE.GE.4)
 
534
 
345
535
      IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED).OR.GOODAMP(SQUAREDSOIND
346
 
     $ EX,LOOPNUM)) THEN
 
536
     $EX,LOOPNUM)) THEN
347
537
        WE(1)=W1
348
538
        WE(2)=W2
349
539
        WE(3)=W3
364
554
          TEMP=1
365
555
          DO J=1,NLOOPLINE
366
556
            PL(I,J)=0.D0
 
557
            IF (DOING_QP) THEN
 
558
              MP_PL(I,J)=0.0E+0_16
 
559
            ENDIF
367
560
            DO K=TEMP,(TEMP+PAIRING(J)-1)
368
561
              PL(I,J)=PL(I,J)-DBLE(W(1+I,WE(K)))
 
562
              IF (DOING_QP) THEN
 
563
                MP_PL(I,J)=MP_PL(I,J)-REAL(MP_W(1+I,WE(K)),KIND=16)
 
564
              ENDIF
369
565
            ENDDO
370
566
            TEMP=TEMP+PAIRING(J)
371
567
          ENDDO
382
578
            EXIT
383
579
          ENDIF
384
580
        ENDDO
385
 
C       Determine it uses qp or not
386
 
        DOING_QP=.FALSE.
387
 
        IF(CTMODE.GE.4)DOING_QP=.TRUE.
388
581
C       Choose the correct loop library
389
582
        CALL ML5_0_CHOOSE_LOOPLIB(LIBINDEX,NLOOPLINE,RANK,COMPLEX_MASS
390
 
     $   ,DOING_QP,I_LIB)
 
583
     $   ,ID,DOING_QP,I_LIB)
391
584
        IF(MLREDUCTIONLIB(I_LIB).EQ.1)THEN
392
 
C         USE CUTTOOLS
393
 
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1,SQUAREDSOIN
394
 
     $     DEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
 
585
C         CutTools is used
 
586
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1
 
587
     $     ,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
395
588
        ELSE
396
 
C         USE TIR
 
589
C         Tensor Integral Reduction is used 
397
590
          CALL ML5_0_TIRLOOP(SQUAREDSOINDEX,LOOPNUM,I_LIB,NLOOPLINE,PL
398
 
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX
399
 
     $     ,LOOPNUM))
 
591
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)
 
592
     $     ,S(SQUAREDSOINDEX,LOOPNUM))
400
593
        ENDIF
401
594
      ELSE
402
595
        LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
406
599
      ENDIF
407
600
      END
408
601
 
409
 
      SUBROUTINE ML5_0_LOOP_4(W1, W2, W3, W4, M1, M2, M3, M4,  RANK
410
 
     $ , SQUAREDSOINDEX, LOOPNUM)
 
602
      SUBROUTINE ML5_0_LOOP_4(W1, W2, W3, W4, M1, M2, M3, M4,  RANK,
 
603
     $  SQUAREDSOINDEX, LOOPNUM)
411
604
      INTEGER    NEXTERNAL
412
605
      PARAMETER (NEXTERNAL=5)
413
606
      INTEGER    NLOOPLINE
433
626
C     LOCAL VARIABLES 
434
627
C     
435
628
      REAL*8 PL(0:3,NLOOPLINE)
 
629
      REAL*16 MP_PL(0:3,NLOOPLINE)
436
630
      COMPLEX*16 M2L(NLOOPLINE)
437
631
      INTEGER PAIRING(NLOOPLINE),WE(4)
438
632
      INTEGER I, J, K, TEMP,I_LIB
459
653
 
460
654
      COMPLEX*16 W(20,NWAVEFUNCS)
461
655
      COMMON/ML5_0_W/W
 
656
      COMPLEX*32 MP_W(20,NWAVEFUNCS)
 
657
      COMMON/ML5_0_MP_W/MP_W
 
658
 
462
659
      REAL*8 LSCALE
463
660
      INTEGER CTMODE
464
661
      COMMON/ML5_0_CT/LSCALE,CTMODE
469
666
C     BEGIN CODE
470
667
C     ----------
471
668
 
 
669
C     Determine it uses qp or not
 
670
      DOING_QP = (CTMODE.GE.4)
 
671
 
472
672
      IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED).OR.GOODAMP(SQUAREDSOIND
473
 
     $ EX,LOOPNUM)) THEN
 
673
     $EX,LOOPNUM)) THEN
474
674
        WE(1)=W1
475
675
        WE(2)=W2
476
676
        WE(3)=W3
490
690
          TEMP=1
491
691
          DO J=1,NLOOPLINE
492
692
            PL(I,J)=0.D0
 
693
            IF (DOING_QP) THEN
 
694
              MP_PL(I,J)=0.0E+0_16
 
695
            ENDIF
493
696
            DO K=TEMP,(TEMP+PAIRING(J)-1)
494
697
              PL(I,J)=PL(I,J)-DBLE(W(1+I,WE(K)))
 
698
              IF (DOING_QP) THEN
 
699
                MP_PL(I,J)=MP_PL(I,J)-REAL(MP_W(1+I,WE(K)),KIND=16)
 
700
              ENDIF
495
701
            ENDDO
496
702
            TEMP=TEMP+PAIRING(J)
497
703
          ENDDO
508
714
            EXIT
509
715
          ENDIF
510
716
        ENDDO
511
 
C       Determine it uses qp or not
512
 
        DOING_QP=.FALSE.
513
 
        IF(CTMODE.GE.4)DOING_QP=.TRUE.
514
717
C       Choose the correct loop library
515
718
        CALL ML5_0_CHOOSE_LOOPLIB(LIBINDEX,NLOOPLINE,RANK,COMPLEX_MASS
516
 
     $   ,DOING_QP,I_LIB)
 
719
     $   ,ID,DOING_QP,I_LIB)
517
720
        IF(MLREDUCTIONLIB(I_LIB).EQ.1)THEN
518
 
C         USE CUTTOOLS
519
 
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1,SQUAREDSOIN
520
 
     $     DEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
 
721
C         CutTools is used
 
722
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1
 
723
     $     ,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
521
724
        ELSE
522
 
C         USE TIR
 
725
C         Tensor Integral Reduction is used 
523
726
          CALL ML5_0_TIRLOOP(SQUAREDSOINDEX,LOOPNUM,I_LIB,NLOOPLINE,PL
524
 
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX
525
 
     $     ,LOOPNUM))
 
727
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)
 
728
     $     ,S(SQUAREDSOINDEX,LOOPNUM))
526
729
        ENDIF
527
730
      ELSE
528
731
        LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
532
735
      ENDIF
533
736
      END
534
737
 
535
 
      SUBROUTINE ML5_0_LOOP_2_3(P1, P2, W1, W2, W3, M1, M2,  RANK
536
 
     $ , SQUAREDSOINDEX, LOOPNUM)
 
738
      SUBROUTINE ML5_0_LOOP_2_3(P1, P2, W1, W2, W3, M1, M2,  RANK,
 
739
     $  SQUAREDSOINDEX, LOOPNUM)
537
740
      INTEGER    NEXTERNAL
538
741
      PARAMETER (NEXTERNAL=5)
539
742
      INTEGER    NLOOPLINE
559
762
C     LOCAL VARIABLES 
560
763
C     
561
764
      REAL*8 PL(0:3,NLOOPLINE)
 
765
      REAL*16 MP_PL(0:3,NLOOPLINE)
562
766
      COMPLEX*16 M2L(NLOOPLINE)
563
767
      INTEGER PAIRING(NLOOPLINE),WE(3)
564
768
      INTEGER I, J, K, TEMP,I_LIB
585
789
 
586
790
      COMPLEX*16 W(20,NWAVEFUNCS)
587
791
      COMMON/ML5_0_W/W
 
792
      COMPLEX*32 MP_W(20,NWAVEFUNCS)
 
793
      COMMON/ML5_0_MP_W/MP_W
 
794
 
588
795
      REAL*8 LSCALE
589
796
      INTEGER CTMODE
590
797
      COMMON/ML5_0_CT/LSCALE,CTMODE
595
802
C     BEGIN CODE
596
803
C     ----------
597
804
 
 
805
C     Determine it uses qp or not
 
806
      DOING_QP = (CTMODE.GE.4)
 
807
 
598
808
      IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED).OR.GOODAMP(SQUAREDSOIND
599
 
     $ EX,LOOPNUM)) THEN
 
809
     $EX,LOOPNUM)) THEN
600
810
        WE(1)=W1
601
811
        WE(2)=W2
602
812
        WE(3)=W3
611
821
          TEMP=1
612
822
          DO J=1,NLOOPLINE
613
823
            PL(I,J)=0.D0
 
824
            IF (DOING_QP) THEN
 
825
              MP_PL(I,J)=0.0E+0_16
 
826
            ENDIF
614
827
            DO K=TEMP,(TEMP+PAIRING(J)-1)
615
828
              PL(I,J)=PL(I,J)-DBLE(W(1+I,WE(K)))
 
829
              IF (DOING_QP) THEN
 
830
                MP_PL(I,J)=MP_PL(I,J)-REAL(MP_W(1+I,WE(K)),KIND=16)
 
831
              ENDIF
616
832
            ENDDO
617
833
            TEMP=TEMP+PAIRING(J)
618
834
          ENDDO
629
845
            EXIT
630
846
          ENDIF
631
847
        ENDDO
632
 
C       Determine it uses qp or not
633
 
        DOING_QP=.FALSE.
634
 
        IF(CTMODE.GE.4)DOING_QP=.TRUE.
635
848
C       Choose the correct loop library
636
849
        CALL ML5_0_CHOOSE_LOOPLIB(LIBINDEX,NLOOPLINE,RANK,COMPLEX_MASS
637
 
     $   ,DOING_QP,I_LIB)
 
850
     $   ,ID,DOING_QP,I_LIB)
638
851
        IF(MLREDUCTIONLIB(I_LIB).EQ.1)THEN
639
 
C         USE CUTTOOLS
640
 
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1,SQUAREDSOIN
641
 
     $     DEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
 
852
C         CutTools is used
 
853
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1
 
854
     $     ,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
642
855
        ELSE
643
 
C         USE TIR
 
856
C         Tensor Integral Reduction is used 
644
857
          CALL ML5_0_TIRLOOP(SQUAREDSOINDEX,LOOPNUM,I_LIB,NLOOPLINE,PL
645
 
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX
646
 
     $     ,LOOPNUM))
 
858
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)
 
859
     $     ,S(SQUAREDSOINDEX,LOOPNUM))
647
860
        ENDIF
648
861
      ELSE
649
862
        LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
653
866
      ENDIF
654
867
      END
655
868
 
656
 
      SUBROUTINE ML5_0_LOOP_2(W1, W2, M1, M2,  RANK, SQUAREDSOINDEX
657
 
     $ , LOOPNUM)
 
869
      SUBROUTINE ML5_0_LOOP_2(W1, W2, M1, M2,  RANK, SQUAREDSOINDEX,
 
870
     $  LOOPNUM)
658
871
      INTEGER    NEXTERNAL
659
872
      PARAMETER (NEXTERNAL=5)
660
873
      INTEGER    NLOOPLINE
680
893
C     LOCAL VARIABLES 
681
894
C     
682
895
      REAL*8 PL(0:3,NLOOPLINE)
 
896
      REAL*16 MP_PL(0:3,NLOOPLINE)
683
897
      COMPLEX*16 M2L(NLOOPLINE)
684
898
      INTEGER PAIRING(NLOOPLINE),WE(2)
685
899
      INTEGER I, J, K, TEMP,I_LIB
706
920
 
707
921
      COMPLEX*16 W(20,NWAVEFUNCS)
708
922
      COMMON/ML5_0_W/W
 
923
      COMPLEX*32 MP_W(20,NWAVEFUNCS)
 
924
      COMMON/ML5_0_MP_W/MP_W
 
925
 
709
926
      REAL*8 LSCALE
710
927
      INTEGER CTMODE
711
928
      COMMON/ML5_0_CT/LSCALE,CTMODE
716
933
C     BEGIN CODE
717
934
C     ----------
718
935
 
 
936
C     Determine it uses qp or not
 
937
      DOING_QP = (CTMODE.GE.4)
 
938
 
719
939
      IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED).OR.GOODAMP(SQUAREDSOIND
720
 
     $ EX,LOOPNUM)) THEN
 
940
     $EX,LOOPNUM)) THEN
721
941
        WE(1)=W1
722
942
        WE(2)=W2
723
943
        M2L(1)=M2**2
733
953
          TEMP=1
734
954
          DO J=1,NLOOPLINE
735
955
            PL(I,J)=0.D0
 
956
            IF (DOING_QP) THEN
 
957
              MP_PL(I,J)=0.0E+0_16
 
958
            ENDIF
736
959
            DO K=TEMP,(TEMP+PAIRING(J)-1)
737
960
              PL(I,J)=PL(I,J)-DBLE(W(1+I,WE(K)))
 
961
              IF (DOING_QP) THEN
 
962
                MP_PL(I,J)=MP_PL(I,J)-REAL(MP_W(1+I,WE(K)),KIND=16)
 
963
              ENDIF
738
964
            ENDDO
739
965
            TEMP=TEMP+PAIRING(J)
740
966
          ENDDO
751
977
            EXIT
752
978
          ENDIF
753
979
        ENDDO
754
 
C       Determine it uses qp or not
755
 
        DOING_QP=.FALSE.
756
 
        IF(CTMODE.GE.4)DOING_QP=.TRUE.
757
980
C       Choose the correct loop library
758
981
        CALL ML5_0_CHOOSE_LOOPLIB(LIBINDEX,NLOOPLINE,RANK,COMPLEX_MASS
759
 
     $   ,DOING_QP,I_LIB)
 
982
     $   ,ID,DOING_QP,I_LIB)
760
983
        IF(MLREDUCTIONLIB(I_LIB).EQ.1)THEN
761
 
C         USE CUTTOOLS
762
 
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1,SQUAREDSOIN
763
 
     $     DEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
 
984
C         CutTools is used
 
985
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1
 
986
     $     ,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
764
987
        ELSE
765
 
C         USE TIR
 
988
C         Tensor Integral Reduction is used 
766
989
          CALL ML5_0_TIRLOOP(SQUAREDSOINDEX,LOOPNUM,I_LIB,NLOOPLINE,PL
767
 
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX
768
 
     $     ,LOOPNUM))
 
990
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)
 
991
     $     ,S(SQUAREDSOINDEX,LOOPNUM))
769
992
        ENDIF
770
993
      ELSE
771
994
        LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
802
1025
C     LOCAL VARIABLES 
803
1026
C     
804
1027
      REAL*8 PL(0:3,NLOOPLINE)
 
1028
      REAL*16 MP_PL(0:3,NLOOPLINE)
805
1029
      COMPLEX*16 M2L(NLOOPLINE)
806
1030
      INTEGER PAIRING(NLOOPLINE),WE(4)
807
1031
      INTEGER I, J, K, TEMP,I_LIB
828
1052
 
829
1053
      COMPLEX*16 W(20,NWAVEFUNCS)
830
1054
      COMMON/ML5_0_W/W
 
1055
      COMPLEX*32 MP_W(20,NWAVEFUNCS)
 
1056
      COMMON/ML5_0_MP_W/MP_W
 
1057
 
831
1058
      REAL*8 LSCALE
832
1059
      INTEGER CTMODE
833
1060
      COMMON/ML5_0_CT/LSCALE,CTMODE
838
1065
C     BEGIN CODE
839
1066
C     ----------
840
1067
 
 
1068
C     Determine it uses qp or not
 
1069
      DOING_QP = (CTMODE.GE.4)
 
1070
 
841
1071
      IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED).OR.GOODAMP(SQUAREDSOIND
842
 
     $ EX,LOOPNUM)) THEN
 
1072
     $EX,LOOPNUM)) THEN
843
1073
        WE(1)=W1
844
1074
        WE(2)=W2
845
1075
        WE(3)=W3
857
1087
          TEMP=1
858
1088
          DO J=1,NLOOPLINE
859
1089
            PL(I,J)=0.D0
 
1090
            IF (DOING_QP) THEN
 
1091
              MP_PL(I,J)=0.0E+0_16
 
1092
            ENDIF
860
1093
            DO K=TEMP,(TEMP+PAIRING(J)-1)
861
1094
              PL(I,J)=PL(I,J)-DBLE(W(1+I,WE(K)))
 
1095
              IF (DOING_QP) THEN
 
1096
                MP_PL(I,J)=MP_PL(I,J)-REAL(MP_W(1+I,WE(K)),KIND=16)
 
1097
              ENDIF
862
1098
            ENDDO
863
1099
            TEMP=TEMP+PAIRING(J)
864
1100
          ENDDO
875
1111
            EXIT
876
1112
          ENDIF
877
1113
        ENDDO
878
 
C       Determine it uses qp or not
879
 
        DOING_QP=.FALSE.
880
 
        IF(CTMODE.GE.4)DOING_QP=.TRUE.
881
1114
C       Choose the correct loop library
882
1115
        CALL ML5_0_CHOOSE_LOOPLIB(LIBINDEX,NLOOPLINE,RANK,COMPLEX_MASS
883
 
     $   ,DOING_QP,I_LIB)
 
1116
     $   ,ID,DOING_QP,I_LIB)
884
1117
        IF(MLREDUCTIONLIB(I_LIB).EQ.1)THEN
885
 
C         USE CUTTOOLS
886
 
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1,SQUAREDSOIN
887
 
     $     DEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
 
1118
C         CutTools is used
 
1119
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1
 
1120
     $     ,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
888
1121
        ELSE
889
 
C         USE TIR
 
1122
C         Tensor Integral Reduction is used 
890
1123
          CALL ML5_0_TIRLOOP(SQUAREDSOINDEX,LOOPNUM,I_LIB,NLOOPLINE,PL
891
 
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX
892
 
     $     ,LOOPNUM))
 
1124
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)
 
1125
     $     ,S(SQUAREDSOINDEX,LOOPNUM))
893
1126
        ENDIF
894
1127
      ELSE
895
1128
        LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
899
1132
      ENDIF
900
1133
      END
901
1134
 
902
 
      SUBROUTINE ML5_0_LOOP_3(W1, W2, W3, M1, M2, M3,  RANK, SQUAREDSOI
903
 
     $ NDEX, LOOPNUM)
 
1135
      SUBROUTINE ML5_0_LOOP_3(W1, W2, W3, M1, M2, M3,  RANK,
 
1136
     $  SQUAREDSOINDEX, LOOPNUM)
904
1137
      INTEGER    NEXTERNAL
905
1138
      PARAMETER (NEXTERNAL=5)
906
1139
      INTEGER    NLOOPLINE
926
1159
C     LOCAL VARIABLES 
927
1160
C     
928
1161
      REAL*8 PL(0:3,NLOOPLINE)
 
1162
      REAL*16 MP_PL(0:3,NLOOPLINE)
929
1163
      COMPLEX*16 M2L(NLOOPLINE)
930
1164
      INTEGER PAIRING(NLOOPLINE),WE(3)
931
1165
      INTEGER I, J, K, TEMP,I_LIB
952
1186
 
953
1187
      COMPLEX*16 W(20,NWAVEFUNCS)
954
1188
      COMMON/ML5_0_W/W
 
1189
      COMPLEX*32 MP_W(20,NWAVEFUNCS)
 
1190
      COMMON/ML5_0_MP_W/MP_W
 
1191
 
955
1192
      REAL*8 LSCALE
956
1193
      INTEGER CTMODE
957
1194
      COMMON/ML5_0_CT/LSCALE,CTMODE
962
1199
C     BEGIN CODE
963
1200
C     ----------
964
1201
 
 
1202
C     Determine it uses qp or not
 
1203
      DOING_QP = (CTMODE.GE.4)
 
1204
 
965
1205
      IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED).OR.GOODAMP(SQUAREDSOIND
966
 
     $ EX,LOOPNUM)) THEN
 
1206
     $EX,LOOPNUM)) THEN
967
1207
        WE(1)=W1
968
1208
        WE(2)=W2
969
1209
        WE(3)=W3
981
1221
          TEMP=1
982
1222
          DO J=1,NLOOPLINE
983
1223
            PL(I,J)=0.D0
 
1224
            IF (DOING_QP) THEN
 
1225
              MP_PL(I,J)=0.0E+0_16
 
1226
            ENDIF
984
1227
            DO K=TEMP,(TEMP+PAIRING(J)-1)
985
1228
              PL(I,J)=PL(I,J)-DBLE(W(1+I,WE(K)))
 
1229
              IF (DOING_QP) THEN
 
1230
                MP_PL(I,J)=MP_PL(I,J)-REAL(MP_W(1+I,WE(K)),KIND=16)
 
1231
              ENDIF
986
1232
            ENDDO
987
1233
            TEMP=TEMP+PAIRING(J)
988
1234
          ENDDO
999
1245
            EXIT
1000
1246
          ENDIF
1001
1247
        ENDDO
1002
 
C       Determine it uses qp or not
1003
 
        DOING_QP=.FALSE.
1004
 
        IF(CTMODE.GE.4)DOING_QP=.TRUE.
1005
1248
C       Choose the correct loop library
1006
1249
        CALL ML5_0_CHOOSE_LOOPLIB(LIBINDEX,NLOOPLINE,RANK,COMPLEX_MASS
1007
 
     $   ,DOING_QP,I_LIB)
 
1250
     $   ,ID,DOING_QP,I_LIB)
1008
1251
        IF(MLREDUCTIONLIB(I_LIB).EQ.1)THEN
1009
 
C         USE CUTTOOLS
1010
 
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1,SQUAREDSOIN
1011
 
     $     DEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
 
1252
C         CutTools is used
 
1253
          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1
 
1254
     $     ,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
1012
1255
        ELSE
1013
 
C         USE TIR
 
1256
C         Tensor Integral Reduction is used 
1014
1257
          CALL ML5_0_TIRLOOP(SQUAREDSOINDEX,LOOPNUM,I_LIB,NLOOPLINE,PL
1015
 
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX
1016
 
     $     ,LOOPNUM))
 
1258
     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)
 
1259
     $     ,S(SQUAREDSOINDEX,LOOPNUM))
1017
1260
        ENDIF
1018
1261
      ELSE
1019
1262
        LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)