| Fortran |
c
c######################################################################
c
subroutine alcheb(pet, pcf, kncf, pfit, kermsg)
c
c this is an ornl:cfadc subroutine to calculate cross sections in
c (cm[2]) versus energy in (ev/amu) or rate coefficients in
c (cm[3]/s) versus maxwellian temperature in (ev) from chebyshev
c polynomial fitting coefficients
c
c these fits are valid only between the limits emin and emax,
c which are coefficients pcf(10) and pcf(11) in the entry data field
c
c pet = collision energy in ev/amu or maxwellian temperature in ev
c
c kermsg = blank if no errors
c
c pfit = cross section in cm[2] or rate coefficient in cm[3]/s
c
c written by h. hunter, cfadc oak ridge national laboratory
c (modified to aladdin calling structure 4/21/88 r.a. hulse)
c
c------------------------------------------------------------------------
c
double precision pet, pcf, pfit
double precision emin, emax, cheb, eminl, emaxl, enl, xnorm
double precision twox, prev, prev2
dimension pcf(11)
character*(*) kermsg
emin = pcf(10)
emax = pcf(11)
if(pet .ge. emin .and. pet .le. emax) then
kermsg = ' '
else
kermsg = 'outside range of fit in alcheb'
return
endif
c
c calculate polynomial using recursion relation
c
k = 9
cheb = pcf(k)
eminl = dlog(emin)
emaxl = dlog(emax)
enl= dlog(pet)
k = k-1
xnorm = (enl-eminl-(emaxl-enl)) / (emaxl-eminl)
twox = 2.0d0 * xnorm
prev2 = 0.0d+00
10 prev = cheb
if(k .ne. 1) then
cheb = pcf(k) + twox*prev - prev2
prev2 = prev
k = k-1
go to 10
endif
cheb = 0.5d0*pcf(1) + xnorm*prev - prev2
pfit = dexp(cheb)
100 return
c
end |