Program f08ztfe ! F08ZTF Example Program Text ! Mark 24 Release. NAG Copyright 2012. ! .. Use Statements .. Use nag_library, Only: dznrm2, nag_wp, zgemv, zggrqf, ztrmv, ztrtrs, & zunmqr, zunmrq ! .. Implicit None Statement .. Implicit None ! .. Parameters .. Complex (Kind=nag_wp), Parameter :: one = (1.0E0_nag_wp,0.0E0_nag_wp) Integer, Parameter :: nb = 64, nin = 5, nout = 6 ! .. Local Scalars .. Real (Kind=nag_wp) :: rnorm Integer :: i, info, lda, ldb, lwork, m, n, p ! .. Local Arrays .. Complex (Kind=nag_wp), Allocatable :: a(:,:), b(:,:), c(:), d(:), & taua(:), taub(:), work(:), x(:) ! .. Intrinsic Procedures .. Intrinsic :: min ! .. Executable Statements .. Write (nout,*) 'F08ZTF Example Program Results' Write (nout,*) ! Skip heading in data file Read (nin,*) Read (nin,*) p, n, m lda = m ldb = p lwork = nb*(p+n) Allocate (a(lda,n),b(ldb,n),c(p),d(m),taua(min(m,n)),taub(min(p, & n)),work(lwork),x(n)) ! Read B, A, C and D from data file Read (nin,*)(b(i,1:n),i=1,p) Read (nin,*)(a(i,1:n),i=1,m) Read (nin,*) c(1:p) Read (nin,*) d(1:m) ! Compute the generalized RQ factorization of (B,A) as ! A = (0 R12)*Q, B = Z*(T11 T12 T13)*Q, where R12, T11 and T22 ! ( 0 T22 T23) ! are upper triangular ! The NAG name equivalent of zggrqf is f08ztf Call zggrqf(m,p,n,a,lda,taua,b,ldb,taub,work,lwork,info) ! Compute (f1) = (Z**H)*c, storing the result in C ! (f2) ! The NAG name equivalent of zunmqr is f08auf Call zunmqr('Left','Conjugate transpose',p,1,min(p,n),b,ldb,taub,c,p, & work,lwork,info) ! Putting Q*x = (y1), solve R12*w = d for w, storing result in D ! (w ) ! The NAG name equivalent of ztrtrs is f07tsf Call ztrtrs('Upper','No transpose','Non-unit',m,1,a(1,n-m+1),lda,d,m, & info) If (info>0) Then Write (nout,*) 'The upper triangular factor, R12, of A is singular, ' Write (nout,*) 'the least squares solution could not be computed' Else ! Form f1 - T1*w, T1 = (T12 T13), in C ! The NAG name equivalent of zgemv is f06saf Call zgemv('No transpose',n-m,m,-one,b(1,n-m+1),ldb,d,1,one,c,1) ! Solve T11*y1 = f1 - T1*w for y1, storing result in C ! The NAG name equivalent of ztrtrs is f07tsf Call ztrtrs('Upper','No transpose','Non-unit',n-m,1,b,ldb,c,n-m,info) If (info>0) Then Write (nout,*) & 'The upper triangular factor, T11, of B is singular, ' Write (nout,*) 'the least squares solution could not be computed' Else ! Copy y into X (first y1, then w) x(1:n-m) = c(1:n-m) x(n-m+1:n) = d(1:m) ! Compute x = (Q**H)*y ! The NAG name equivalent of zunmrq is f08cxf Call zunmrq('Left','Conjugate transpose',n,1,m,a,lda,taua,x,n,work, & lwork,info) ! Putting w = (y2), form f2 - T22*y2 - T23*y3 ! (y3) ! T22*y2 ! The NAG name equivalent of ztrmv is f06sff Call ztrmv('Upper','No transpose','Non-unit',min(p,n)-n+m, & b(n-m+1,n-m+1),ldb,d,1) ! f2 - T22*y2 Do i = 1, min(p,n) - n + m c(n-m+i) = c(n-m+i) - d(i) End Do If (p