16 double complex :: mominv6(15), masses2_6(0:5), masses2_6b(0:5)
17 double complex :: mominv5(10), masses2_5(0:4)
18 double complex :: mominv4(6), masses2_4(0:3)
19 double complex :: mominv3(3), masses2_3(0:2)
20 double complex :: mominv2(1), masses2_2(0:1)
21 double complex,
allocatable :: t1coeff(:),t1coeffuv(:)
22 double complex,
allocatable :: t2coeff(:),t2coeffuv(:)
23 double complex,
allocatable :: t3coeff(:),t3coeffuv(:)
24 double complex,
allocatable :: t4coeff(:),t4coeffuv(:)
25 double complex,
allocatable :: t5coeff(:),t5coeffuv(:)
26 double complex,
allocatable :: t6coeff(:),t6coeffuv(:)
27 double precision,
allocatable :: t1err(:),t2err(:),t3err(:),t4err(:)
28 double precision,
allocatable :: t5err(:),t6err(:)
29 real :: starttime,endtime
30 integer :: rank2,rank3,rank4,rank5,rank6,nmax,rmax
31 integer :: i,k,nevents,nloops,icache,mycase
39 call init_cll(nmax,rmax,
"output_cll")
43 call initcachesystem_cll(icache,nmax)
52 allocate(t2coeff(getnc_cll(2,rank2)))
53 allocate(t2coeffuv(getnc_cll(2,rank2)))
54 allocate(t2err(0:rank2))
55 allocate(t3coeff(getnc_cll(3,rank3)))
56 allocate(t3coeffuv(getnc_cll(3,rank3)))
57 allocate(t3err(0:rank3))
58 allocate(t4coeff(getnc_cll(4,rank4)))
59 allocate(t4coeffuv(getnc_cll(4,rank4)))
60 allocate(t4err(0:rank4))
61 allocate(t5coeff(getnc_cll(5,rank5)))
62 allocate(t5coeffuv(getnc_cll(5,rank5)))
63 allocate(t5err(0:rank5))
64 allocate(t6coeff(getnc_cll(6,rank6)))
65 allocate(t6coeffuv(getnc_cll(6,rank6)))
66 allocate(t6err(0:rank6))
69 masses2_6(0) = dcmplx(1d2,-1d0)
70 masses2_6(1) = dcmplx(1d2,-1d0)
71 masses2_6(2) = dcmplx(0d0,0d0)
72 masses2_6(3) = dcmplx(1d1,0d0)
73 masses2_6(4) = dcmplx(1d1,0d0)
74 masses2_6(5) = dcmplx(0d0,0d0)
86 write(*,
'(/(a))')
' *** Demonstration of cache system ***'
87 write(*,
'(/(a))')
' NOTE: This example merely illustrates the usage of the cache system.'
88 write(*,
'((a))')
' The results are not representative for the performance'
89 write(*,
'((a))')
' of Collier, neither wrt to speed nor precision.'
90 write(*,
'(/a21,i10/)')
' Number of "Events": ',nevents
91 write(*,
'(/a39,i4/)')
' Number of "identical calls/integral": ',nloops
99 call initcachesystem_cll(icache,6)
101 write(*,
'(/(a))')
' COLI branch with cache: '
103 elseif (mycase.eq.2)
then
107 call switchoffcachesystem_cll
109 write(*,
'(/(a))')
' COLI branch without cache: '
111 else if (mycase.eq.3)
then
115 call initcachesystem_cll(icache,6)
117 write(*,
'(/(a))')
' DD branch with cache: '
119 else if (mycase.eq.4)
then
123 call switchoffcachesystem_cll
125 write(*,
'(/(a))')
' DD branch without cache '
129 call cpu_time(starttime)
133 eventloop:
do i=1,nevents
139 call initevent_cll(icache)
144 testloop:
do k=1,nloops
146 call tn_cll(t6coeff,t6coeffuv,mominv6,masses2_6,6,rank6,t6err)
151 call tn_cll(t5coeff,t5coeffuv,mominv5,masses2_5,5,rank5,t5err)
157 call tn_cll(t3coeff,t3coeffuv,mominv3,masses2_3,3,rank3,t3err)
158 call tn_cll(t4coeff,t4coeffuv,mominv4,masses2_4,4,rank4,t4err)
161 masses2_6b(0) = (0d0,0d0)
162 masses2_6b(1:5) = masses2_6(1:5)
163 call tn_cll(t6coeff,t6coeffuv,mominv6,masses2_6b,6,rank6,t6err)
173 call tn_cll(t2coeff,t2coeffuv,mominv2,masses2_2,2,rank2,t2err)
174 call tn_cll(t5coeff,t5coeffuv,mominv5,masses2_5,5,rank5,t5err)
175 call tn_cll(t4coeff,t4coeffuv,mominv4,masses2_4,4,rank4,t4err)
181 call cpu_time(endtime)
183 write(*,
'(/a29,f10.1/)')
' Needed CPU time in seconds: ',endtime-starttime
189 deallocate(t2coeff,t2coeffuv,t2err)
190 deallocate(t3coeff,t3coeffuv,t3err)
191 deallocate(t4coeff,t4coeffuv,t4err)
192 deallocate(t5coeff,t5coeffuv,t5err)
193 deallocate(t6coeff,t6coeffuv,t6err)
201 double precision p(1:6,0:3)
202 double complex MomInv6(1:15)
203 double precision :: s(0:5,0:5),fact=47d0/49d0
204 integer :: N=6,k,i,newmom,ievent
209 if (ievent.eq.1)
then
212 p(1,0) = -0.9500000000000000d+02
213 p(1,1) = 0.0000000000000000d+00
214 p(1,2) = 0.0000000000000000d+00
215 p(1,3) = 0.9500000000000000d+02
217 p(2,0) = -0.9500000000000000d+02
218 p(2,1) = 0.0000000000000000d+00
219 p(2,2) = 0.0000000000000000d+00
220 p(2,3) = -0.9500000000000000d+02
222 p(3,0) = 0.2046111001757171d+02
223 p(3,1) = 0.1057734233089455d+02
224 p(3,2) = -0.2324961261504543d+01
225 p(3,3) = 0.1736005205921753d+02
227 p(4,0) = 0.3558305163378869d+01
228 p(4,1) = 0.1436222934374051d+01
229 p(4,2) = -0.2174258125294355d+01
230 p(4,3) = -0.2423097382091398d+01
232 p(5,0) = 0.8154540918019539d+02
233 p(5,1) = -0.5230395944682889d+02
234 p(5,2) = 0.3083642435466509d+02
235 p(5,3) = 0.5443403822581044d+02
237 p(6,0) = 0.8443517563885433d+02
238 p(6,1) = 0.4029039418156027d+02
239 p(6,2) = -0.2633720496786619d+02
240 p(6,3) = -0.6937099290293661d+02
243 else if (ievent.eq.50)
then
246 p(1,0) = -250.0000000000000d0
247 p(1,1) = 0.0000000000000000d+00
248 p(1,2) = 0.0000000000000000d+00
249 p(1,3) = -250.0000000000000d0
251 p(2,0) = -250.0000000000000d0
252 p(2,1) = 0.0000000000000000d+00
253 p(2,2) = 0.0000000000000000d+00
254 p(2,3) = 250.0000000000000d0
256 p(3,0) = 239.91697591218110d0
257 p(3,1) = 31.640107125154238d0
258 p(3,2) = -116.12928664566658d0
259 p(3,3) = -207.54047252312847d0
261 p(4,0) = 10.276659578244022d0
262 p(4,1) = 1.2584629023924521d0
263 p(4,2) = 10.131487190720421d0
264 p(4,3) = -1.1742957526504370d0
266 p(5,0) = 181.18800109198216d0
267 p(5,1) = -35.098898474333737d0
268 p(5,2) = 109.26326872331617d0
269 p(5,3) = 140.20947604742324d0
271 p(6,0) = 68.618363417592718d0
272 p(6,1) = 2.2003284467870472d0
273 p(6,2) = -3.2654692683700119d0
274 p(6,3) = 68.505292228355671d0
276 else if (ievent.eq.100)
then
278 p(1,0) = -100.0000000000000d0
279 p(1,1) = 0.000000000000000d0
280 p(1,2) = 0.000000000000000d0
281 p(1,3) = -100.0000000000000d0
283 p(2,0) = -100.0000000000000d0
284 p(2,1) = 0.000000000000000d0
285 p(2,2) = 0.000000000000000d0
286 p(2,3) = 100.0000000000000d0
288 p(5,0) = 42.13045937773281d0
289 p(5,1) = 38.71045260055064d0
290 p(5,2) = 16.54764683624944d0
291 p(5,3) = 1.628450497276325d0
293 p(6,0) = 73.31699208541762d0
294 p(6,1) = -66.17895647220249d0
295 p(6,2) = 31.53024134467396d0
296 p(6,3) = -1.253367244715983d0
298 p(4,0) = 69.60509150033207d0
299 p(4,1) = 34.42307953372518d0
300 p(4,2) = -60.31848303205667d0
301 p(4,3) = 4.647683605402023d0
303 p(3,0) = 14.94745703651754d0
304 p(3,1) = -6.954575662073326d0
305 p(3,2) = 12.24059485113326d0
306 p(3,3) = -5.022766857962365d0
307 else if (ievent.eq.150)
then
309 p(1,0) = -100.0000000000000d0
310 p(1,1) = 0.0000000000000000d+00
311 p(1,2) = 0.0000000000000000d+00
312 p(1,3) = -100.0000000000000d0
314 p(2,0) = -100.0000000000000d0
315 p(2,1) = 0.0000000000000000d+00
316 p(2,2) = 0.0000000000000000d+00
317 p(2,3) = 100.0000000000000d0
319 p(3,0) = 72.02905523601244d0
320 p(3,1) = 32.33193950844802d0
321 p(3,2) = 64.10893732125734d0
322 p(3,3) = 5.733641195059672d0
324 p(4,0) = 25.07519219892066d0
325 p(4,1) = 2.906699172283076d0
326 p(4,2) = -24.65281899347465d0
327 p(4,3) = 3.543286526606821d0
329 p(5,0) = 60.19609558319355d0
330 p(5,1) = 2.282219923475870d0
331 p(5,2) = -59.72688022443283d0
332 p(5,3) = -7.145710205299344d0
334 p(6,0) = 42.69965698187335d0
335 p(6,1) = -37.52085860420697d0
336 p(6,2) = 20.27076189665014d0
337 p(6,3) = -2.131217516367149d0
338 else if (ievent.eq.200)
then
340 p(1,0) = -100.0000000000000d0
341 p(1,1) = 0.0000000000000000d+00
342 p(1,2) = 0.0000000000000000d+00
343 p(1,3) = -100.0000000000000d0
345 p(2,0) = -100.0000000000000d0
346 p(2,1) = 0.0000000000000000d+00
347 p(2,2) = 0.0000000000000000d+00
348 p(2,3) = 100.0000000000000d0
350 p(3,0) = 69.63925440091566d0
351 p(3,1) = -68.33683461638704d0
352 p(3,2) = -3.803463998624776d0
353 p(3,3) = -12.85443307731553d0
355 p(4,0) = 34.49634785978589d0
356 p(4,1) = 21.54290146454439d0
357 p(4,2) = -15.60187514686951d0
358 p(4,3) = 21.96549348532329d0
360 p(5,0) = 66.20628864972379d0
361 p(5,1) = 65.60907499667245d0
362 p(5,2) = -3.500460884404577d0
363 p(5,3) = -8.152834381263276d0
365 p(6,0) = 29.65810908957467d0
366 p(6,1) = -18.81514184482980d0
367 p(6,2) = 22.90580002989887d0
368 p(6,3) = -0.9582260267444891d0
370 else if (ievent.eq.250)
then
372 p( 1, 0) = -159.2798637127597488d0
373 p( 1, 1) = 0.0000000000000000d0
374 p( 1, 2) = 0.0000000000000000d0
375 p( 1, 3) = -159.2798637127597488d0
376 p( 2, 0) = -237.6261878491357891d0
377 p( 2, 1) = 0.0000000000000000d0
378 p( 2, 2) = 0.0000000000000000d0
379 p( 2, 3) = 237.6261878491357891d0
380 p( 3, 0) = 133.6018149124976446d0
381 p( 3, 1) = 60.5134748218625731d0
382 p( 3, 2) = 112.6312302087531094d0
383 p( 3, 3) = -38.7526811273741174d0
384 p( 4, 0) = 69.9256648527542950d0
385 p( 4, 1) = 66.6463990867762135d0
386 p( 4, 2) = 7.6538437283820429d0
387 p( 4, 3) = -19.7300473909055754d0
388 p( 5, 0) = 85.7314709474216841d0
389 p( 5, 1) = -38.6625297933730749d0
390 p( 5, 2) = -63.4432116127613241d0
391 p( 5, 3) = -42.7791164126088006d0
392 p( 6, 0) = 107.6471008492219283d0
393 p( 6, 1) = -88.4973441152657188d0
394 p( 6, 2) = -56.8418623243738281d0
395 p( 6, 3) = 22.9155207945124531d0
400 if (newmom.eq.1)
then
404 s(modulo(k+1,n),k) = p(k+1,0)**2
406 s(modulo(k+1,n),k) = s(modulo(k+1,n),k) - p(k+1,i)**2
408 if(abs(s(modulo(k+1,n),k)).lt.1d-14* abs(p(k+1,0))**2) s(modulo(k+1,n),k) = 0d0
410 s(k,modulo(k+1,n)) = s(modulo(k+1,n),k)
414 s(modulo(k+2,n),k) = (p(modulo(k+1,n)+1,0) + p(k+1,0))**2
416 s(modulo(k+2,n),k) = s(modulo(k+2,n),k) &
417 - (p(modulo(k+1,n)+1,i) + p(k+1,i))**2
419 s(k,modulo(k+2,n)) = s(modulo(k+2,n),k)
422 s(modulo(k+3,n),k) = (p(modulo(k+2,n)+1,0) + p(modulo(k+1,n)+1,0) + p(k+1,0))**2
424 s(modulo(k+3,n),k) = s(modulo(k+3,n),k) &
425 - (p(modulo(k+2,n)+1,i) + p(modulo(k+1,n)+1,i) + p(k+1,i))**2
427 s(k,modulo(k+3,n)) = s(modulo(k+3,n),k)
431 mominv6(k) = s(modulo(k,n),k-1)
434 mominv6(k+n) = s(modulo(k+1,n),k-1)
437 mominv6(k+2*n) = s(modulo(k+2,n),k-1)
443 fact = 4*fact*(1-fact)
446 mominv6(k) = s(modulo(k,n),k-1)*fact
449 mominv6(k+n) = s(modulo(k+1,n),k-1)*fact
452 mominv6(k+2*n) = s(modulo(k+2,n),k-1)*fact
470 double complex :: masses(0:n-1),
submasses(0:n-2)
496 integer,
intent(in) :: n,k
497 integer :: i,cnt,moms,limit
498 double complex :: mominv(n*(n-1)/2),
submominv((n-1)*(n-2)/2)
506 if ((k.ge.i).and.(i.ge.k-moms+1).or.(n+k.ge.i).and.(i.ge.n+k-moms+1))
then
517 if (mod(n,2).eq.1)
then
519 if(k.le.(n-1)/2)
then
526 if ((k.ge.i).and.(i.ge.k-moms+1).or.(n+k.ge.i).and.(i.ge.n+k-moms+1))
then
527 submominv(cnt) = mominv(n*moms+i-(n-1)/2)
540 if ((k.ge.i).and.(i.ge.k-moms+1).or.(n+k.ge.i).and.(i.ge.n+k-moms+1))
then