PROGRAM f04zcfe ! F04ZCF Example Program Text ! Mark 23 Release. NAG Copyright 2011. ! .. Use Statements .. USE nag_library, ONLY : f04zcf, f06ubf, nag_wp, zgbtrf, zgbtrs ! .. Implicit None Statement .. IMPLICIT NONE ! .. Parameters .. INTEGER, PARAMETER :: nin = 5, nout = 6 ! .. Local Scalars .. REAL (KIND=nag_wp) :: anorm, cond, estnrm INTEGER :: i, icase, ifail, info, j, k, kl, ku, & lda, ldx, n, nrhs ! .. Local Arrays .. COMPLEX (KIND=nag_wp), ALLOCATABLE :: a(:,:), work(:), x(:,:) REAL (KIND=nag_wp) :: rwork(1) INTEGER, ALLOCATABLE :: ipiv(:) ! .. Intrinsic Functions .. INTRINSIC max, min ! .. Executable Statements .. WRITE (nout,*) 'F04ZCF Example Program Results' ! Skip heading in data file READ (nin,*) READ (nin,*) n, kl, ku, nrhs lda = 2*kl + ku + 1 ldx = n ALLOCATE (a(lda,n),work(n),x(ldx,nrhs),ipiv(n)) k = kl + ku + 1 READ (nin,*) ((a(k+i-j,j),j=max(i-kl,1),min(i+ku,n)),i=1,n) ! First compute the 1-norm of A. anorm = f06ubf('1-norm',n,kl,ku,a(kl+1,1),lda,rwork) WRITE (nout,*) WRITE (nout,99999) 'Computed norm of A =', anorm ! Next estimate the 1-norm of inverse(A). ! Factorise A into P*L*U. ! The NAG name equivalent of zgbtrf is f07brf CALL zgbtrf(n,n,kl,ku,a,lda,ipiv,info) icase = 0 LOOP: DO ! ifail: behaviour on error exit ! =0 for hard exit, =1 for quiet-soft, =-1 for noisy-soft ifail = 0 CALL f04zcf(icase,n,x,estnrm,work,ifail) IF (icase/=0) THEN ! The NAG name equivalent of the backsolve routine zgbtrs is f07bsf IF (icase==1) THEN ! Return X := inv(A)*X by solving A*Y = X, overwriting ! Y on X. CALL zgbtrs('No transpose',n,kl,ku,nrhs,a,lda,ipiv,x,ldx,info) ELSE IF (icase==2) THEN ! Return X := conjg(inv(A)')*X by solving conjg(A')*Y ! = X, overwriting Y on X. CALL zgbtrs('Conjugate transpose',n,kl,ku,nrhs,a,lda,ipiv,x, & ldx,info) END IF ! Continue until icase is returned as 0. ELSE WRITE (nout,99999) 'Estimated norm of inverse(A) =', estnrm cond = anorm*estnrm WRITE (nout,99998) 'Estimated condition number of A =', cond EXIT LOOP END IF END DO LOOP 99999 FORMAT (1X,A,F8.4) 99998 FORMAT (1X,A,F6.1) END PROGRAM f04zcfe