FUNCTION mysin(x)
  IMPLICIT NONE
  DOUBLE PRECISION, INTENT (in) :: x
  DOUBLE PRECISION :: mysin

  DOUBLE PRECISION :: sum, term
  INTEGER :: k

  sum = 0.0d0
  term = x
  DO k = 1,11,2
     sum = sum + term
     term = term * (- x * x / (k+1) / (k+2))
  END DO
  mysin = sum
END FUNCTION mysin

PROGRAM sin35
  IMPLICIT NONE

  DOUBLE PRECISION :: x
  DOUBLE PRECISION, PARAMETER :: pi = 3.141592653589793d0
  DOUBLE PRECISION, PARAMETER :: deg2rad = pi / 180.0d0

  INTERFACE
     FUNCTION mysin(x)
       IMPLICIT NONE
       DOUBLE PRECISION, INTENT (in) :: x
       DOUBLE PRECISION :: mysin
     END FUNCTION mysin
  END INTERFACE

  x = 35.0d0 * deg2rad

  PRINT *, mysin(x), sin(x)

END PROGRAM sin35
