JHUGen MELA  JHUGen v7.5.6, MELA v2.4.2
Matrix element calculations as used in JHUGen.
Functions/Subroutines
reductionab Module Reference

Functions/Subroutines

subroutine calca (A, Auv, m02, rmax, Aerr)
 
subroutine calcb (B, Buv, p10, m02, m12, rmax, id, Berr)
 
subroutine calcbred (B, Buv, p10, m02, m12, rmax, id, Berr)
 
subroutine calcdb (DB, DBuv, p10, m02, m12, rmax, id, DBerr)
 

Function/Subroutine Documentation

◆ calca()

subroutine reductionab::calca ( double complex, dimension(0:rmax/2), intent(out)  A,
double complex, dimension(0:rmax/2), intent(out)  Auv,
double complex, intent(in)  m02,
integer, intent(in)  rmax,
double precision, dimension(0:rmax), intent(out), optional  Aerr 
)

Definition at line 53 of file reductionAB.F90.

53 
54  integer, intent(in) :: rmax
55  double complex, intent(in) :: m02
56  double complex, intent(out) :: A(0:rmax/2), Auv(0:rmax/2)
57  double precision, optional, intent(out) :: Aerr(0:rmax)
58  double complex :: A0_coli, mm02, elimminf2_coli
59  double complex :: fac, ccc
60  integer :: n0
61 
62  if(present(aerr)) then
63  aerr = 0d0
64  end if
65  auv(0) = m02
66  a(0) = a0_coli(m02)
67 
68  ccc = 2d0
69  mm02 = elimminf2_coli(m02)
70  do n0=1,rmax/2
71  fac = m02/(2*(n0+1))
72  ccc = ccc*fac
73  auv(n0) = fac*auv(n0-1)
74  a(n0) = fac*(a(n0-1)+ccc)
75  end do
76  if(present(aerr)) then
77  aerr = abs(a(0))*acc_def_b
78  end if
79 

◆ calcb()

subroutine reductionab::calcb ( double complex, dimension(0:rmax,0:rmax), intent(inout)  B,
double complex, dimension(0:rmax,0:rmax), intent(inout)  Buv,
double complex, intent(in)  p10,
double complex, intent(in)  m02,
double complex, intent(in)  m12,
integer, intent(in)  rmax,
integer, intent(in)  id,
double precision, dimension(0:rmax), intent(out), optional  Berr 
)

Definition at line 92 of file reductionAB.F90.

92 
93  integer, intent(in) :: rmax, id
94  double complex, intent(in) :: p10,m02,m12
95  double complex, intent(inout) :: B(0:rmax,0:rmax), Buv(0:rmax,0:rmax)
96  double precision, optional, intent(out) :: Berr(0:rmax)
97  double complex, allocatable :: Baux(:,:), Buvaux(:,:),fct(:)
98  double precision, allocatable :: Berraux(:)
99  double complex :: x(3)
100  integer :: rank,switch,cnt,n0,n1,r
101  logical :: nocalc,wrica
102 
103  if (use_cache_system) then
104  if ((ncache.gt.0).and.(ncache.le.ncache_max)) then
105 ! if (use_cache(ncache).ge.2) then
106  x(1)=p10
107  x(2)=m02
108  x(3)=m12
109  rank = rmax
110  switch = 0
111 
112  allocate(fct(2*ncoefsg(rmax,2)+rmax+1))
113  call readcache(fct,2*ncoefsg(rmax,2)+rmax+1,x,3,1,id,2,rank,nocalc,wrica)
114 
115  if(nocalc)then
116  cnt = 0
117  do r=0,rmax
118  do n0=0,r
119  n1 = r-n0
120 
121  cnt = cnt+1
122  b(n0,n1) = fct(cnt)
123  cnt = cnt+1
124  buv(n0,n1) = fct(cnt)
125 
126  end do
127  cnt = cnt+1
128  if(present(berr))then
129  berr(r) = real(fct(cnt))
130  end if
131  end do
132  return
133  endif
134 
135 
136  if(rank.eq.rmax) then
137 
138  allocate(berraux(0:rank))
139  call calcbred(b,buv,p10,m02,m12,rank,id,berraux)
140  if(present(berr)) berr=berraux
141 
142  if (wrica) then
143  cnt = 0
144  do r=0,rank
145  do n0=0,r
146  n1 = r-n0
147 
148  cnt = cnt+1
149  fct(cnt) = b(n0,n1)
150  cnt = cnt+1
151  fct(cnt) = buv(n0,n1)
152 
153  end do
154  cnt = cnt+1
155  fct(cnt) = berraux(r)
156  end do
157 
158  call writecache(fct,2*ncoefsg(rank,2)+rank+1,id,2,rank)
159 
160  end if
161 
162  return
163 
164 
165  else
166  allocate(baux(0:rank,0:rank))
167  allocate(buvaux(0:rank,0:rank))
168  allocate(berraux(0:rank))
169 
170  call calcbred(baux,buvaux,p10,m02,m12,rank,id,berraux)
171 
172  if (wrica) then
173  cnt = 0
174  deallocate(fct)
175  allocate(fct(2*ncoefsg(rank,2)+rank+1))
176  do r=0,rank
177  do n0=0,r
178  n1 = r-n0
179 
180  cnt = cnt+1
181  fct(cnt) = baux(n0,n1)
182  cnt = cnt+1
183  fct(cnt) = buvaux(n0,n1)
184 
185  end do
186  cnt = cnt+1
187  fct(cnt) = berraux(r)
188  end do
189 
190  call writecache(fct,2*ncoefsg(rank,2)+rank+1,id,2,rank)
191 
192  end if
193 
194  b(0:rmax,0:rmax) = baux(0:rmax,0:rmax)
195  buv(0:rmax,0:rmax) = buvaux(0:rmax,0:rmax)
196  if(present(berr))then
197  berr = berraux(0:rmax)
198  end if
199  return
200 
201  end if
202 ! end if
203  end if
204  end if
205 
206 
207  if(present(berr)) then
208  call calcbred(b,buv,p10,m02,m12,rmax,id,berr)
209  else
210  allocate(berraux(0:rmax))
211  call calcbred(b,buv,p10,m02,m12,rmax,id,berraux)
212  deallocate(berraux)
213  end if
214 

◆ calcbred()

subroutine reductionab::calcbred ( double complex, dimension(0:rmax,0:rmax), intent(inout)  B,
double complex, dimension(0:rmax,0:rmax), intent(inout)  Buv,
double complex, intent(in)  p10,
double complex, intent(in)  m02,
double complex, intent(in)  m12,
integer, intent(in)  rmax,
integer, intent(in)  id,
double precision, dimension(0:rmax), intent(out)  Berr 
)

Definition at line 227 of file reductionAB.F90.

227 
228  integer, intent(in) :: rmax, id
229  double complex, intent(in) :: p10,m02,m12
230  double complex, intent(inout) :: B(0:rmax,0:rmax), Buv(0:rmax,0:rmax)
231  double precision, intent(out) :: Berr(0:rmax)
232 ! double complex, allocatable :: A(:), Auv(:)
233  double complex :: A(0:rmax-1), Auv(0:rmax-1)
234  double complex :: q2, mm02, mm12, f1
235  double complex :: Bn_coli,elimminf2_coli
236  integer :: rmaxA,n0,n1,r,sgn,k,rid,r0
237  logical :: finarg
238  double complex, parameter :: cd0 = dcmplx(0d0,0d0)
239 
240  r0 = 0
241 
242  ! r=0
243  if (r0.eq.0) then
244  buv(0,0) = 1d0
245  b(0,0) = bn_coli(0,p10,m02,m12)
246  end if
247 
248  ! calculate B(0,n1)
249  sgn = (-1)**r0
250  do n1=max(r0,1),rmax
251  sgn = -sgn
252  buv(0,n1) = sgn/(n1+1d0)
253  b(0,n1) = bn_coli(n1,p10,m02,m12)
254  end do
255 
256 
257  ! calculate B(n0,n1) for n0=/=0
258  if (rmax.ge.1) then
259 
260  call elminf2iv_coli(p10,m02,m12,q2,mm02,mm12,finarg)
261 
262 ! q2 = elimminf2_coli(p10)
263 ! mm12 = elimminf2_coli(m12)
264 ! mm02 = elimminf2_coli(m02)
265 ! finarg=q2.ne.cd0.or.mm12.ne.cd0.or.mm02.ne.cd0
266 
267  ! at least one of the invariants p10,m02,m12 different from zero
268  if (finarg) then
269 
270  ! calculation of A-functions
271  rmaxa = max(rmax-1,0)
272 ! allocate(A(0:rmaxA))
273 ! allocate(Auv(0:rmaxA))
274  call calca(a,auv,mm12,2*rmaxa)
275 
276  ! p10 non-negligible
277  if (abs(q2)/(abs(q2)+abs(mm02+mm12)).gt.dprec_coli) then
278 
279  f1 = q2-mm12+mm02
280 
281 ! do n0=1,rmax/2
282 ! k=rmax-2*n0
283  do r=2,2*rmax
284  do n0=max(1,r-rmax),r/2
285  k=r-2*n0
286  do n1=max(r0-2*n0,0),k
287  buv(n0,n1) = ((-1)**n1*auv(n0-1) + 2*mm02*buv(n0-1,n1) &
288  + f1*buv(n0-1,n1+1))/(2*n0+n1+1)/2
289  b(n0,n1) = ((-1)**n1*a(n0-1) + 2*mm02*b(n0-1,n1) &
290  + f1*b(n0-1,n1+1) + 4*buv(n0,n1))/(2*n0+n1+1)/2
291 ! write(*,*) 'CalcBred B(n0,n1)',n0,n1,B(n0,n1)
292  end do
293  end do
294  end do
295 
296  ! p10 negligible
297  else
298 
299 ! do n0=1,rmax/2
300 ! k=rmax-2*n0
301  do r=2,2*rmax
302  do n0=max(1,r-rmax),r/2
303  k=r-2*n0
304  do n1=max(r0-2*n0,0),k
305 
306  buv(n0,n1) = ((-1)**n1*auv(n0-1) + mm02*buv(n0-1,n1))/(n0+n1+1)/2
307  b(n0,n1) = ((-1)**n1*a(n0-1) + mm02*b(n0-1,n1) &
308  + 2*buv(n0,n1))/(n0+n1+1)/2
309  end do
310  end do
311  end do
312 
313  end if
314 
315  ! p10=m02=m12=0
316  else
317 
318  buv(1:,:) = 0d0
319  b(1:,:) = 0d0
320 
321  end if
322  end if
323 
324  berr(0) = abs(b(0,0))
325  do r=1,rmax
326  berr(r) = max(berr(r-1),abs(b(0,r)))
327  end do
328  berr = berr*acc_def_b
329 

◆ calcdb()

subroutine reductionab::calcdb ( double complex, dimension(0:rmax/2,0:rmax), intent(inout)  DB,
double complex, dimension(0:rmax/2,0:rmax), intent(inout)  DBuv,
double complex, intent(in)  p10,
double complex, intent(in)  m02,
double complex, intent(in)  m12,
integer, intent(in)  rmax,
integer, intent(in)  id,
double precision, dimension(0:rmax), intent(out), optional  DBerr 
)

Definition at line 341 of file reductionAB.F90.

341 
342  integer, intent(in) :: rmax, id
343  double complex, intent(in) :: p10,m02,m12
344  double complex, intent(inout) :: DB(0:rmax/2,0:rmax), DBuv(0:rmax/2,0:rmax)
345  double precision, optional, intent(out) :: DBerr(0:rmax)
346  double precision :: Berr(0:rmax)
347  double complex :: B(0:rmax,0:rmax), Buv(0:rmax,0:rmax)
348 ! double complex, allocatable :: A(:), Auv(:)
349  double complex :: A(0:rmax-1), Auv(0:rmax-1)
350  double complex :: q2, mm02, mm12, f1
351  double complex :: DBn_coli,elimminf2_coli
352  integer :: rmaxA,n0,n1,r,sgn,k,rid
353  logical :: finarg
354  double complex, parameter :: cd0 = dcmplx(0d0,0d0)
355 
356 ! write(*,*) 'CalcDB in',p10,m02,m12,rmax,id
357 
358  ! r=0
359  dbuv(0,0) = 0d0
360  db(0,0) = dbn_coli(0,p10,m02,m12)
361 
362  call calcbred(b,buv,p10,m02,m12,rmax,id,berr)
363 
364  ! calculate DB(0,n1)
365  do n1=1,rmax
366  sgn = -sgn
367  dbuv(0,n1) = 0d0
368  db(0,n1) = dbn_coli(n1,p10,m02,m12)
369  end do
370 
371 
372  ! calculate DB(n0,n1) for n0=/=0
373  if (rmax.ge.2) then
374 
375  call elminf2iv_coli(p10,m02,m12,q2,mm02,mm12,finarg)
376 
377 ! q2 = elimminf2_coli(p10)
378 ! mm12 = elimminf2_coli(m12)
379 ! mm02 = elimminf2_coli(m02)
380 ! finarg=q2.ne.cd0.or.mm12.ne.cd0.or.mm02.ne.0d0
381 
382  ! at least one of the invariants p10,m02,m12 different from zero
383  if (q2.ne.0d0.or.mm02.ne.0d0.or.mm12.ne.0d0) then
384 
385  ! calculation of A-functions
386  rmaxa = max(rmax-1,0)
387 ! allocate(A(0:rmaxA))
388 ! allocate(Auv(0:rmaxA))
389  call calca(a,auv,mm12,2*rmaxa)
390 
391  f1 = q2-mm12+mm02
392 
393  do r=2,rmax
394  do n0=max(1,r-rmax),r/2
395  k=r-2*n0
396  do n1=0,k
397  dbuv(n0,n1) = (buv(n0-1,n1+1) + 2*mm02*dbuv(n0-1,n1) &
398  + f1*dbuv(n0-1,n1+1))/(2*n0+n1+1)/2
399 
400  db(n0,n1) = b(n0-1,n1+1) + 4*dbuv(n0,n1)
401  if(mm02.ne.0d0) then
402  db(n0,n1) = db(n0,n1) + 2*mm02*db(n0-1,n1)
403  end if
404  if(f1.ne.0d0) then
405  db(n0,n1) = db(n0,n1) + f1*db(n0-1,n1+1)
406  end if
407 
408  db(n0,n1) = db(n0,n1) /(2*(2*n0+n1+1))
409 
410  end do
411  end do
412  end do
413 
414  ! p10=m02=m12=0
415  else
416 
417  dbuv(1:,:) = 0d0
418 
419  do r=2,rmax
420  do n0=max(1,r-rmax),r/2
421  k=r-2*n0
422  do n1=0,k
423  if (.not.ir_rational_terms_cll) then
424  dbuv(n0,n1) = buv(n0-1,n1+1)/(2*n0+n1+1)/2
425  end if
426  db(n0,n1) = (b(n0-1,n1+1) + 4*dbuv(n0,n1))/(2*n0+n1+1)/2
427  end do
428  end do
429  end do
430 
431  end if
432  end if
433 
434  dberr(0) = abs(db(0,0))
435  do r=1,rmax
436  dberr(r) = max(dberr(r-1),abs(db(0,r)))
437  end do
438  dberr = dberr*acc_def_b
439