PROGRAM f11ddfe

!      F11DDF Example Program Text

!      Mark 23 Release. NAG Copyright 2011.

!      .. Use Statements ..
       USE nag_library, ONLY : f11bdf, f11bef, f11bff, f11ddf, f11xaf, nag_wp
!      .. Implicit None Statement ..
       IMPLICIT NONE
!      .. Parameters ..
       INTEGER, PARAMETER              :: nin = 5, nout = 6
!      .. Local Scalars ..
       REAL (KIND=nag_wp)              :: anorm, omega, sigmax, stplhs,        &
                                          stprhs, tol
       INTEGER                         :: i, ifail, irevcm, iterm, itn, la,    &
                                          liwork, lwneed, lwork, m, maxitn,    &
                                          monit, n, nnz
       CHARACTER (1)                   :: ckddf, ckxaf, norm, precon, trans,   &
                                          weight
       CHARACTER (8)                   :: method
!      .. Local Arrays ..
       REAL (KIND=nag_wp), ALLOCATABLE :: a(:), b(:), rdiag(:), wgt(:),        &
                                          work(:), x(:)
       INTEGER, ALLOCATABLE            :: icol(:), irow(:), iwork(:)
!      .. Intrinsic Functions ..
       INTRINSIC                          max
!      .. Executable Statements ..
       WRITE (nout,*) 'F11DDF Example Program Results'
       WRITE (nout,*)
!      Skip heading in data file
       READ (nin,*)

!      Read algorithmic parameters
       READ (nin,*) n, m
       READ (nin,*) nnz
       la = 3*nnz
       lwork = max(n*(m+3)+m*(m+5)+101,7*n+100,(2*n+m)*(m+2)+n+100,10*n+100)
       liwork = 2*n + 1
       ALLOCATE (a(la),b(n),rdiag(n),wgt(n),work(lwork),x(n),icol(la), &
          irow(la),iwork(liwork))
       READ (nin,*) method
       READ (nin,*) precon, norm, iterm
       READ (nin,*) tol, maxitn
       READ (nin,*) anorm, sigmax
       READ (nin,*) omega

!      Read the matrix A

       DO i = 1, nnz
          READ (nin,*) a(i), irow(i), icol(i)
       END DO

!      Read right-hand side vector b and initial approximate solution x

       READ (nin,*) b(1:n)
       READ (nin,*) x(1:n)

!      Call F11BDF to initialize solver

       weight = 'N'
       monit = 0

!      ifail: behaviour on error exit
!             =0 for hard exit, =1 for quiet-soft, =-1 for noisy-soft
       ifail = 0
       CALL f11bdf(method,precon,norm,weight,iterm,n,m,tol,maxitn,anorm, &
          sigmax,monit,lwneed,work,lwork,ifail)

!      Calculate reciprocal diagonal matrix elements if necessary

       IF (precon=='P' .OR. precon=='p') THEN

          iwork(1:n) = 0

          DO i = 1, nnz
             IF (irow(i)==icol(i)) THEN
                iwork(irow(i)) = iwork(irow(i)) + 1
                IF (a(i)/=0.E0_nag_wp) THEN
                   rdiag(irow(i)) = 1.E0_nag_wp/a(i)
                ELSE
                   WRITE (nout,*) 'Matrix has a zero diagonal element'
                   GO TO 20
                END IF
             END IF
          END DO

          DO i = 1, n
             IF (iwork(i)==0) THEN
                WRITE (nout,*) 'Matrix has a missing diagonal element'
                GO TO 20
             END IF
             IF (iwork(i)>=2) THEN
                WRITE (nout,*) 'Matrix has a multiple diagonal element'
                GO TO 20
             END IF
          END DO

       END IF

!      Call F11BEF to solve the linear system

       irevcm = 0
       ckxaf = 'C'
       ckddf = 'C'

LOOP:  DO
          ifail = 0
          CALL f11bef(irevcm,x,b,wgt,work,lwork,ifail)

          SELECT CASE (irevcm)
          CASE (1)
!            Compute matrix-vector product
             trans = 'N'

             CALL f11xaf(trans,n,nnz,a,irow,icol,ckxaf,x,b,ifail)

             ckxaf = 'N'
          CASE (-1)
!            Compute transposed matrix-vector product
             trans = 'T'

             CALL f11xaf(trans,n,nnz,a,irow,icol,ckxaf,x,b,ifail)

             ckxaf = 'N'
          CASE (2)

!            SSOR preconditioning
             trans = 'N'

             CALL f11ddf(trans,n,nnz,a,irow,icol,rdiag,omega,ckddf,x,b,iwork, &
                ifail)

             ckddf = 'N'
          CASE (4)
!            Termination

             ifail = 0
             CALL f11bff(itn,stplhs,stprhs,anorm,sigmax,work,lwork,ifail)

             WRITE (nout,'(A,I10,A)') ' Converged in', itn, ' iterations'
             WRITE (nout,'(A,1P,E16.3)') ' Matrix norm =', anorm
             WRITE (nout,'(A,1P,E16.3)') ' Final residual norm =', stplhs
             WRITE (nout,*)

!            Output x
             WRITE (nout,*) '           X'
             WRITE (nout,'(1X,1P,E16.4)') x(1:n)
             EXIT LOOP
          CASE DEFAULT
             EXIT LOOP
          END SELECT
       END DO LOOP

20     CONTINUE

    END PROGRAM f11ddfe