Drag Function Subroutines (csnx, frdrag)
SUBROUTINE csnx(x,cosx,sinx,nmax)
C
C: Iterative computation of cos(nx), sin(nx)
C
C cos(nx) = cos[(n-1)x+x] = cos(n-1)x*cosx - sin(n-1)x*sinx
C sin(nx) = sin[(n-1)x+x] = sin(n-1)x*cosx + cos(n-1)x*sinx
C
C Parameter description:
C
C x ...... variable
C cosx ... vector of cosnx
C sinx ... vector of sinnx
C nmax ... max. evaluation degree
C
C Author: Zongping Chen, GFZ/D-PAF 95/03/30
C with some modifications by Andrew Sinclair, RGO
C
parameter (n0=30)
implicit DOUBLE PRECISION (A-H,O-Z)
dimension cosx(n0),sinx(n0)
if (nmax.gt.n0) then
print *,' !! Fatal !! SUBR CSNX: parameter N0 too small',
+ ' dimensioned. stop!'
stop ' ++ Error ++ '
endif
cosx(1)=cos(x)
sinx(1)=sin(x)
if (nmax.eq.1) goto 11
do 10 ii=2,nmax
cosx(ii)=cosx(ii-1)*cosx(1)-sinx(ii-1)*sinx(1)
sinx(ii)=sinx(ii-1)*cosx(1)+cosx(ii-1)*sinx(1)
10 continue
11 continue
return
end
Subroutine frdrag(t,x1,x2,x3,nmax,v)
C
C: Evaluation of the drag function
C
C Drag Function:
C v = x1 + x2*sum[(-1)**k coskX/k**2] + x3*sum[(-1)**k sinkX/k]
C where X = 2*pi*(t-0.5)
C
C Parameter description:
C
C t ........... time in days from epoch
C x1,x2,x3 .... DRAG FRCO coefficients
C nmax ........ degree of the Fourier expansion
C v ........... Value after the Fourier expansion
C
C Author: Zongping Chen, GFZ/D-PAF 95/03/30
C with some modifications by Andrew Sinclair, RGO
C
parameter (n0=30)
implicit DOUBLE PRECISION (A-H,O-Z)
dimension cosx(n0),sinx(n0)
if (nmax.gt.n0) then
print *,' !! Fatal !! SUBR FRDRAG: parameter N0 too small',
+ ' dimensioned. stop!'
stop ' ++ Error ++ '
endif
pi=4.d0*atan(1.d0)
gx=2.d0*pi*(t-0.5d0)
call csnx (gx,cosx,sinx,nmax)
k=-1
sumr2i=0.d0
sumr3i=0.d0
do 20 ii=1,nmax
cosxii=cosx(ii)*dble(k)/dble(ii*ii)
sinxii=sinx(ii)*dble(k)/dble(ii)
sumr2i=sumr2i+cosxii
sumr3i=sumr3i+sinxii
k=-k
20 continue
v=x1+x2*sumr2i+x3*sumr3i
return
end
