LAPACK 3.3.0
|
00001 SUBROUTINE SSYSV( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK, 00002 $ LWORK, INFO ) 00003 * 00004 * -- LAPACK driver routine (version 3.2.2) -- 00005 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00006 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00007 * May 2010 00008 * 00009 * .. Scalar Arguments .. 00010 CHARACTER UPLO 00011 INTEGER INFO, LDA, LDB, LWORK, N, NRHS 00012 * .. 00013 * .. Array Arguments .. 00014 INTEGER IPIV( * ) 00015 REAL A( LDA, * ), B( LDB, * ), WORK( * ) 00016 * .. 00017 * 00018 * Purpose 00019 * ======= 00020 * 00021 * SSYSV computes the solution to a real system of linear equations 00022 * A * X = B, 00023 * where A is an N-by-N symmetric matrix and X and B are N-by-NRHS 00024 * matrices. 00025 * 00026 * The diagonal pivoting method is used to factor A as 00027 * A = U * D * U**T, if UPLO = 'U', or 00028 * A = L * D * L**T, if UPLO = 'L', 00029 * where U (or L) is a product of permutation and unit upper (lower) 00030 * triangular matrices, and D is symmetric and block diagonal with 00031 * 1-by-1 and 2-by-2 diagonal blocks. The factored form of A is then 00032 * used to solve the system of equations A * X = B. 00033 * 00034 * Arguments 00035 * ========= 00036 * 00037 * UPLO (input) CHARACTER*1 00038 * = 'U': Upper triangle of A is stored; 00039 * = 'L': Lower triangle of A is stored. 00040 * 00041 * N (input) INTEGER 00042 * The number of linear equations, i.e., the order of the 00043 * matrix A. N >= 0. 00044 * 00045 * NRHS (input) INTEGER 00046 * The number of right hand sides, i.e., the number of columns 00047 * of the matrix B. NRHS >= 0. 00048 * 00049 * A (input/output) REAL array, dimension (LDA,N) 00050 * On entry, the symmetric matrix A. If UPLO = 'U', the leading 00051 * N-by-N upper triangular part of A contains the upper 00052 * triangular part of the matrix A, and the strictly lower 00053 * triangular part of A is not referenced. If UPLO = 'L', the 00054 * leading N-by-N lower triangular part of A contains the lower 00055 * triangular part of the matrix A, and the strictly upper 00056 * triangular part of A is not referenced. 00057 * 00058 * On exit, if INFO = 0, the block diagonal matrix D and the 00059 * multipliers used to obtain the factor U or L from the 00060 * factorization A = U*D*U**T or A = L*D*L**T as computed by 00061 * SSYTRF. 00062 * 00063 * LDA (input) INTEGER 00064 * The leading dimension of the array A. LDA >= max(1,N). 00065 * 00066 * IPIV (output) INTEGER array, dimension (N) 00067 * Details of the interchanges and the block structure of D, as 00068 * determined by SSYTRF. If IPIV(k) > 0, then rows and columns 00069 * k and IPIV(k) were interchanged, and D(k,k) is a 1-by-1 00070 * diagonal block. If UPLO = 'U' and IPIV(k) = IPIV(k-1) < 0, 00071 * then rows and columns k-1 and -IPIV(k) were interchanged and 00072 * D(k-1:k,k-1:k) is a 2-by-2 diagonal block. If UPLO = 'L' and 00073 * IPIV(k) = IPIV(k+1) < 0, then rows and columns k+1 and 00074 * -IPIV(k) were interchanged and D(k:k+1,k:k+1) is a 2-by-2 00075 * diagonal block. 00076 * 00077 * B (input/output) REAL array, dimension (LDB,NRHS) 00078 * On entry, the N-by-NRHS right hand side matrix B. 00079 * On exit, if INFO = 0, the N-by-NRHS solution matrix X. 00080 * 00081 * LDB (input) INTEGER 00082 * The leading dimension of the array B. LDB >= max(1,N). 00083 * 00084 * WORK (workspace/output) REAL array, dimension (MAX(1,LWORK)) 00085 * On exit, if INFO = 0, WORK(1) returns the optimal LWORK. 00086 * 00087 * LWORK (input) INTEGER 00088 * The length of WORK. LWORK >= 1, and for best performance 00089 * LWORK >= max(1,N*NB), where NB is the optimal blocksize for 00090 * SSYTRF. 00091 * 00092 * If LWORK = -1, then a workspace query is assumed; the routine 00093 * only calculates the optimal size of the WORK array, returns 00094 * this value as the first entry of the WORK array, and no error 00095 * message related to LWORK is issued by XERBLA. 00096 * 00097 * INFO (output) INTEGER 00098 * = 0: successful exit 00099 * < 0: if INFO = -i, the i-th argument had an illegal value 00100 * > 0: if INFO = i, D(i,i) is exactly zero. The factorization 00101 * has been completed, but the block diagonal matrix D is 00102 * exactly singular, so the solution could not be computed. 00103 * 00104 * ===================================================================== 00105 * 00106 * .. Local Scalars .. 00107 LOGICAL LQUERY 00108 INTEGER LWKOPT, NB 00109 * .. 00110 * .. External Functions .. 00111 LOGICAL LSAME 00112 INTEGER ILAENV 00113 EXTERNAL ILAENV, LSAME 00114 * .. 00115 * .. External Subroutines .. 00116 EXTERNAL SSYTRF, SSYTRS2, XERBLA 00117 * .. 00118 * .. Intrinsic Functions .. 00119 INTRINSIC MAX 00120 * .. 00121 * .. Executable Statements .. 00122 * 00123 * Test the input parameters. 00124 * 00125 INFO = 0 00126 LQUERY = ( LWORK.EQ.-1 ) 00127 IF( .NOT.LSAME( UPLO, 'U' ) .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN 00128 INFO = -1 00129 ELSE IF( N.LT.0 ) THEN 00130 INFO = -2 00131 ELSE IF( NRHS.LT.0 ) THEN 00132 INFO = -3 00133 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 00134 INFO = -5 00135 ELSE IF( LDB.LT.MAX( 1, N ) ) THEN 00136 INFO = -8 00137 ELSE IF( LWORK.LT.1 .AND. .NOT.LQUERY ) THEN 00138 INFO = -10 00139 END IF 00140 * 00141 IF( INFO.EQ.0 ) THEN 00142 IF( N.EQ.0 ) THEN 00143 LWKOPT = 1 00144 ELSE 00145 NB = ILAENV( 1, 'SSYTRF', UPLO, N, -1, -1, -1 ) 00146 LWKOPT = N*NB 00147 END IF 00148 WORK( 1 ) = LWKOPT 00149 END IF 00150 * 00151 IF( INFO.NE.0 ) THEN 00152 CALL XERBLA( 'SSYSV ', -INFO ) 00153 RETURN 00154 ELSE IF( LQUERY ) THEN 00155 RETURN 00156 END IF 00157 * 00158 * Compute the factorization A = U*D*U' or A = L*D*L'. 00159 * 00160 CALL SSYTRF( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO ) 00161 IF( INFO.EQ.0 ) THEN 00162 * 00163 * Solve the system A*X = B, overwriting B with X. 00164 * 00165 CALL SSYTRS2( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK, INFO ) 00166 * 00167 END IF 00168 * 00169 WORK( 1 ) = LWKOPT 00170 * 00171 RETURN 00172 * 00173 * End of SSYSV 00174 * 00175 END