Program f04zdfe ! f04zdf Example Program Text ! Mark 24 Release. NAG Copyright 2012. ! .. Use Statements .. Use nag_library, Only: f04zdf, f06uaf, nag_wp, zgetrf, zgetrs ! .. Implicit None Statement .. Implicit None ! .. Parameters .. Integer, Parameter :: nin = 5, nout = 6 ! .. Local Scalars .. Real (Kind=nag_wp) :: cond, nrma, nrminv Integer :: i, ifail, irevcm, lda, ldx, ldy, m, & n, seed, t ! .. Local Arrays .. Complex (Kind=nag_wp), Allocatable :: a(:,:), work(:), x(:,:), y(:,:) Real (Kind=nag_wp), Allocatable :: rwork(:) Real (Kind=nag_wp) :: workr(1) Integer, Allocatable :: ipiv(:), iwork(:) ! .. Executable Statements .. Write (nout,*) 'F04ZDF Example Program Results' Write (nout,*) ! Skip heading in data file Read (nin,*) Read (nin,*) m, n, t lda = m ldx = n ldy = m Allocate (a(lda,n)) Allocate (x(ldx,t)) Allocate (y(ldy,t)) Allocate (work(m*t)) Allocate (rwork(2*n)) Allocate (iwork(2*n+5*t+20)) Allocate (ipiv(n)) ! Read A from data file Read (nin,*)(a(i,1:n),i=1,m) ! Compute 1-norm of A nrma = f06uaf('1',m,n,a,lda,workr) Write (nout,99999) 'The norm of A is: ', nrma ! Estimate the norm of A^(-1) without forming A^(-1) irevcm = 0 ifail = 0 seed = 652 ! Perform an LU factorization so that A=LU where L and U are lower ! and upper triangular. Call zgetrf(m,n,a,lda,ipiv,ifail) loop: Do Call f04zdf(irevcm,m,n,x,ldx,y,ldy,nrminv,t,seed,work,rwork,iwork, & ifail) If (irevcm/=0) Then If (irevcm==1) Then ! Compute Y = inv(A)*X Call zgetrs('N',n,t,a,lda,ipiv,x,ldx,ifail) ! X was overwritten by ZGETRS, so set Y=X y(1:n,1:t) = x(1:n,1:t) Else ! Compute X = herm(inv(A))*Y Call zgetrs('C',n,t,a,lda,ipiv,y,ldy,ifail) ! Y was overwritten by ZGETRS, so set X=Y x(1:n,1:t) = y(1:n,1:t) End If Else Write (nout,99999) 'The estimated norm of inverse(A) is: ', nrminv ! Compute and print the estimated condition number cond = nrminv*nrma Write (nout,*) Write (nout,99999) 'Estimated condition number of A: ', cond Write (nout,*) Exit loop End If End Do loop 99999 Format (1X,A,F6.2) End Program f04zdfe