octave-lyh: libcruft/lapack/clalsd.f annotate

author	David Bateman <dbateman@free.fr>
date	Sun, 27 Apr 2008 22:34:17 +0200
parents
children

rev	line source
7789 82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	1 SUBROUTINE CLALSD( UPLO, SMLSIZ, N, NRHS, D, E, B, LDB, RCOND,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	2 $ RANK, WORK, RWORK, IWORK, INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	3 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	4 * -- LAPACK routine (version 3.1) --
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	5 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	6 * November 2006
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	7 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	8 * .. Scalar Arguments ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	9 CHARACTER UPLO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	10 INTEGER INFO, LDB, N, NRHS, RANK, SMLSIZ
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	11 REAL RCOND
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	12 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	13 * .. Array Arguments ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	14 INTEGER IWORK( * )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	15 REAL D( * ), E( * ), RWORK( * )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	16 COMPLEX B( LDB, * ), WORK( * )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	17 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	18 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	19 * Purpose
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	20 * =======
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	21 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	22 * CLALSD uses the singular value decomposition of A to solve the least
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	23 * squares problem of finding X to minimize the Euclidean norm of each
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	24 * column of A*X-B, where A is N-by-N upper bidiagonal, and X and B
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	25 * are N-by-NRHS. The solution X overwrites B.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	26 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	27 * The singular values of A smaller than RCOND times the largest
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	28 * singular value are treated as zero in solving the least squares
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	29 * problem; in this case a minimum norm solution is returned.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	30 * The actual singular values are returned in D in ascending order.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	31 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	32 * This code makes very mild assumptions about floating point
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	33 * arithmetic. It will work on machines with a guard digit in
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	34 * add/subtract, or on those binary machines without guard digits
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	35 * which subtract like the Cray XMP, Cray YMP, Cray C 90, or Cray 2.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	36 * It could conceivably fail on hexadecimal or decimal machines
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	37 * without guard digits, but we know of none.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	38 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	39 * Arguments
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	40 * =========
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	41 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	42 * UPLO (input) CHARACTER*1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	43 * = 'U': D and E define an upper bidiagonal matrix.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	44 * = 'L': D and E define a lower bidiagonal matrix.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	45 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	46 * SMLSIZ (input) INTEGER
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	47 * The maximum size of the subproblems at the bottom of the
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	48 * computation tree.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	49 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	50 * N (input) INTEGER
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	51 * The dimension of the bidiagonal matrix. N >= 0.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	52 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	53 * NRHS (input) INTEGER
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	54 * The number of columns of B. NRHS must be at least 1.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	55 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	56 * D (input/output) REAL array, dimension (N)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	57 * On entry D contains the main diagonal of the bidiagonal
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	58 * matrix. On exit, if INFO = 0, D contains its singular values.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	59 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	60 * E (input/output) REAL array, dimension (N-1)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	61 * Contains the super-diagonal entries of the bidiagonal matrix.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	62 * On exit, E has been destroyed.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	63 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	64 * B (input/output) COMPLEX array, dimension (LDB,NRHS)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	65 * On input, B contains the right hand sides of the least
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	66 * squares problem. On output, B contains the solution X.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	67 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	68 * LDB (input) INTEGER
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	69 * The leading dimension of B in the calling subprogram.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	70 * LDB must be at least max(1,N).
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	71 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	72 * RCOND (input) REAL
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	73 * The singular values of A less than or equal to RCOND times
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	74 * the largest singular value are treated as zero in solving
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	75 * the least squares problem. If RCOND is negative,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	76 * machine precision is used instead.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	77 * For example, if diag(S)*X=B were the least squares problem,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	78 * where diag(S) is a diagonal matrix of singular values, the
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	79 * solution would be X(i) = B(i) / S(i) if S(i) is greater than
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	80 * RCOND*max(S), and X(i) = 0 if S(i) is less than or equal to
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	81 * RCOND*max(S).
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	82 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	83 * RANK (output) INTEGER
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	84 * The number of singular values of A greater than RCOND times
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	85 * the largest singular value.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	86 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	87 * WORK (workspace) COMPLEX array, dimension (N * NRHS).
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	88 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	89 * RWORK (workspace) REAL array, dimension at least
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	90 * (9N + 2NSMLSIZ + 8NNLVL + 3SMLSIZNRHS + (SMLSIZ+1)*2),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	91 * where
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	92 * NLVL = MAX( 0, INT( LOG_2( MIN( M,N )/(SMLSIZ+1) ) ) + 1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	93 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	94 * IWORK (workspace) INTEGER array, dimension (3NNLVL + 11*N).
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	95 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	96 * INFO (output) INTEGER
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	97 * = 0: successful exit.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	98 * < 0: if INFO = -i, the i-th argument had an illegal value.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	99 * > 0: The algorithm failed to compute an singular value while
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	100 * working on the submatrix lying in rows and columns
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	101 * INFO/(N+1) through MOD(INFO,N+1).
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	102 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	103 * Further Details
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	104 * ===============
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	105 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	106 * Based on contributions by
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	107 * Ming Gu and Ren-Cang Li, Computer Science Division, University of
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	108 * California at Berkeley, USA
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	109 * Osni Marques, LBNL/NERSC, USA
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	110 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	111 * =====================================================================
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	112 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	113 * .. Parameters ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	114 REAL ZERO, ONE, TWO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	115 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	116 COMPLEX CZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	117 PARAMETER ( CZERO = ( 0.0E0, 0.0E0 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	118 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	119 * .. Local Scalars ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	120 INTEGER BX, BXST, C, DIFL, DIFR, GIVCOL, GIVNUM,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	121 $ GIVPTR, I, ICMPQ1, ICMPQ2, IRWB, IRWIB, IRWRB,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	122 $ IRWU, IRWVT, IRWWRK, IWK, J, JCOL, JIMAG,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	123 $ JREAL, JROW, K, NLVL, NM1, NRWORK, NSIZE, NSUB,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	124 $ PERM, POLES, S, SIZEI, SMLSZP, SQRE, ST, ST1,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	125 $ U, VT, Z
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	126 REAL CS, EPS, ORGNRM, R, RCND, SN, TOL
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	127 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	128 * .. External Functions ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	129 INTEGER ISAMAX
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	130 REAL SLAMCH, SLANST
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	131 EXTERNAL ISAMAX, SLAMCH, SLANST
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	132 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	133 * .. External Subroutines ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	134 EXTERNAL CCOPY, CLACPY, CLALSA, CLASCL, CLASET, CSROT,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	135 $ SGEMM, SLARTG, SLASCL, SLASDA, SLASDQ, SLASET,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	136 $ SLASRT, XERBLA
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	137 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	138 * .. Intrinsic Functions ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	139 INTRINSIC ABS, AIMAG, CMPLX, INT, LOG, REAL, SIGN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	140 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	141 * .. Executable Statements ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	142 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	143 * Test the input parameters.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	144 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	145 INFO = 0
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	146 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	147 IF( N.LT.0 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	148 INFO = -3
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	149 ELSE IF( NRHS.LT.1 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	150 INFO = -4
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	151 ELSE IF( ( LDB.LT.1 ) .OR. ( LDB.LT.N ) ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	152 INFO = -8
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	153 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	154 IF( INFO.NE.0 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	155 CALL XERBLA( 'CLALSD', -INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	156 RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	157 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	158 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	159 EPS = SLAMCH( 'Epsilon' )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	160 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	161 * Set up the tolerance.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	162 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	163 IF( ( RCOND.LE.ZERO ) .OR. ( RCOND.GE.ONE ) ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	164 RCND = EPS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	165 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	166 RCND = RCOND
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	167 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	168 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	169 RANK = 0
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	170 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	171 * Quick return if possible.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	172 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	173 IF( N.EQ.0 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	174 RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	175 ELSE IF( N.EQ.1 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	176 IF( D( 1 ).EQ.ZERO ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	177 CALL CLASET( 'A', 1, NRHS, CZERO, CZERO, B, LDB )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	178 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	179 RANK = 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	180 CALL CLASCL( 'G', 0, 0, D( 1 ), ONE, 1, NRHS, B, LDB, INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	181 D( 1 ) = ABS( D( 1 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	182 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	183 RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	184 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	185 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	186 * Rotate the matrix if it is lower bidiagonal.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	187 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	188 IF( UPLO.EQ.'L' ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	189 DO 10 I = 1, N - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	190 CALL SLARTG( D( I ), E( I ), CS, SN, R )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	191 D( I ) = R
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	192 E( I ) = SN*D( I+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	193 D( I+1 ) = CS*D( I+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	194 IF( NRHS.EQ.1 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	195 CALL CSROT( 1, B( I, 1 ), 1, B( I+1, 1 ), 1, CS, SN )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	196 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	197 RWORK( I*2-1 ) = CS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	198 RWORK( I*2 ) = SN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	199 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	200 10 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	201 IF( NRHS.GT.1 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	202 DO 30 I = 1, NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	203 DO 20 J = 1, N - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	204 CS = RWORK( J*2-1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	205 SN = RWORK( J*2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	206 CALL CSROT( 1, B( J, I ), 1, B( J+1, I ), 1, CS, SN )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	207 20 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	208 30 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	209 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	210 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	211 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	212 * Scale.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	213 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	214 NM1 = N - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	215 ORGNRM = SLANST( 'M', N, D, E )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	216 IF( ORGNRM.EQ.ZERO ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	217 CALL CLASET( 'A', N, NRHS, CZERO, CZERO, B, LDB )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	218 RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	219 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	220 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	221 CALL SLASCL( 'G', 0, 0, ORGNRM, ONE, N, 1, D, N, INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	222 CALL SLASCL( 'G', 0, 0, ORGNRM, ONE, NM1, 1, E, NM1, INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	223 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	224 * If N is smaller than the minimum divide size SMLSIZ, then solve
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	225 * the problem with another solver.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	226 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	227 IF( N.LE.SMLSIZ ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	228 IRWU = 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	229 IRWVT = IRWU + N*N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	230 IRWWRK = IRWVT + N*N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	231 IRWRB = IRWWRK
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	232 IRWIB = IRWRB + N*NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	233 IRWB = IRWIB + N*NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	234 CALL SLASET( 'A', N, N, ZERO, ONE, RWORK( IRWU ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	235 CALL SLASET( 'A', N, N, ZERO, ONE, RWORK( IRWVT ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	236 CALL SLASDQ( 'U', 0, N, N, N, 0, D, E, RWORK( IRWVT ), N,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	237 $ RWORK( IRWU ), N, RWORK( IRWWRK ), 1,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	238 $ RWORK( IRWWRK ), INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	239 IF( INFO.NE.0 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	240 RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	241 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	242 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	243 * In the real version, B is passed to SLASDQ and multiplied
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	244 * internally by Q'. Here B is complex and that product is
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	245 * computed below in two steps (real and imaginary parts).
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	246 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	247 J = IRWB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	248 DO 50 JCOL = 1, NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	249 DO 40 JROW = 1, N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	250 J = J + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	251 RWORK( J ) = REAL( B( JROW, JCOL ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	252 40 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	253 50 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	254 CALL SGEMM( 'T', 'N', N, NRHS, N, ONE, RWORK( IRWU ), N,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	255 $ RWORK( IRWB ), N, ZERO, RWORK( IRWRB ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	256 J = IRWB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	257 DO 70 JCOL = 1, NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	258 DO 60 JROW = 1, N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	259 J = J + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	260 RWORK( J ) = AIMAG( B( JROW, JCOL ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	261 60 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	262 70 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	263 CALL SGEMM( 'T', 'N', N, NRHS, N, ONE, RWORK( IRWU ), N,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	264 $ RWORK( IRWB ), N, ZERO, RWORK( IRWIB ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	265 JREAL = IRWRB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	266 JIMAG = IRWIB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	267 DO 90 JCOL = 1, NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	268 DO 80 JROW = 1, N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	269 JREAL = JREAL + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	270 JIMAG = JIMAG + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	271 B( JROW, JCOL ) = CMPLX( RWORK( JREAL ), RWORK( JIMAG ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	272 80 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	273 90 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	274 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	275 TOL = RCND*ABS( D( ISAMAX( N, D, 1 ) ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	276 DO 100 I = 1, N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	277 IF( D( I ).LE.TOL ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	278 CALL CLASET( 'A', 1, NRHS, CZERO, CZERO, B( I, 1 ), LDB )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	279 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	280 CALL CLASCL( 'G', 0, 0, D( I ), ONE, 1, NRHS, B( I, 1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	281 $ LDB, INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	282 RANK = RANK + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	283 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	284 100 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	285 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	286 * Since B is complex, the following call to SGEMM is performed
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	287 * in two steps (real and imaginary parts). That is for V * B
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	288 * (in the real version of the code V' is stored in WORK).
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	289 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	290 * CALL SGEMM( 'T', 'N', N, NRHS, N, ONE, WORK, N, B, LDB, ZERO,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	291 * $ WORK( NWORK ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	292 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	293 J = IRWB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	294 DO 120 JCOL = 1, NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	295 DO 110 JROW = 1, N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	296 J = J + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	297 RWORK( J ) = REAL( B( JROW, JCOL ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	298 110 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	299 120 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	300 CALL SGEMM( 'T', 'N', N, NRHS, N, ONE, RWORK( IRWVT ), N,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	301 $ RWORK( IRWB ), N, ZERO, RWORK( IRWRB ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	302 J = IRWB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	303 DO 140 JCOL = 1, NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	304 DO 130 JROW = 1, N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	305 J = J + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	306 RWORK( J ) = AIMAG( B( JROW, JCOL ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	307 130 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	308 140 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	309 CALL SGEMM( 'T', 'N', N, NRHS, N, ONE, RWORK( IRWVT ), N,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	310 $ RWORK( IRWB ), N, ZERO, RWORK( IRWIB ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	311 JREAL = IRWRB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	312 JIMAG = IRWIB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	313 DO 160 JCOL = 1, NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	314 DO 150 JROW = 1, N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	315 JREAL = JREAL + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	316 JIMAG = JIMAG + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	317 B( JROW, JCOL ) = CMPLX( RWORK( JREAL ), RWORK( JIMAG ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	318 150 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	319 160 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	320 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	321 * Unscale.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	322 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	323 CALL SLASCL( 'G', 0, 0, ONE, ORGNRM, N, 1, D, N, INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	324 CALL SLASRT( 'D', N, D, INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	325 CALL CLASCL( 'G', 0, 0, ORGNRM, ONE, N, NRHS, B, LDB, INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	326 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	327 RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	328 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	329 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	330 * Book-keeping and setting up some constants.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	331 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	332 NLVL = INT( LOG( REAL( N ) / REAL( SMLSIZ+1 ) ) / LOG( TWO ) ) + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	333 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	334 SMLSZP = SMLSIZ + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	335 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	336 U = 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	337 VT = 1 + SMLSIZ*N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	338 DIFL = VT + SMLSZP*N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	339 DIFR = DIFL + NLVL*N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	340 Z = DIFR + NLVLN2
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	341 C = Z + NLVL*N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	342 S = C + N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	343 POLES = S + N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	344 GIVNUM = POLES + 2NLVLN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	345 NRWORK = GIVNUM + 2NLVLN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	346 BX = 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	347 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	348 IRWRB = NRWORK
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	349 IRWIB = IRWRB + SMLSIZ*NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	350 IRWB = IRWIB + SMLSIZ*NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	351 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	352 SIZEI = 1 + N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	353 K = SIZEI + N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	354 GIVPTR = K + N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	355 PERM = GIVPTR + N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	356 GIVCOL = PERM + NLVL*N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	357 IWK = GIVCOL + NLVLN2
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	358 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	359 ST = 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	360 SQRE = 0
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	361 ICMPQ1 = 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	362 ICMPQ2 = 0
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	363 NSUB = 0
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	364 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	365 DO 170 I = 1, N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	366 IF( ABS( D( I ) ).LT.EPS ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	367 D( I ) = SIGN( EPS, D( I ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	368 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	369 170 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	370 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	371 DO 240 I = 1, NM1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	372 IF( ( ABS( E( I ) ).LT.EPS ) .OR. ( I.EQ.NM1 ) ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	373 NSUB = NSUB + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	374 IWORK( NSUB ) = ST
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	375 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	376 * Subproblem found. First determine its size and then
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	377 * apply divide and conquer on it.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	378 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	379 IF( I.LT.NM1 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	380 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	381 * A subproblem with E(I) small for I < NM1.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	382 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	383 NSIZE = I - ST + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	384 IWORK( SIZEI+NSUB-1 ) = NSIZE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	385 ELSE IF( ABS( E( I ) ).GE.EPS ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	386 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	387 * A subproblem with E(NM1) not too small but I = NM1.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	388 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	389 NSIZE = N - ST + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	390 IWORK( SIZEI+NSUB-1 ) = NSIZE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	391 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	392 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	393 * A subproblem with E(NM1) small. This implies an
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	394 * 1-by-1 subproblem at D(N), which is not solved
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	395 * explicitly.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	396 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	397 NSIZE = I - ST + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	398 IWORK( SIZEI+NSUB-1 ) = NSIZE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	399 NSUB = NSUB + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	400 IWORK( NSUB ) = N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	401 IWORK( SIZEI+NSUB-1 ) = 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	402 CALL CCOPY( NRHS, B( N, 1 ), LDB, WORK( BX+NM1 ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	403 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	404 ST1 = ST - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	405 IF( NSIZE.EQ.1 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	406 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	407 * This is a 1-by-1 subproblem and is not solved
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	408 * explicitly.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	409 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	410 CALL CCOPY( NRHS, B( ST, 1 ), LDB, WORK( BX+ST1 ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	411 ELSE IF( NSIZE.LE.SMLSIZ ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	412 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	413 * This is a small subproblem and is solved by SLASDQ.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	414 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	415 CALL SLASET( 'A', NSIZE, NSIZE, ZERO, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	416 $ RWORK( VT+ST1 ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	417 CALL SLASET( 'A', NSIZE, NSIZE, ZERO, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	418 $ RWORK( U+ST1 ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	419 CALL SLASDQ( 'U', 0, NSIZE, NSIZE, NSIZE, 0, D( ST ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	420 $ E( ST ), RWORK( VT+ST1 ), N, RWORK( U+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	421 $ N, RWORK( NRWORK ), 1, RWORK( NRWORK ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	422 $ INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	423 IF( INFO.NE.0 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	424 RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	425 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	426 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	427 * In the real version, B is passed to SLASDQ and multiplied
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	428 * internally by Q'. Here B is complex and that product is
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	429 * computed below in two steps (real and imaginary parts).
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	430 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	431 J = IRWB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	432 DO 190 JCOL = 1, NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	433 DO 180 JROW = ST, ST + NSIZE - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	434 J = J + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	435 RWORK( J ) = REAL( B( JROW, JCOL ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	436 180 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	437 190 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	438 CALL SGEMM( 'T', 'N', NSIZE, NRHS, NSIZE, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	439 $ RWORK( U+ST1 ), N, RWORK( IRWB ), NSIZE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	440 $ ZERO, RWORK( IRWRB ), NSIZE )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	441 J = IRWB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	442 DO 210 JCOL = 1, NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	443 DO 200 JROW = ST, ST + NSIZE - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	444 J = J + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	445 RWORK( J ) = AIMAG( B( JROW, JCOL ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	446 200 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	447 210 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	448 CALL SGEMM( 'T', 'N', NSIZE, NRHS, NSIZE, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	449 $ RWORK( U+ST1 ), N, RWORK( IRWB ), NSIZE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	450 $ ZERO, RWORK( IRWIB ), NSIZE )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	451 JREAL = IRWRB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	452 JIMAG = IRWIB - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	453 DO 230 JCOL = 1, NRHS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	454 DO 220 JROW = ST, ST + NSIZE - 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	455 JREAL = JREAL + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	456 JIMAG = JIMAG + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	457 B( JROW, JCOL ) = CMPLX( RWORK( JREAL ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	458 $ RWORK( JIMAG ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	459 220 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	460 230 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	461 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	462 CALL CLACPY( 'A', NSIZE, NRHS, B( ST, 1 ), LDB,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	463 $ WORK( BX+ST1 ), N )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	464 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	465 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	466 * A large problem. Solve it using divide and conquer.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	467 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	468 CALL SLASDA( ICMPQ1, SMLSIZ, NSIZE, SQRE, D( ST ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	469 $ E( ST ), RWORK( U+ST1 ), N, RWORK( VT+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	470 $ IWORK( K+ST1 ), RWORK( DIFL+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	471 $ RWORK( DIFR+ST1 ), RWORK( Z+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	472 $ RWORK( POLES+ST1 ), IWORK( GIVPTR+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	473 $ IWORK( GIVCOL+ST1 ), N, IWORK( PERM+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	474 $ RWORK( GIVNUM+ST1 ), RWORK( C+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	475 $ RWORK( S+ST1 ), RWORK( NRWORK ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	476 $ IWORK( IWK ), INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	477 IF( INFO.NE.0 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	478 RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	479 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	480 BXST = BX + ST1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	481 CALL CLALSA( ICMPQ2, SMLSIZ, NSIZE, NRHS, B( ST, 1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	482 $ LDB, WORK( BXST ), N, RWORK( U+ST1 ), N,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	483 $ RWORK( VT+ST1 ), IWORK( K+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	484 $ RWORK( DIFL+ST1 ), RWORK( DIFR+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	485 $ RWORK( Z+ST1 ), RWORK( POLES+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	486 $ IWORK( GIVPTR+ST1 ), IWORK( GIVCOL+ST1 ), N,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	487 $ IWORK( PERM+ST1 ), RWORK( GIVNUM+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	488 $ RWORK( C+ST1 ), RWORK( S+ST1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	489 $ RWORK( NRWORK ), IWORK( IWK ), INFO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	490 IF( INFO.NE.0 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	491 RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	492 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	493 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	494 ST = I + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	495 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	496 240 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	497 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	498 * Apply the singular values and treat the tiny ones as zero.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	499 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	500 TOL = RCND*ABS( D( ISAMAX( N, D, 1 ) ) )

7789

82be108cc558 First attempt at single precision tyeps