c_______________________________________________________________________________
c
c     Este programa utiliza la subrutina MNewton() para calcular la raiz
c     de la funcion
c
c               f(x) = cos(x) - x
c
c     Los valores de la función f(x) y de su derivada f'(x) se
c     proporcionan mediante una subrutina definida por el usuario
c
c               SP1(x,fx,dfx,fxdfx),
c
c     que para cada valor de x calcula
c
c               fx = f(x), dfx = f'(x), fxdfx = f(x)/f'(x).                                                                                                                                                                
c 
c_______________________________________________________________________________

      implicit real*8 (a-h,o-z)

      external SP1             ! nombre de la subrutina que calcula f(x) y f'(x)

	write(6,'( a,$)') ' Valor Inicial  -> x      = '
      read(5,*) x
	write(6,'( a,$)') ' Err. Relativo  -> r      = '
      read(5,*) r
	write(6,'( a,$)') ' Err. Absoluto  -> e      = '
      read(5,*) e
	write(6,'( a,$)') ' No. Max. Iter. -> mxiter = '
      read(5,*) mxiter
      write(6,'()')

      call MNewton(SP1,x,fx,dfx,r,e,mxiter,niter,.TRUE.)

      if (niter.le.mxiter) then
        write(6,100) niter,x,fx,dfx
  100   format(///' SOLUCION: '//
     .       ' niter = ',i6,' -> x = ',d23.16,'; f (x) = ',d23.16,/
     .       '         ',6x,'        ',23x   ,'  f''(x) = ',d23.16) ! ''
      else
        write(6,'(a)') ' Error: el metodo no ha convergido.'
      endif

      read(5,'()')

      end


c____________________________________________________________________________SP1
c
c     Subroutine SP1(x,fx,dfx,fxdfx)
c
c     DESCRIPCION:
c
c               Siendo
c                         f(x) = cos(x) - x,
c
c               para cada valor de x esta subrutina calcula 
c
c                         fx = f(x), dfx = f'(x), fxdfx = f(x)/f'(x).
c
c_______________________________________________________________________________
                                                                                                                                                         
      Subroutine SP1(x,fx,dfx,fxdfx)
      implicit real*8 (a-h,o-z)

      fx    =  cos(x)-x
      dfx   = -sin(x)-1
      fxdfx =  fx/dfx

      return
      end


c_______________________________________________________________________________
c
c     Subrutina MNewton(S,x,fx,dfx,r,e,mxiter,niter,Snoopy)
c
c     DESCRIPCION:
c
c               Esta subrutina calcula la raiz de una funcion cualquiera f(x).
c               Se utiliza el Método de Newton.
c
c               Los valores de la función f(x) y de su derivada f'(x) se
C               proporcionan mediante una subrutina definida por el usuario
c
c                         S(x,fx,dfx,dfofx),
c
c               que para cada valor de x calcula
c
c                         fx = f(x), dfx = f'(x), fxdfx = f(x)/f'(x).  
c
c     ARGUMENTOS:
c
c        S      -> entrada : nombre de la subrutina definida por el usuario                                                                                      
c        x      -> entrada : aproximacion inicial
c               -> salida  : solucion aproximada de la ecuacion f(x) = 0
c        fx     -> salida  : valor de f(x) en la solucion aproximada x
c        dfx    -> salida  : valor de f'(x) en la solucion aproximada x
c        r      -> entrada : tolerancia de la sol. aprox. en error relativo 
c        e      -> salida  : tolerancia de la sol. aprox. en error absoluto 
c        mxiter -> entrada : no. maximo de iteraciones permitidas
c        niter  -> salida  : no. de iteraciones realizadas
c                            [niter > mxiter => el proceso no ha convergido]
c        Snoopy -> entrada : indica si hay que imprimir (.TRUE.) o no (.FALSE.)
c                            el resultado de cada iteracion 
c_______________________________________________________________________________

      Subroutine MNewton(S,x,fx,dfx,r,e,mxiter,niter,Snoopy)
      implicit real*8 (a-h,o-z)
      logical Snoopy

      niter=0
      call S(x,fx,dfx,fxdfx)
      if (Snoopy) write(6,100) niter,x,fx,dfx
  100 format(' niter = ',i6,' -> x = ',d23.16,'; f (x) = ',d23.16,/
     .       '         ',6x,'        ',23x   ,'  f''(x) = ',d23.16) ! ''

      do while (niter.lt.mxiter)
        ex=abs(x)*r
        if (e.gt.ex) ex=e
        dx=-fxdfx      ! metodo de Newton (segundo orden)
        x=x+dx
        niter=niter+1
        call S(x,fx,dfx,fxdfx)
        if (Snoopy) write(6,100) niter,x,fx,dfx
        if (abs(dx).lt.ex) return
      enddo

      niter=mxiter+1   ! se indica asi que el proceso no ha convergido
      return

      end


