 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE CLATSP( UPLO, N, X, ISEED ) 00002 * 00003 * -- LAPACK auxiliary test routine (version 3.1) -- 00004 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00005 * November 2006 00006 * 00007 * .. Scalar Arguments .. 00008 CHARACTER UPLO 00009 INTEGER N 00010 * .. 00011 * .. Array Arguments .. 00012 INTEGER ISEED( * ) 00013 COMPLEX X( * ) 00014 * .. 00015 * 00016 * Purpose 00017 * ======= 00018 * 00019 * CLATSP generates a special test matrix for the complex symmetric 00020 * (indefinite) factorization for packed matrices. The pivot blocks of 00021 * the generated matrix will be in the following order: 00022 * 2x2 pivot block, non diagonalizable 00023 * 1x1 pivot block 00024 * 2x2 pivot block, diagonalizable 00025 * (cycle repeats) 00026 * A row interchange is required for each non-diagonalizable 2x2 block. 00027 * 00028 * Arguments 00029 * ========= 00030 * 00031 * UPLO (input) CHARACTER 00032 * Specifies whether the generated matrix is to be upper or 00033 * lower triangular. 00034 * = 'U': Upper triangular 00035 * = 'L': Lower triangular 00036 * 00037 * N (input) INTEGER 00038 * The dimension of the matrix to be generated. 00039 * 00040 * X (output) COMPLEX array, dimension (N*(N+1)/2) 00041 * The generated matrix in packed storage format. The matrix 00042 * consists of 3x3 and 2x2 diagonal blocks which result in the 00043 * pivot sequence given above. The matrix outside these 00044 * diagonal blocks is zero. 00045 * 00046 * ISEED (input/output) INTEGER array, dimension (4) 00047 * On entry, the seed for the random number generator. The last 00048 * of the four integers must be odd. (modified on exit) 00049 * 00050 * ===================================================================== 00051 * 00052 * .. Parameters .. 00053 COMPLEX EYE 00054 PARAMETER ( EYE = ( 0.0, 1.0 ) ) 00055 * .. 00056 * .. Local Scalars .. 00057 INTEGER J, JJ, N5 00058 REAL ALPHA, ALPHA3, BETA 00059 COMPLEX A, B, C, R 00060 * .. 00061 * .. External Functions .. 00062 COMPLEX CLARND 00063 EXTERNAL CLARND 00064 * .. 00065 * .. Intrinsic Functions .. 00066 INTRINSIC ABS, SQRT 00067 * .. 00068 * .. Executable Statements .. 00069 * 00070 * Initialize constants 00071 * 00072 ALPHA = ( 1.+SQRT( 17. ) ) / 8. 00073 BETA = ALPHA - 1. / 1000. 00074 ALPHA3 = ALPHA*ALPHA*ALPHA 00075 * 00076 * Fill the matrix with zeros. 00077 * 00078 DO 10 J = 1, N*( N+1 ) / 2 00079 X( J ) = 0.0 00080 10 CONTINUE 00081 * 00082 * UPLO = 'U': Upper triangular storage 00083 * 00084 IF( UPLO.EQ.'U' ) THEN 00085 N5 = N / 5 00086 N5 = N - 5*N5 + 1 00087 * 00088 JJ = N*( N+1 ) / 2 00089 DO 20 J = N, N5, -5 00090 A = ALPHA3*CLARND( 5, ISEED ) 00091 B = CLARND( 5, ISEED ) / ALPHA 00092 C = A - 2.*B*EYE 00093 R = C / BETA 00094 X( JJ ) = A 00095 X( JJ-2 ) = B 00096 JJ = JJ - J 00097 X( JJ ) = CLARND( 2, ISEED ) 00098 X( JJ-1 ) = R 00099 JJ = JJ - ( J-1 ) 00100 X( JJ ) = C 00101 JJ = JJ - ( J-2 ) 00102 X( JJ ) = CLARND( 2, ISEED ) 00103 JJ = JJ - ( J-3 ) 00104 X( JJ ) = CLARND( 2, ISEED ) 00105 IF( ABS( X( JJ+( J-3 ) ) ).GT.ABS( X( JJ ) ) ) THEN 00106 X( JJ+( J-4 ) ) = 2.0*X( JJ+( J-3 ) ) 00107 ELSE 00108 X( JJ+( J-4 ) ) = 2.0*X( JJ ) 00109 END IF 00110 JJ = JJ - ( J-4 ) 00111 20 CONTINUE 00112 * 00113 * Clean-up for N not a multiple of 5. 00114 * 00115 J = N5 - 1 00116 IF( J.GT.2 ) THEN 00117 A = ALPHA3*CLARND( 5, ISEED ) 00118 B = CLARND( 5, ISEED ) / ALPHA 00119 C = A - 2.*B*EYE 00120 R = C / BETA 00121 X( JJ ) = A 00122 X( JJ-2 ) = B 00123 JJ = JJ - J 00124 X( JJ ) = CLARND( 2, ISEED ) 00125 X( JJ-1 ) = R 00126 JJ = JJ - ( J-1 ) 00127 X( JJ ) = C 00128 JJ = JJ - ( J-2 ) 00129 J = J - 3 00130 END IF 00131 IF( J.GT.1 ) THEN 00132 X( JJ ) = CLARND( 2, ISEED ) 00133 X( JJ-J ) = CLARND( 2, ISEED ) 00134 IF( ABS( X( JJ ) ).GT.ABS( X( JJ-J ) ) ) THEN 00135 X( JJ-1 ) = 2.0*X( JJ ) 00136 ELSE 00137 X( JJ-1 ) = 2.0*X( JJ-J ) 00138 END IF 00139 JJ = JJ - J - ( J-1 ) 00140 J = J - 2 00141 ELSE IF( J.EQ.1 ) THEN 00142 X( JJ ) = CLARND( 2, ISEED ) 00143 J = J - 1 00144 END IF 00145 * 00146 * UPLO = 'L': Lower triangular storage 00147 * 00148 ELSE 00149 N5 = N / 5 00150 N5 = N5*5 00151 * 00152 JJ = 1 00153 DO 30 J = 1, N5, 5 00154 A = ALPHA3*CLARND( 5, ISEED ) 00155 B = CLARND( 5, ISEED ) / ALPHA 00156 C = A - 2.*B*EYE 00157 R = C / BETA 00158 X( JJ ) = A 00159 X( JJ+2 ) = B 00160 JJ = JJ + ( N-J+1 ) 00161 X( JJ ) = CLARND( 2, ISEED ) 00162 X( JJ+1 ) = R 00163 JJ = JJ + ( N-J ) 00164 X( JJ ) = C 00165 JJ = JJ + ( N-J-1 ) 00166 X( JJ ) = CLARND( 2, ISEED ) 00167 JJ = JJ + ( N-J-2 ) 00168 X( JJ ) = CLARND( 2, ISEED ) 00169 IF( ABS( X( JJ-( N-J-2 ) ) ).GT.ABS( X( JJ ) ) ) THEN 00170 X( JJ-( N-J-2 )+1 ) = 2.0*X( JJ-( N-J-2 ) ) 00171 ELSE 00172 X( JJ-( N-J-2 )+1 ) = 2.0*X( JJ ) 00173 END IF 00174 JJ = JJ + ( N-J-3 ) 00175 30 CONTINUE 00176 * 00177 * Clean-up for N not a multiple of 5. 00178 * 00179 J = N5 + 1 00180 IF( J.LT.N-1 ) THEN 00181 A = ALPHA3*CLARND( 5, ISEED ) 00182 B = CLARND( 5, ISEED ) / ALPHA 00183 C = A - 2.*B*EYE 00184 R = C / BETA 00185 X( JJ ) = A 00186 X( JJ+2 ) = B 00187 JJ = JJ + ( N-J+1 ) 00188 X( JJ ) = CLARND( 2, ISEED ) 00189 X( JJ+1 ) = R 00190 JJ = JJ + ( N-J ) 00191 X( JJ ) = C 00192 JJ = JJ + ( N-J-1 ) 00193 J = J + 3 00194 END IF 00195 IF( J.LT.N ) THEN 00196 X( JJ ) = CLARND( 2, ISEED ) 00197 X( JJ+( N-J+1 ) ) = CLARND( 2, ISEED ) 00198 IF( ABS( X( JJ ) ).GT.ABS( X( JJ+( N-J+1 ) ) ) ) THEN 00199 X( JJ+1 ) = 2.0*X( JJ ) 00200 ELSE 00201 X( JJ+1 ) = 2.0*X( JJ+( N-J+1 ) ) 00202 END IF 00203 JJ = JJ + ( N-J+1 ) + ( N-J ) 00204 J = J + 2 00205 ELSE IF( J.EQ.N ) THEN 00206 X( JJ ) = CLARND( 2, ISEED ) 00207 JJ = JJ + ( N-J+1 ) 00208 J = J + 1 00209 END IF 00210 END IF 00211 * 00212 RETURN 00213 * 00214 * End of CLATSP 00215 * 00216 END
1.7.3 
