1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
|
subroutine splint(r,db,db2,spin,start,n,ri,dbi)
c splint calculates the cubic spline interpolation of the density
c together with subroutine spline.
c the main ideas are "stolen" from "numerical recipes", cambridge
c university press (1992).
c input r : position coordinate (array)
c input db : spin up and down density. (array)
c input db2 : second derivative of db (array)
c input spin : 1 if spin up, -1 if spin down
c input start : starting point of spline = r(start)
c input n : index, r(n) = r_max
c input ri : position of interpolation
c output dbi : interpolated value of db
implicit logical (a-z)
integer spin,spinn,n,klo,khi,k,start
double precision r,db,db2,ri,dbi,h,a,b
dimension r(561),db(2,561),db2(2,561)
spinn = spin
if (spinn.eq.-1) spinn = 2
klo = start
khi = n
10 if (khi-klo.gt.1) then
k = (khi+klo)/2
if (r(k).gt.ri) then
khi = k
else
klo = k
endif
goto 10
endif
h = r(khi)-r(klo)
a = (r(khi)-ri)/h
b = (ri-r(klo))/h
dbi = a*db(spinn,klo)+b*db(spinn,khi)+((a**3-a)*
j db2(spinn,klo)+(b**3-b)*db2(spinn,khi))*(h**2)/6.d0
return
end
|