cfuncs.f

Go to the documentation of this file.

      SUBROUTINE cfuncs(Z,C)
*
*
*       C-functions for two-body iteration.
*       -----------------------------------
*
*       Coded by Seppo Mikkola (11/99).
*       Reference: New Astron. 3, 309, 1998.
*       ------------------------------------
*
      IMPLICIT REAL*8 (a-h,o-z)
      parameter(c6=1.d0/6.d0,c132=1.d0/132.d0,c56=1.d0/56.d0,
     &   c30=1.d0/30d0,c24=1.d0/24.d0,c156=1.d0/156.d0,
     &   c90=1.d0/90.d0,c110=1.d0/110.d0,c16=1.d0/16.d0,c8=1.d0/8.d0,
     &   c72=1.d0/72.d0,c42=1.d0/42.d0,c120=1.d0/120.d0,u=1.d0)
      REAL*8 c(5)
      data large/0/
*
*
*       Reduce the argument by 4 until < 0.07.
      h=z
      DO 1 k=0,20
          kk = k
          IF (abs(h).lt. 0.07) go to 5
          h=0.25d0*h
    1 CONTINUE
*
      if (large.lt.100) then
          large=large+1
          WRITE (6,3) large, k, z, h
    3     FORMAT (' WARNING!    CFUNCS    # k Z h   ',2i4,1p,2e10.2)
      end if
*
*       Form the two highest c-functions.
    5 c(4)=(u-h*(u-h*(u-h*c90/(u+h*c132))*c56)*c30)*c24
      c(5)=(u-h*(u-h*(u-h*c110/(u+h*c156))*c72)*c42)*c120
*
*       Obtain the lower c-function by recursion.
      DO 4 i=1,kk
          c3=c6-h*c(5)
          c2=.5d0-h*c(4)
          c(5)=(c(5)+c(4)+c2*c3)*c16
          c(4)=c3*(2.d0-h*c3)*c8
          h=4.d0*h
    4 CONTINUE
*
      c(3)=c6-z*c(5)
      c(2)=.5d0-z*c(4)
      c(1)=1.d0-z*c(3)
*
      RETURN
*
      END