!   D02GBF Example Program Text
!   Mark 24 Release. NAG Copyright 2012.

    Module d02gbfe_mod

!     Data for D02GBF example program

!     .. Use Statements ..
      Use nag_library, Only: nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Real (Kind=nag_wp), Parameter        :: one = 1.0_nag_wp
      Real (Kind=nag_wp), Parameter        :: zero = 0.0_nag_wp
      Integer, Parameter                   :: iset = 1, n = 2, nin = 5, nout = 6
!     .. Local Scalars ..
      Real (Kind=nag_wp)                   :: eps
    Contains
      Subroutine fcnf(x,f)

!       .. Scalar Arguments ..
        Real (Kind=nag_wp), Intent (In)      :: x
!       .. Array Arguments ..
        Real (Kind=nag_wp), Intent (Out)     :: f(*)
!       .. Executable Statements ..
        f(1:2) = 0.0E0_nag_wp
        f(3) = 1.0E0_nag_wp
        f(4) = -1.0E0_nag_wp/eps
        Return
      End Subroutine fcnf
      Subroutine fcng(x,g)

!       .. Scalar Arguments ..
        Real (Kind=nag_wp), Intent (In)      :: x
!       .. Array Arguments ..
        Real (Kind=nag_wp), Intent (Out)     :: g(*)
!       .. Executable Statements ..
        g(1:2) = 0.0E0_nag_wp
        Return
      End Subroutine fcng
    End Module d02gbfe_mod
    Program d02gbfe

!     D02GBF Example Main Program

!     .. Use Statements ..
      Use nag_library, Only: d02gbf, nag_wp, x04abf
      Use d02gbfe_mod, Only: eps, fcnf, fcng, iset, n, nin, nout, one, zero
!     .. Implicit None Statement ..
      Implicit None
!     .. Local Scalars ..
      Real (Kind=nag_wp)                   :: a, b, tol
      Integer                              :: i, ifail, j, liw, lw, mnp, np,   &
                                              outchn
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable      :: c(:,:), d(:,:), gam(:), w(:),    &
                                              x(:), y(:,:)
      Integer, Allocatable                 :: iw(:)
!     .. Executable Statements ..
      Write (nout,*) 'D02GBF Example Program Results'
!     Skip heading in data file
      Read (nin,*)
!     mnp: maximum permitted number of mesh points.
      Read (nin,*) mnp
      liw = mnp*(2*n+1) + n
      lw = mnp*(3*n*n+5*n+2) + 3*n*n + 5*n
      Allocate (iw(liw),c(n,n),d(n,n),gam(n),w(lw),x(mnp),y(n,mnp))
!     tol: positive absolute error tolerance
!     np : determines whether a default or user-supplied mesh is used. 
!     a  : left-hand boundary point, b: right-hand boundary point.
      Read (nin,*) tol
      Read (nin,*) np
      Read (nin,*) a, b
      outchn = nout
      Call x04abf(iset,outchn)
      gam(1:n) = zero
      c(1:n,1:n) = zero
      d(1:n,1:n) = zero
      c(1,1) = one
      d(2,1) = one
      gam(2) = one
loop: Do i = 1, 2
        eps = 10.0E0_nag_wp**(-i)
        Write (nout,*)
!       ifail: behaviour on error exit
!              =1 for quiet-soft exit
!       * Set ifail to 111 to obtain monitoring information *
        ifail = 1
        Call d02gbf(a,b,n,tol,fcnf,fcng,c,d,gam,mnp,x,y,np,w,lw,iw,liw,ifail)

        If (ifail>=0) Write (nout,99999) 'Problem with epsilon = ', eps
        If (ifail==0) Then
          Write (nout,99998) np
          Write (nout,*) '       X(I)     Y(1,I)'
          Write (nout,99997)(x(j),y(1,j),j=1,np)
        Else
          Write (nout,99996) ifail
          Exit loop
        End If
      End Do loop

99999 Format (1X,A,E10.2)
99998 Format (/1X,'Approximate solution on final mesh of ',I2,' points')
99997 Format (1X,2F11.4)
99996 Format (1X/1X,' ** D02GBF returned with IFAIL = ',I5)
    End Program d02gbfe
