c===========================================================
c     Program illustrating "catastrophic" loss of precision
c     resulting from the subtraction of two nearly equal 
c     floating point values.
c===========================================================
      program         catprec

      implicit        none

      real*8          x  
      parameter     ( x = 0.2d0 )

      integer         i
      real*8          h,     dsinx

      write(*,*) '      h      d(sin) approx   '//
     &           'd(sin) exact    d(sin) err'
      write(*,*)

      h = 0.5d0
      do i = 1 , 16
c-----------------------------------------------------------
c        Algebraically, in the limit h -> 0, dsinx should
c        approach cos(x), but sin(x+h) -> sin(x) so 
c        catastrophic loss of precision occurs.
c-----------------------------------------------------------
         dsinx = (sin(x+h) - sin(x)) / h
         write(*,1000) h,  dsinx,  cos(x),   dsinx - cos(x)
1000     format(1P,E12.3,2E16.8,E12.3)
         h = 0.125d0 * h
      end do

      stop

      end