LAPACK  3.8.0
LAPACK: Linear Algebra PACKage

◆ cpot01()

subroutine cpot01 ( character  UPLO,
integer  N,
complex, dimension( lda, * )  A,
integer  LDA,
complex, dimension( ldafac, * )  AFAC,
integer  LDAFAC,
real, dimension( * )  RWORK,
real  RESID 
)

CPOT01

Purpose:
 CPOT01 reconstructs a Hermitian positive definite matrix  A  from
 its L*L' or U'*U factorization and computes the residual
    norm( L*L' - A ) / ( N * norm(A) * EPS ) or
    norm( U'*U - A ) / ( N * norm(A) * EPS ),
 where EPS is the machine epsilon, L' is the conjugate transpose of L,
 and U' is the conjugate transpose of U.
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 (LDA,N)
          The original Hermitian matrix A.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N)
[in,out]AFAC
          AFAC is COMPLEX array, dimension (LDAFAC,N)
          On entry, the factor L or U from the L*L' or U'*U
          factorization of A.
          Overwritten with the reconstructed matrix, and then with the
          difference L*L' - A (or U'*U - A).
[in]LDAFAC
          LDAFAC is INTEGER
          The leading dimension of the array AFAC.  LDAFAC >= max(1,N).
[out]RWORK
          RWORK is REAL array, dimension (N)
[out]RESID
          RESID is REAL
          If UPLO = 'L', norm(L*L' - A) / ( N * norm(A) * EPS )
          If UPLO = 'U', norm(U'*U - A) / ( N * norm(A) * EPS )
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
December 2016

Definition at line 108 of file cpot01.f.

108 *
109 * -- LAPACK test routine (version 3.7.0) --
110 * -- LAPACK is a software package provided by Univ. of Tennessee, --
111 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
112 * December 2016
113 *
114 * .. Scalar Arguments ..
115  CHARACTER uplo
116  INTEGER lda, ldafac, n
117  REAL resid
118 * ..
119 * .. Array Arguments ..
120  REAL rwork( * )
121  COMPLEX a( lda, * ), afac( ldafac, * )
122 * ..
123 *
124 * =====================================================================
125 *
126 * .. Parameters ..
127  REAL zero, one
128  parameter( zero = 0.0e+0, one = 1.0e+0 )
129 * ..
130 * .. Local Scalars ..
131  INTEGER i, j, k
132  REAL anorm, eps, tr
133  COMPLEX tc
134 * ..
135 * .. External Functions ..
136  LOGICAL lsame
137  REAL clanhe, slamch
138  COMPLEX cdotc
139  EXTERNAL lsame, clanhe, slamch, cdotc
140 * ..
141 * .. External Subroutines ..
142  EXTERNAL cher, cscal, ctrmv
143 * ..
144 * .. Intrinsic Functions ..
145  INTRINSIC aimag, real
146 * ..
147 * .. Executable Statements ..
148 *
149 * Quick exit if N = 0.
150 *
151  IF( n.LE.0 ) THEN
152  resid = zero
153  RETURN
154  END IF
155 *
156 * Exit with RESID = 1/EPS if ANORM = 0.
157 *
158  eps = slamch( 'Epsilon' )
159  anorm = clanhe( '1', uplo, n, a, lda, rwork )
160  IF( anorm.LE.zero ) THEN
161  resid = one / eps
162  RETURN
163  END IF
164 *
165 * Check the imaginary parts of the diagonal elements and return with
166 * an error code if any are nonzero.
167 *
168  DO 10 j = 1, n
169  IF( aimag( afac( j, j ) ).NE.zero ) THEN
170  resid = one / eps
171  RETURN
172  END IF
173  10 CONTINUE
174 *
175 * Compute the product U'*U, overwriting U.
176 *
177  IF( lsame( uplo, 'U' ) ) THEN
178  DO 20 k = n, 1, -1
179 *
180 * Compute the (K,K) element of the result.
181 *
182  tr = cdotc( k, afac( 1, k ), 1, afac( 1, k ), 1 )
183  afac( k, k ) = tr
184 *
185 * Compute the rest of column K.
186 *
187  CALL ctrmv( 'Upper', 'Conjugate', 'Non-unit', k-1, afac,
188  $ ldafac, afac( 1, k ), 1 )
189 *
190  20 CONTINUE
191 *
192 * Compute the product L*L', overwriting L.
193 *
194  ELSE
195  DO 30 k = n, 1, -1
196 *
197 * Add a multiple of column K of the factor L to each of
198 * columns K+1 through N.
199 *
200  IF( k+1.LE.n )
201  $ CALL cher( 'Lower', n-k, one, afac( k+1, k ), 1,
202  $ afac( k+1, k+1 ), ldafac )
203 *
204 * Scale column K by the diagonal element.
205 *
206  tc = afac( k, k )
207  CALL cscal( n-k+1, tc, afac( k, k ), 1 )
208 *
209  30 CONTINUE
210  END IF
211 *
212 * Compute the difference L*L' - A (or U'*U - A).
213 *
214  IF( lsame( uplo, 'U' ) ) THEN
215  DO 50 j = 1, n
216  DO 40 i = 1, j - 1
217  afac( i, j ) = afac( i, j ) - a( i, j )
218  40 CONTINUE
219  afac( j, j ) = afac( j, j ) - REAL( A( J, J ) )
220  50 CONTINUE
221  ELSE
222  DO 70 j = 1, n
223  afac( j, j ) = afac( j, j ) - REAL( A( J, J ) )
224  DO 60 i = j + 1, n
225  afac( i, j ) = afac( i, j ) - a( i, j )
226  60 CONTINUE
227  70 CONTINUE
228  END IF
229 *
230 * Compute norm( L*U - A ) / ( N * norm(A) * EPS )
231 *
232  resid = clanhe( '1', uplo, n, afac, ldafac, rwork )
233 *
234  resid = ( ( resid / REAL( N ) ) / anorm ) / eps
235 *
236  RETURN
237 *
238 * End of CPOT01
239 *
subroutine cher(UPLO, N, ALPHA, X, INCX, A, LDA)
CHER
Definition: cher.f:137
subroutine cscal(N, CA, CX, INCX)
CSCAL
Definition: cscal.f:80
complex function cdotc(N, CX, INCX, CY, INCY)
CDOTC
Definition: cdotc.f:85
subroutine ctrmv(UPLO, TRANS, DIAG, N, A, LDA, X, INCX)
CTRMV
Definition: ctrmv.f:149
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
real function clanhe(NORM, UPLO, N, A, LDA, WORK)
CLANHE returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a complex Hermitian matrix.
Definition: clanhe.f:126
Here is the call graph for this function:
Here is the caller graph for this function: