00001 /* zspsv.f -- translated by f2c (version 20061008). 00002 You must link the resulting object file with libf2c: 00003 on Microsoft Windows system, link with libf2c.lib; 00004 on Linux or Unix systems, link with .../path/to/libf2c.a -lm 00005 or, if you install libf2c.a in a standard place, with -lf2c -lm 00006 -- in that order, at the end of the command line, as in 00007 cc *.o -lf2c -lm 00008 Source for libf2c is in /netlib/f2c/libf2c.zip, e.g., 00009 00010 http://www.netlib.org/f2c/libf2c.zip 00011 */ 00012 00013 #include "f2c.h" 00014 #include "blaswrap.h" 00015 00016 /* Subroutine */ int zspsv_(char *uplo, integer *n, integer *nrhs, 00017 doublecomplex *ap, integer *ipiv, doublecomplex *b, integer *ldb, 00018 integer *info) 00019 { 00020 /* System generated locals */ 00021 integer b_dim1, b_offset, i__1; 00022 00023 /* Local variables */ 00024 extern logical lsame_(char *, char *); 00025 extern /* Subroutine */ int xerbla_(char *, integer *), zsptrf_( 00026 char *, integer *, doublecomplex *, integer *, integer *), 00027 zsptrs_(char *, integer *, integer *, doublecomplex *, integer *, 00028 doublecomplex *, integer *, integer *); 00029 00030 00031 /* -- LAPACK driver routine (version 3.2) -- */ 00032 /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */ 00033 /* November 2006 */ 00034 00035 /* .. Scalar Arguments .. */ 00036 /* .. */ 00037 /* .. Array Arguments .. */ 00038 /* .. */ 00039 00040 /* Purpose */ 00041 /* ======= */ 00042 00043 /* ZSPSV computes the solution to a complex system of linear equations */ 00044 /* A * X = B, */ 00045 /* where A is an N-by-N symmetric matrix stored in packed format and X */ 00046 /* and B are N-by-NRHS matrices. */ 00047 00048 /* The diagonal pivoting method is used to factor A as */ 00049 /* A = U * D * U**T, if UPLO = 'U', or */ 00050 /* A = L * D * L**T, if UPLO = 'L', */ 00051 /* where U (or L) is a product of permutation and unit upper (lower) */ 00052 /* triangular matrices, D is symmetric and block diagonal with 1-by-1 */ 00053 /* and 2-by-2 diagonal blocks. The factored form of A is then used to */ 00054 /* solve the system of equations A * X = B. */ 00055 00056 /* Arguments */ 00057 /* ========= */ 00058 00059 /* UPLO (input) CHARACTER*1 */ 00060 /* = 'U': Upper triangle of A is stored; */ 00061 /* = 'L': Lower triangle of A is stored. */ 00062 00063 /* N (input) INTEGER */ 00064 /* The number of linear equations, i.e., the order of the */ 00065 /* matrix A. N >= 0. */ 00066 00067 /* NRHS (input) INTEGER */ 00068 /* The number of right hand sides, i.e., the number of columns */ 00069 /* of the matrix B. NRHS >= 0. */ 00070 00071 /* AP (input/output) COMPLEX*16 array, dimension (N*(N+1)/2) */ 00072 /* On entry, the upper or lower triangle of the symmetric matrix */ 00073 /* A, packed columnwise in a linear array. The j-th column of A */ 00074 /* is stored in the array AP as follows: */ 00075 /* if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j; */ 00076 /* if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n. */ 00077 /* See below for further details. */ 00078 00079 /* On exit, the block diagonal matrix D and the multipliers used */ 00080 /* to obtain the factor U or L from the factorization */ 00081 /* A = U*D*U**T or A = L*D*L**T as computed by ZSPTRF, stored as */ 00082 /* a packed triangular matrix in the same storage format as A. */ 00083 00084 /* IPIV (output) INTEGER array, dimension (N) */ 00085 /* Details of the interchanges and the block structure of D, as */ 00086 /* determined by ZSPTRF. If IPIV(k) > 0, then rows and columns */ 00087 /* k and IPIV(k) were interchanged, and D(k,k) is a 1-by-1 */ 00088 /* diagonal block. If UPLO = 'U' and IPIV(k) = IPIV(k-1) < 0, */ 00089 /* then rows and columns k-1 and -IPIV(k) were interchanged and */ 00090 /* D(k-1:k,k-1:k) is a 2-by-2 diagonal block. If UPLO = 'L' and */ 00091 /* IPIV(k) = IPIV(k+1) < 0, then rows and columns k+1 and */ 00092 /* -IPIV(k) were interchanged and D(k:k+1,k:k+1) is a 2-by-2 */ 00093 /* diagonal block. */ 00094 00095 /* B (input/output) COMPLEX*16 array, dimension (LDB,NRHS) */ 00096 /* On entry, the N-by-NRHS right hand side matrix B. */ 00097 /* On exit, if INFO = 0, the N-by-NRHS solution matrix X. */ 00098 00099 /* LDB (input) INTEGER */ 00100 /* The leading dimension of the array B. LDB >= max(1,N). */ 00101 00102 /* INFO (output) INTEGER */ 00103 /* = 0: successful exit */ 00104 /* < 0: if INFO = -i, the i-th argument had an illegal value */ 00105 /* > 0: if INFO = i, D(i,i) is exactly zero. The factorization */ 00106 /* has been completed, but the block diagonal matrix D is */ 00107 /* exactly singular, so the solution could not be */ 00108 /* computed. */ 00109 00110 /* Further Details */ 00111 /* =============== */ 00112 00113 /* The packed storage scheme is illustrated by the following example */ 00114 /* when N = 4, UPLO = 'U': */ 00115 00116 /* Two-dimensional storage of the symmetric matrix A: */ 00117 00118 /* a11 a12 a13 a14 */ 00119 /* a22 a23 a24 */ 00120 /* a33 a34 (aij = aji) */ 00121 /* a44 */ 00122 00123 /* Packed storage of the upper triangle of A: */ 00124 00125 /* AP = [ a11, a12, a22, a13, a23, a33, a14, a24, a34, a44 ] */ 00126 00127 /* ===================================================================== */ 00128 00129 /* .. External Functions .. */ 00130 /* .. */ 00131 /* .. External Subroutines .. */ 00132 /* .. */ 00133 /* .. Intrinsic Functions .. */ 00134 /* .. */ 00135 /* .. Executable Statements .. */ 00136 00137 /* Test the input parameters. */ 00138 00139 /* Parameter adjustments */ 00140 --ap; 00141 --ipiv; 00142 b_dim1 = *ldb; 00143 b_offset = 1 + b_dim1; 00144 b -= b_offset; 00145 00146 /* Function Body */ 00147 *info = 0; 00148 if (! lsame_(uplo, "U") && ! lsame_(uplo, "L")) { 00149 *info = -1; 00150 } else if (*n < 0) { 00151 *info = -2; 00152 } else if (*nrhs < 0) { 00153 *info = -3; 00154 } else if (*ldb < max(1,*n)) { 00155 *info = -7; 00156 } 00157 if (*info != 0) { 00158 i__1 = -(*info); 00159 xerbla_("ZSPSV ", &i__1); 00160 return 0; 00161 } 00162 00163 /* Compute the factorization A = U*D*U' or A = L*D*L'. */ 00164 00165 zsptrf_(uplo, n, &ap[1], &ipiv[1], info); 00166 if (*info == 0) { 00167 00168 /* Solve the system A*X = B, overwriting B with X. */ 00169 00170 zsptrs_(uplo, n, nrhs, &ap[1], &ipiv[1], &b[b_offset], ldb, info); 00171 00172 } 00173 return 0; 00174 00175 /* End of ZSPSV */ 00176 00177 } /* zspsv_ */