function rjay(zetfac,gnoo) * * Function to determine value of J for each element * implicit double precision(a-h,o-z) if(gnoo.ge.-20.0 .and. gnoo.lt.0.d0)then zfac = zetfac rjnum = (0.0151/zfac)+(0.242/zfac**2.25)+(1.21/zfac**4.55) rjden = 1.d0 + 0.0371*dexp(0.906*gnoo)*(1.d0+0.928*zfac) rjay = rjnum*dexp(gnoo)/rjden elseif(gnoo .le. 10.0)then den = 1.d0+(-0.0120*gnoo)+(0.029*gnoo*gnoo)+(-0.00104*gnoo**3) zfac = zetfac/den rjnum = (0.0123/zfac)+(0.266/zfac**2.25)+(1.30/zfac**4.55) rjden = 1.d0 + 0.117*dexp(0.897*gnoo)*(1.d0+0.177*zfac) rjay = rjnum*dexp(gnoo)/rjden endif return end ************************************************************************