Program e02agfe

!     E02AGF Example Program Text

!     Mark 25 Release. NAG Copyright 2014.

!     .. Use Statements ..
      Use nag_library, Only: e02agf, e02akf, f16dnf, nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
      Real (Kind=nag_wp)               :: fiti, xmax, xmin
      Integer                          :: i, ifail, ipmax, k, kplus1, la, lda, &
                                          liwrk, lwrk, lyf, m, mf, n, np1
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable  :: a(:,:), s(:), w(:), wrk(:), x(:),    &
                                          xf(:), y(:), yf(:)
      Integer, Allocatable             :: ip(:), iwrk(:)
!     .. Intrinsic Procedures ..
      Intrinsic                        :: max, sum
!     .. Executable Statements ..
      Write (nout,*) 'E02AGF Example Program Results'

!     Skip heading in data file
      Read (nin,*)

      Read (nin,*) mf
      liwrk = 2*mf + 2
      Allocate (ip(mf),xf(mf),iwrk(liwrk))

      Read (nin,*) ip(1:mf)

!     Get max(IP) for later use

      Call f16dnf(mf,ip,1,i,ipmax)

      Read (nin,*) xf(1:mf)

      lyf = mf + sum(ip(1:mf))
      Allocate (yf(lyf))

      Read (nin,*) yf(1:lyf)

      Read (nin,*) m
      Allocate (x(m),y(m),w(m))

      Read (nin,*)(x(i),y(i),w(i),i=1,m)

      Read (nin,*) k, xmin, xmax
      kplus1 = k + 1
      n = lyf
      np1 = n + 1
      lwrk = max(4*m+3*kplus1,8*n+5*ipmax+mf+10) + 2*n + 2
      lda = kplus1
      Allocate (wrk(lwrk),a(lda,kplus1),s(kplus1))

      ifail = 0
      Call e02agf(m,kplus1,lda,xmin,xmax,x,y,w,mf,xf,yf,lyf,ip,a,s,np1,wrk, &
        lwrk,iwrk,liwrk,ifail)

      Write (nout,*)
      Write (nout,*) 'Degree  RMS residual'
      Write (nout,99999)(i,s(i+1),i=np1-1,k)
      Write (nout,*)
      Write (nout,99996) 'Details of the fit of degree ', k
      Write (nout,*)
      Write (nout,*) '  Index   Coefficient'

      Do i = 1, kplus1
        Write (nout,99997) i - 1, a(kplus1,i)
      End Do

      Write (nout,*)
      Write (nout,*) '     I      X(I)       Y(I)       Fit     Residual'

      Do i = 1, m
        la = lda*kplus1 - k

        ifail = 0
        Call e02akf(kplus1,xmin,xmax,a(kplus1,1),lda,la,x(i),fiti,ifail)

        Write (nout,99998) i, x(i), y(i), fiti, fiti - y(i)
      End Do

99999 Format (1X,I4,1P,E15.2)
99998 Format (1X,I6,3F11.4,E11.2)
99997 Format (1X,I6,F11.4)
99996 Format (1X,A,I2)
    End Program e02agfe