LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ cppt03()

subroutine cppt03 ( character  uplo,
integer  n,
complex, dimension( * )  a,
complex, dimension( * )  ainv,
complex, dimension( ldwork, * )  work,
integer  ldwork,
real, dimension( * )  rwork,
real  rcond,
real  resid 
)

CPPT03

Purpose:
 CPPT03 computes the residual for a Hermitian packed matrix times its
 inverse:
    norm( I - A*AINV ) / ( N * norm(A) * norm(AINV) * EPS ),
 where EPS is the machine epsilon.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          Specifies whether the upper or lower triangular part of the
          Hermitian matrix A is stored:
          = 'U':  Upper triangular
          = 'L':  Lower triangular
[in]N
          N is INTEGER
          The number of rows and columns of the matrix A.  N >= 0.
[in]A
          A is COMPLEX array, dimension (N*(N+1)/2)
          The original Hermitian matrix A, stored as a packed
          triangular matrix.
[in]AINV
          AINV is COMPLEX array, dimension (N*(N+1)/2)
          The (Hermitian) inverse of the matrix A, stored as a packed
          triangular matrix.
[out]WORK
          WORK is COMPLEX array, dimension (LDWORK,N)
[in]LDWORK
          LDWORK is INTEGER
          The leading dimension of the array WORK.  LDWORK >= max(1,N).
[out]RWORK
          RWORK is REAL array, dimension (N)
[out]RCOND
          RCOND is REAL
          The reciprocal of the condition number of A, computed as
          ( 1/norm(A) ) / norm(AINV).
[out]RESID
          RESID is REAL
          norm(I - A*AINV) / ( N * norm(A) * norm(AINV) * EPS )
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 108 of file cppt03.f.

110*
111* -- LAPACK test routine --
112* -- LAPACK is a software package provided by Univ. of Tennessee, --
113* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
114*
115* .. Scalar Arguments ..
116 CHARACTER UPLO
117 INTEGER LDWORK, N
118 REAL RCOND, RESID
119* ..
120* .. Array Arguments ..
121 REAL RWORK( * )
122 COMPLEX A( * ), AINV( * ), WORK( LDWORK, * )
123* ..
124*
125* =====================================================================
126*
127* .. Parameters ..
128 REAL ZERO, ONE
129 parameter( zero = 0.0e+0, one = 1.0e+0 )
130 COMPLEX CZERO, CONE
131 parameter( czero = ( 0.0e+0, 0.0e+0 ),
132 $ cone = ( 1.0e+0, 0.0e+0 ) )
133* ..
134* .. Local Scalars ..
135 INTEGER I, J, JJ
136 REAL AINVNM, ANORM, EPS
137* ..
138* .. External Functions ..
139 LOGICAL LSAME
140 REAL CLANGE, CLANHP, SLAMCH
141 EXTERNAL lsame, clange, clanhp, slamch
142* ..
143* .. Intrinsic Functions ..
144 INTRINSIC conjg, real
145* ..
146* .. External Subroutines ..
147 EXTERNAL ccopy, chpmv
148* ..
149* .. Executable Statements ..
150*
151* Quick exit if N = 0.
152*
153 IF( n.LE.0 ) THEN
154 rcond = one
155 resid = zero
156 RETURN
157 END IF
158*
159* Exit with RESID = 1/EPS if ANORM = 0 or AINVNM = 0.
160*
161 eps = slamch( 'Epsilon' )
162 anorm = clanhp( '1', uplo, n, a, rwork )
163 ainvnm = clanhp( '1', uplo, n, ainv, rwork )
164 IF( anorm.LE.zero .OR. ainvnm.LE.zero ) THEN
165 rcond = zero
166 resid = one / eps
167 RETURN
168 END IF
169 rcond = ( one/anorm ) / ainvnm
170*
171* UPLO = 'U':
172* Copy the leading N-1 x N-1 submatrix of AINV to WORK(1:N,2:N) and
173* expand it to a full matrix, then multiply by A one column at a
174* time, moving the result one column to the left.
175*
176 IF( lsame( uplo, 'U' ) ) THEN
177*
178* Copy AINV
179*
180 jj = 1
181 DO 20 j = 1, n - 1
182 CALL ccopy( j, ainv( jj ), 1, work( 1, j+1 ), 1 )
183 DO 10 i = 1, j - 1
184 work( j, i+1 ) = conjg( ainv( jj+i-1 ) )
185 10 CONTINUE
186 jj = jj + j
187 20 CONTINUE
188 jj = ( ( n-1 )*n ) / 2 + 1
189 DO 30 i = 1, n - 1
190 work( n, i+1 ) = conjg( ainv( jj+i-1 ) )
191 30 CONTINUE
192*
193* Multiply by A
194*
195 DO 40 j = 1, n - 1
196 CALL chpmv( 'Upper', n, -cone, a, work( 1, j+1 ), 1, czero,
197 $ work( 1, j ), 1 )
198 40 CONTINUE
199 CALL chpmv( 'Upper', n, -cone, a, ainv( jj ), 1, czero,
200 $ work( 1, n ), 1 )
201*
202* UPLO = 'L':
203* Copy the trailing N-1 x N-1 submatrix of AINV to WORK(1:N,1:N-1)
204* and multiply by A, moving each column to the right.
205*
206 ELSE
207*
208* Copy AINV
209*
210 DO 50 i = 1, n - 1
211 work( 1, i ) = conjg( ainv( i+1 ) )
212 50 CONTINUE
213 jj = n + 1
214 DO 70 j = 2, n
215 CALL ccopy( n-j+1, ainv( jj ), 1, work( j, j-1 ), 1 )
216 DO 60 i = 1, n - j
217 work( j, j+i-1 ) = conjg( ainv( jj+i ) )
218 60 CONTINUE
219 jj = jj + n - j + 1
220 70 CONTINUE
221*
222* Multiply by A
223*
224 DO 80 j = n, 2, -1
225 CALL chpmv( 'Lower', n, -cone, a, work( 1, j-1 ), 1, czero,
226 $ work( 1, j ), 1 )
227 80 CONTINUE
228 CALL chpmv( 'Lower', n, -cone, a, ainv( 1 ), 1, czero,
229 $ work( 1, 1 ), 1 )
230*
231 END IF
232*
233* Add the identity matrix to WORK .
234*
235 DO 90 i = 1, n
236 work( i, i ) = work( i, i ) + cone
237 90 CONTINUE
238*
239* Compute norm(I - A*AINV) / (N * norm(A) * norm(AINV) * EPS)
240*
241 resid = clange( '1', n, n, work, ldwork, rwork )
242*
243 resid = ( ( resid*rcond )/eps ) / real( n )
244*
245 RETURN
246*
247* End of CPPT03
248*
subroutine ccopy(n, cx, incx, cy, incy)
CCOPY
Definition ccopy.f:81
subroutine chpmv(uplo, n, alpha, ap, x, incx, beta, y, incy)
CHPMV
Definition chpmv.f:149
real function slamch(cmach)
SLAMCH
Definition slamch.f:68
real function clange(norm, m, n, a, lda, work)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition clange.f:115
real function clanhp(norm, uplo, n, ap, work)
CLANHP returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition clanhp.f:117
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
Here is the call graph for this function:
Here is the caller graph for this function: