!********************************************************
!*       Program to demonstrate arcsine recursion       *
!* ---------------------------------------------------- *
!*   Reference: BASIC Scientific Subroutines, Vol. II   *
!*   By F.R. Ruckdeschel, BYTE/McGRAWW-HILL, 1981 [1].  *
!*                                                      *
!*              F90 Version by J.-P. Moreau, Paris.     *
!*                       (www.jpmoreau.fr)              *
!* ---------------------------------------------------- *
!* SAMPLE RUN:                                          *
!*                                                      *
!*   X       ARCSIN(X)       STEPS        ERROR         *
!*  -----------------------------------------------     *
!*  0.00     0.0000000         0       0.0000000000     *
!*  0.05     0.0500209         7      -0.0000000013     *
!*  0.10     0.1001674         8      -0.0000000025     *
!*  0.15     0.1505683         9      -0.0000000021     *
!*  0.20     0.2013579        10      -0.0000000013     *
!*  0.25     0.2526803        10      -0.0000000025     *
!*  0.30     0.3046927        11      -0.0000000011     *
!*  0.35     0.3575711        11      -0.0000000017     *
!*  0.40     0.4115168        11      -0.0000000025     *
!*  0.45     0.4667653        12      -0.0000000009     *  
!*  0.50     0.5235988        12      -0.0000000012     *
!*  0.55     0.5823642        12      -0.0000000016     *
!*  0.60     0.6435011        12      -0.0000000021     *
!*  0.65     0.7075844        13      -0.0000000007     *
!*  0.70     0.7753975        13      -0.0000000008     *
!*  0.75     0.8480621        13      -0.0000000010     *
!*  0.80     0.9272952        13      -0.0000000012     *
!*  0.85     1.0159853        13      -0.0000000014     *
!*  0.90     1.1197695        14      -0.0000000004     *
!*  0.95     1.2532359        14      -0.0000000004     *
!*  1.00     1.5707963         0       0.0000000000     *
!*                                                      *
!********************************************************
PROGRAM Arcsin

real*8   e,x,y
integer  i, m

  e=1.d-8
  write(*,*) ' X       ARCSIN(X)       STEPS        ERROR    '
  write(*,*) '-----------------------------------------------'
  x=0.d0
  do i = 1, 21
    call Arc_sin(m,e,x,y)
    write(*,50) x,y,m, (dsin(y)-x)
    x = x + 0.05d0
  end do
  stop
50 format(' ',f4.2,'     ',f9.7,'        ',i2,'      ',f13.10)
END

!*************************************************
!*       Arcsin(x) recursion subroutine          *
!* Input is x (-1<x<1), output is y=arcsin(x),   *
!* convergence criteria is e.                    *
!* --------------------------------------------- *
!* Reference: Computational Analysis by Henrici. *
!*************************************************
Subroutine Arc_sin(m,e,x,y)  
! Labels used: 100, 200, 300
  integer m
  real*8 e,x,y
  real*8 pi,u0,u1,u2
  m=0; pi = 4.d0*datan(1.d0)
  ! Guard against failure
  if (e <= 0) return
  if (x.ne.0) goto 100
  y=0.d0
  return
  ! Check range
100 if (dabs(x) > 1.d0) goto 300
  u0=x*dsqrt(1.d0-x*x)
  u1=x
200 u2=u1*dsqrt(2.d0*u1/(u1+u0))
  y=u2
  m=m+1
  if (dabs(u2-u1)<e) return
  u0=u1; u1=u2;
  goto 200
300 if (dabs(x-1.d0) < 1e-10) y=pi/2.d0
  return
end

! End of file arcsin.f90