diff options
author | jofret | 2009-04-28 07:17:00 +0000 |
---|---|---|
committer | jofret | 2009-04-28 07:17:00 +0000 |
commit | 8c8d2f518968ce7057eec6aa5cd5aec8faab861a (patch) | |
tree | 3dd1788b71d6a3ce2b73d2d475a3133580e17530 /src/lib/lapack/dsptrd.f | |
parent | 9f652ffc16a310ac6641a9766c5b9e2671e0e9cb (diff) | |
download | scilab2c-8c8d2f518968ce7057eec6aa5cd5aec8faab861a.tar.gz scilab2c-8c8d2f518968ce7057eec6aa5cd5aec8faab861a.tar.bz2 scilab2c-8c8d2f518968ce7057eec6aa5cd5aec8faab861a.zip |
Moving lapack to right place
Diffstat (limited to 'src/lib/lapack/dsptrd.f')
-rw-r--r-- | src/lib/lapack/dsptrd.f | 228 |
1 files changed, 0 insertions, 228 deletions
diff --git a/src/lib/lapack/dsptrd.f b/src/lib/lapack/dsptrd.f deleted file mode 100644 index 6d3390e3..00000000 --- a/src/lib/lapack/dsptrd.f +++ /dev/null @@ -1,228 +0,0 @@ - SUBROUTINE DSPTRD( UPLO, N, AP, D, E, TAU, INFO ) -* -* -- LAPACK routine (version 3.1) -- -* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. -* November 2006 -* -* .. Scalar Arguments .. - CHARACTER UPLO - INTEGER INFO, N -* .. -* .. Array Arguments .. - DOUBLE PRECISION AP( * ), D( * ), E( * ), TAU( * ) -* .. -* -* Purpose -* ======= -* -* DSPTRD reduces a real symmetric matrix A stored in packed form to -* symmetric tridiagonal form T by an orthogonal similarity -* transformation: Q**T * A * Q = T. -* -* Arguments -* ========= -* -* UPLO (input) CHARACTER*1 -* = 'U': Upper triangle of A is stored; -* = 'L': Lower triangle of A is stored. -* -* N (input) INTEGER -* The order of the matrix A. N >= 0. -* -* AP (input/output) DOUBLE PRECISION array, dimension (N*(N+1)/2) -* On entry, the upper or lower triangle of the symmetric matrix -* A, packed columnwise in a linear array. The j-th column of A -* is stored in the array AP as follows: -* if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j; -* if UPLO = 'L', AP(i + (j-1)*(2*n-j)/2) = A(i,j) for j<=i<=n. -* On exit, if UPLO = 'U', the diagonal and first superdiagonal -* of A are overwritten by the corresponding elements of the -* tridiagonal matrix T, and the elements above the first -* superdiagonal, with the array TAU, represent the orthogonal -* matrix Q as a product of elementary reflectors; if UPLO -* = 'L', the diagonal and first subdiagonal of A are over- -* written by the corresponding elements of the tridiagonal -* matrix T, and the elements below the first subdiagonal, with -* the array TAU, represent the orthogonal matrix Q as a product -* of elementary reflectors. See Further Details. -* -* D (output) DOUBLE PRECISION array, dimension (N) -* The diagonal elements of the tridiagonal matrix T: -* D(i) = A(i,i). -* -* E (output) DOUBLE PRECISION array, dimension (N-1) -* The off-diagonal elements of the tridiagonal matrix T: -* E(i) = A(i,i+1) if UPLO = 'U', E(i) = A(i+1,i) if UPLO = 'L'. -* -* TAU (output) DOUBLE PRECISION array, dimension (N-1) -* The scalar factors of the elementary reflectors (see Further -* Details). -* -* INFO (output) INTEGER -* = 0: successful exit -* < 0: if INFO = -i, the i-th argument had an illegal value -* -* Further Details -* =============== -* -* If UPLO = 'U', the matrix Q is represented as a product of elementary -* reflectors -* -* Q = H(n-1) . . . H(2) H(1). -* -* Each H(i) has the form -* -* H(i) = I - tau * v * v' -* -* where tau is a real scalar, and v is a real vector with -* v(i+1:n) = 0 and v(i) = 1; v(1:i-1) is stored on exit in AP, -* overwriting A(1:i-1,i+1), and tau is stored in TAU(i). -* -* If UPLO = 'L', the matrix Q is represented as a product of elementary -* reflectors -* -* Q = H(1) H(2) . . . H(n-1). -* -* Each H(i) has the form -* -* H(i) = I - tau * v * v' -* -* where tau is a real scalar, and v is a real vector with -* v(1:i) = 0 and v(i+1) = 1; v(i+2:n) is stored on exit in AP, -* overwriting A(i+2:n,i), and tau is stored in TAU(i). -* -* ===================================================================== -* -* .. Parameters .. - DOUBLE PRECISION ONE, ZERO, HALF - PARAMETER ( ONE = 1.0D0, ZERO = 0.0D0, - $ HALF = 1.0D0 / 2.0D0 ) -* .. -* .. Local Scalars .. - LOGICAL UPPER - INTEGER I, I1, I1I1, II - DOUBLE PRECISION ALPHA, TAUI -* .. -* .. External Subroutines .. - EXTERNAL DAXPY, DLARFG, DSPMV, DSPR2, XERBLA -* .. -* .. External Functions .. - LOGICAL LSAME - DOUBLE PRECISION DDOT - EXTERNAL LSAME, DDOT -* .. -* .. Executable Statements .. -* -* Test the input parameters -* - INFO = 0 - UPPER = LSAME( UPLO, 'U' ) - IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN - INFO = -1 - ELSE IF( N.LT.0 ) THEN - INFO = -2 - END IF - IF( INFO.NE.0 ) THEN - CALL XERBLA( 'DSPTRD', -INFO ) - RETURN - END IF -* -* Quick return if possible -* - IF( N.LE.0 ) - $ RETURN -* - IF( UPPER ) THEN -* -* Reduce the upper triangle of A. -* I1 is the index in AP of A(1,I+1). -* - I1 = N*( N-1 ) / 2 + 1 - DO 10 I = N - 1, 1, -1 -* -* Generate elementary reflector H(i) = I - tau * v * v' -* to annihilate A(1:i-1,i+1) -* - CALL DLARFG( I, AP( I1+I-1 ), AP( I1 ), 1, TAUI ) - E( I ) = AP( I1+I-1 ) -* - IF( TAUI.NE.ZERO ) THEN -* -* Apply H(i) from both sides to A(1:i,1:i) -* - AP( I1+I-1 ) = ONE -* -* Compute y := tau * A * v storing y in TAU(1:i) -* - CALL DSPMV( UPLO, I, TAUI, AP, AP( I1 ), 1, ZERO, TAU, - $ 1 ) -* -* Compute w := y - 1/2 * tau * (y'*v) * v -* - ALPHA = -HALF*TAUI*DDOT( I, TAU, 1, AP( I1 ), 1 ) - CALL DAXPY( I, ALPHA, AP( I1 ), 1, TAU, 1 ) -* -* Apply the transformation as a rank-2 update: -* A := A - v * w' - w * v' -* - CALL DSPR2( UPLO, I, -ONE, AP( I1 ), 1, TAU, 1, AP ) -* - AP( I1+I-1 ) = E( I ) - END IF - D( I+1 ) = AP( I1+I ) - TAU( I ) = TAUI - I1 = I1 - I - 10 CONTINUE - D( 1 ) = AP( 1 ) - ELSE -* -* Reduce the lower triangle of A. II is the index in AP of -* A(i,i) and I1I1 is the index of A(i+1,i+1). -* - II = 1 - DO 20 I = 1, N - 1 - I1I1 = II + N - I + 1 -* -* Generate elementary reflector H(i) = I - tau * v * v' -* to annihilate A(i+2:n,i) -* - CALL DLARFG( N-I, AP( II+1 ), AP( II+2 ), 1, TAUI ) - E( I ) = AP( II+1 ) -* - IF( TAUI.NE.ZERO ) THEN -* -* Apply H(i) from both sides to A(i+1:n,i+1:n) -* - AP( II+1 ) = ONE -* -* Compute y := tau * A * v storing y in TAU(i:n-1) -* - CALL DSPMV( UPLO, N-I, TAUI, AP( I1I1 ), AP( II+1 ), 1, - $ ZERO, TAU( I ), 1 ) -* -* Compute w := y - 1/2 * tau * (y'*v) * v -* - ALPHA = -HALF*TAUI*DDOT( N-I, TAU( I ), 1, AP( II+1 ), - $ 1 ) - CALL DAXPY( N-I, ALPHA, AP( II+1 ), 1, TAU( I ), 1 ) -* -* Apply the transformation as a rank-2 update: -* A := A - v * w' - w * v' -* - CALL DSPR2( UPLO, N-I, -ONE, AP( II+1 ), 1, TAU( I ), 1, - $ AP( I1I1 ) ) -* - AP( II+1 ) = E( I ) - END IF - D( I ) = AP( II ) - TAU( I ) = TAUI - II = I1I1 - 20 CONTINUE - D( N ) = AP( II ) - END IF -* - RETURN -* -* End of DSPTRD -* - END |