1 subroutine sp0001(gout) 2c $Id$ 3c ***** special fast routine for -p- loop for 0001 **** 4 implicit none 5 double precision gout(*) 6c 7 double precision a0, b0, c0, a1, b1, c1, a2 8 common/tabint/a0(333),b0(333),c0(334),a1(333),b1(333),c1(334), 9 + a2(4000) 10c 11#include "auxvar.fh" 12#include "miscg.fh" 13#include "ginf.fh" 14#include "pgeom.fh" 15#include "shllfo.fh" 16#include "spgeom.fh" 17#include "qgeom.fh" 18#include "maxc.fh" 19c 20 double precision qq, theta, cq, dq 21 integer n 22 integer k, l, i 23 double precision sixty, tenm12, gout1, gout2, gout3, gout4 24 double precision q0, q1, ecd, qqq, qqtest 25 double precision dq00, dq01, dq10, dq11, v, aqx, aqz 26 double precision qperp, qperp2, cosp, sinp, p, q, pqab 27 double precision h0000, h0001, h0003, g0001, g0002, g0003 28 double precision theta2, theta3, theta4, u, g, t1, t2, t3 29c 30 data sixty,tenm12/60.0d0,1.0d-12/ 31c 32 gout1 = 0.0d0 33 gout2 = 0.0d0 34 gout3 = 0.0d0 35 gout4 = 0.0d0 36 do 940 k = 1,ngc 37 gc = cg(k) 38 do 940 l = 1,ngd 39 gd = dg(l) 40 gcd = gc+gd 41 ecd = 1.0d0/gcd 42 cq = gd*ecd*rcd 43 dq = cq-rcd 44 qqq = cq*dq*gcd 45 if (qqq+sixty) 480,500,500 46 480 v = 0.0d0 47 go to 520 48 500 v = dexp(qqq)*ecd 49 520 qqtest = cmaxc(k)*cmaxd(l)*v 50 if (qqtest-error1) 560,560,540 51 540 ismlq = 0 52 go to 600 53 560 if (qqtest-error2) 940,940,580 54 580 ismlq = 1 55 600 sc = csc(k) 56 sd = csd(l) 57 pc = cpc(k) 58 pd = cpd(l) 59 dq00 = sc*sd*v 60 dq01 = sc*pd*v 61 dq10 = pc*sd*v 62 dq11 = pc*pd*v 63 aqx = acx+sing*cq 64 aqz = acz+cosg*cq 65 qperp2 = aqx*aqx+acy2 66 qperp = dsqrt(qperp2) 67 if (qperp-tenm12) 640,640,620 68 620 cosp = -aqx/qperp 69 sinp = -acy/qperp 70 go to 660 71 640 cosp = 1.0d0 72 sinp = 0.0d0 73 660 h0000 = 0.d0 74 h0001 = 0.d0 75 h0003 = 0.d0 76 do 180 i = 1,ngangb 77 isml = ismlq+ismlp(i) 78 if (isml .ge. 2) go to 180 79 auxvar = var(isml+1) 80 pqab = aqz-app(i) 81 g = 1.d0/(ep(i)+ecd) 82 p = (pqab*pqab+qperp2)*g 83 if (p .le. auxvar) go to 140 84 q0 = dp00p(i)*dsqrt(0.7853981625d0/(p*(gp(i)+gcd))) 85 q1 = 0.5d0*q0/p 86 go to 160 87 140 q = dp00p(i)/dsqrt(gp(i)+gcd) 88 qq = p*12.5d0 89 n = idint(qq) 90 theta = qq- dble(n) 91 theta2 = theta*(theta-1.d0) 92 theta3 = theta2*(theta-2.d0) 93 theta4 = theta2*(theta+1.d0) 94 q0 = (a0(n+1)+theta*b0(n+1)-theta3*c0(n+1)+theta4*c0(n+2))*q 95 q1 = (a1(n+1)+theta*b1(n+1)-theta3*c1(n+1)+theta4*c1(n+2))*q 96 160 u = g*q1 97 h0000 = h0000+q0 98 h0001 = h0001+u 99 h0003 = h0003-u*pqab 100 180 continue 101 h0001 = h0001*ecd*qperp 102 h0003 = h0003*ecd 103 p = dq*h0000 104 g0001 = h0001*cosp+p*sing 105 g0002 = h0001*sinp 106 g0003 = h0003+p*cosg 107 gout1 = gout1+dq00*h0000 108 gout2 = gout2+dq01*g0001 109 gout3 = gout3+dq01*g0002 110 gout4 = gout4+dq01*g0003 111 940 continue 112 gout(1) = gout1 113 t1 = gout2 114 t2 = gout3 115 t3 = gout4 116 gout(2) = p11*t1+p21*t2+p31*t3 117 gout(3) = p12*t1+p22*t2+p32*t3 118 gout(4) = p13*t1+p23*t2+p33*t3 119 return 120 end 121