octave-nkf: libcruft/lapack/slaqr5.f annotate

annotate libcruft/lapack/slaqr5.f @ 9409:ed392e13c117

remove assertion in force_variable

author	Jaroslav Hajek <highegg@gmail.com>
date	Tue, 30 Jun 2009 08:18:26 +0200
parents	82be108cc558
children

rev	line source
7789 82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	1 SUBROUTINE SLAQR5( WANTT, WANTZ, KACC22, N, KTOP, KBOT, NSHFTS,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	2 $ SR, SI, H, LDH, ILOZ, IHIZ, Z, LDZ, V, LDV, U,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	3 $ LDU, NV, WV, LDWV, NH, WH, LDWH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	4 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	5 * -- LAPACK auxiliary routine (version 3.1) --
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	6 * 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	7 * November 2006
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	8 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	9 * .. Scalar Arguments ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	10 INTEGER IHIZ, ILOZ, KACC22, KBOT, KTOP, LDH, LDU, LDV,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	11 $ LDWH, LDWV, LDZ, N, NH, NSHFTS, NV
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	12 LOGICAL WANTT, WANTZ
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	13 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	14 * .. Array Arguments ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	15 REAL H( LDH, * ), SI( * ), SR( * ), U( LDU, * ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	16 $ V( LDV, * ), WH( LDWH, * ), WV( LDWV, * ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	17 $ Z( LDZ, * )
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 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	20 * This auxiliary subroutine called by SLAQR0 performs a
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	21 * single small-bulge multi-shift QR sweep.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	22 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	23 * WANTT (input) logical scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	24 * WANTT = .true. if the quasi-triangular Schur factor
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	25 * is being computed. WANTT is set to .false. otherwise.
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 * WANTZ (input) logical scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	28 * WANTZ = .true. if the orthogonal Schur factor is being
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	29 * computed. WANTZ is set to .false. otherwise.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	30 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	31 * KACC22 (input) integer with value 0, 1, or 2.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	32 * Specifies the computation mode of far-from-diagonal
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	33 * orthogonal updates.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	34 * = 0: SLAQR5 does not accumulate reflections and does not
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	35 * use matrix-matrix multiply to update far-from-diagonal
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	36 * matrix entries.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	37 * = 1: SLAQR5 accumulates reflections and uses matrix-matrix
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	38 * multiply to update the far-from-diagonal matrix entries.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	39 * = 2: SLAQR5 accumulates reflections, uses matrix-matrix
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	40 * multiply to update the far-from-diagonal matrix entries,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	41 * and takes advantage of 2-by-2 block structure during
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	42 * matrix multiplies.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	43 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	44 * N (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	45 * N is the order of the Hessenberg matrix H upon which this
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	46 * subroutine operates.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	47 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	48 * KTOP (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	49 * KBOT (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	50 * These are the first and last rows and columns of an
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	51 * isolated diagonal block upon which the QR sweep is to be
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	52 * applied. It is assumed without a check that
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	53 * either KTOP = 1 or H(KTOP,KTOP-1) = 0
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	54 * and
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	55 * either KBOT = N or H(KBOT+1,KBOT) = 0.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	56 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	57 * NSHFTS (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	58 * NSHFTS gives the number of simultaneous shifts. NSHFTS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	59 * must be positive and even.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	60 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	61 * SR (input) REAL array of size (NSHFTS)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	62 * SI (input) REAL array of size (NSHFTS)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	63 * SR contains the real parts and SI contains the imaginary
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	64 * parts of the NSHFTS shifts of origin that define the
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	65 * multi-shift QR sweep.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	66 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	67 * H (input/output) REAL array of size (LDH,N)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	68 * On input H contains a Hessenberg matrix. On output a
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	69 * multi-shift QR sweep with shifts SR(J)+i*SI(J) is applied
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	70 * to the isolated diagonal block in rows and columns KTOP
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	71 * through KBOT.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	72 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	73 * LDH (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	74 * LDH is the leading dimension of H just as declared in the
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	75 * calling procedure. LDH.GE.MAX(1,N).
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	76 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	77 * ILOZ (input) INTEGER
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	78 * IHIZ (input) INTEGER
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	79 * Specify the rows of Z to which transformations must be
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	80 * applied if WANTZ is .TRUE.. 1 .LE. ILOZ .LE. IHIZ .LE. N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	81 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	82 * Z (input/output) REAL array of size (LDZ,IHI)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	83 * If WANTZ = .TRUE., then the QR Sweep orthogonal
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	84 * similarity transformation is accumulated into
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	85 * Z(ILOZ:IHIZ,ILO:IHI) from the right.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	86 * If WANTZ = .FALSE., then Z is unreferenced.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	87 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	88 * LDZ (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	89 * LDA is the leading dimension of Z just as declared in
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	90 * the calling procedure. LDZ.GE.N.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	91 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	92 * V (workspace) REAL array of size (LDV,NSHFTS/2)
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 * LDV (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	95 * LDV is the leading dimension of V as declared in the
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	96 * calling procedure. LDV.GE.3.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	97 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	98 * U (workspace) REAL array of size
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	99 * (LDU,3*NSHFTS-3)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	100 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	101 * LDU (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	102 * LDU is the leading dimension of U just as declared in the
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	103 * in the calling subroutine. LDU.GE.3*NSHFTS-3.
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 * NH (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	106 * NH is the number of columns in array WH available for
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	107 * workspace. NH.GE.1.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	108 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	109 * WH (workspace) REAL array of size (LDWH,NH)
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 * LDWH (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	112 * Leading dimension of WH just as declared in the
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	113 * calling procedure. LDWH.GE.3*NSHFTS-3.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	114 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	115 * NV (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	116 * NV is the number of rows in WV agailable for workspace.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	117 * NV.GE.1.
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 * WV (workspace) REAL array of size
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	120 * (LDWV,3*NSHFTS-3)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	121 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	122 * LDWV (input) integer scalar
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	123 * LDWV is the leading dimension of WV as declared in the
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	124 * in the calling subroutine. LDWV.GE.NV.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	125 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	126 *
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 * Based on contributions by
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	129 * Karen Braman and Ralph Byers, Department of Mathematics,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	130 * University of Kansas, USA
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	131 *
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 * Reference:
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	134 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	135 * K. Braman, R. Byers and R. Mathias, The Multi-Shift QR
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	136 * Algorithm Part I: Maintaining Well Focused Shifts, and
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	137 * Level 3 Performance, SIAM Journal of Matrix Analysis,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	138 * volume 23, pages 929--947, 2002.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	139 *
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 * .. Parameters ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	142 REAL ZERO, ONE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	143 PARAMETER ( ZERO = 0.0e0, ONE = 1.0e0 )
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 * .. Local Scalars ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	146 REAL ALPHA, BETA, H11, H12, H21, H22, REFSUM,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	147 $ SAFMAX, SAFMIN, SCL, SMLNUM, SWAP, TST1, TST2,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	148 $ ULP
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	149 INTEGER I, I2, I4, INCOL, J, J2, J4, JBOT, JCOL, JLEN,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	150 $ JROW, JTOP, K, K1, KDU, KMS, KNZ, KRCOL, KZS,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	151 $ M, M22, MBOT, MEND, MSTART, MTOP, NBMPS, NDCOL,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	152 $ NS, NU
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	153 LOGICAL ACCUM, BLK22, BMP22
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	154 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	155 * .. External Functions ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	156 REAL SLAMCH
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	157 EXTERNAL SLAMCH
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 * .. Intrinsic Functions ..
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 INTRINSIC ABS, MAX, MIN, MOD, REAL
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 * .. Local Arrays ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	164 REAL VT( 3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	165 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	166 * .. External Subroutines ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	167 EXTERNAL SGEMM, SLABAD, SLACPY, SLAQR1, SLARFG, SLASET,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	168 $ STRMM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	169 * ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	170 * .. Executable Statements ..
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	171 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	172 * ==== If there are no shifts, then there is nothing to do. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	173 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	174 IF( NSHFTS.LT.2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	175 $ RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	176 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	177 * ==== If the active block is empty or 1-by-1, then there
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	178 * . is nothing to do. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	179 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	180 IF( KTOP.GE.KBOT )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	181 $ RETURN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	182 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	183 * ==== Shuffle shifts into pairs of real shifts and pairs
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	184 * . of complex conjugate shifts assuming complex
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	185 * . conjugate shifts are already adjacent to one
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	186 * . another. ====
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 DO 10 I = 1, NSHFTS - 2, 2
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	189 IF( SI( I ).NE.-SI( I+1 ) ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	190 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	191 SWAP = SR( I )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	192 SR( I ) = SR( I+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	193 SR( I+1 ) = SR( I+2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	194 SR( I+2 ) = SWAP
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	195 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	196 SWAP = SI( I )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	197 SI( I ) = SI( I+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	198 SI( I+1 ) = SI( I+2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	199 SI( I+2 ) = SWAP
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	200 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	201 10 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	202 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	203 * ==== NSHFTS is supposed to be even, but if is odd,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	204 * . then simply reduce it by one. The shuffle above
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	205 * . ensures that the dropped shift is real and that
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	206 * . the remaining shifts are paired. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	207 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	208 NS = NSHFTS - MOD( NSHFTS, 2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	209 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	210 * ==== Machine constants for deflation ====
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 SAFMIN = SLAMCH( 'SAFE MINIMUM' )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	213 SAFMAX = ONE / SAFMIN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	214 CALL SLABAD( SAFMIN, SAFMAX )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	215 ULP = SLAMCH( 'PRECISION' )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	216 SMLNUM = SAFMIN*( REAL( N ) / ULP )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	217 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	218 * ==== Use accumulated reflections to update far-from-diagonal
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	219 * . entries ? ====
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 ACCUM = ( KACC22.EQ.1 ) .OR. ( KACC22.EQ.2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	222 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	223 * ==== If so, exploit the 2-by-2 block structure? ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	224 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	225 BLK22 = ( NS.GT.2 ) .AND. ( KACC22.EQ.2 )
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 * ==== clear trash ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	228 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	229 IF( KTOP+2.LE.KBOT )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	230 $ H( KTOP+2, KTOP ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	231 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	232 * ==== NBMPS = number of 2-shift bulges in the chain ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	233 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	234 NBMPS = NS / 2
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	235 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	236 * ==== KDU = width of slab ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	237 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	238 KDU = 6*NBMPS - 3
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	239 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	240 * ==== Create and chase chains of NBMPS bulges ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	241 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	242 DO 220 INCOL = 3( 1-NBMPS ) + KTOP - 1, KBOT - 2, 3NBMPS - 2
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	243 NDCOL = INCOL + KDU
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	244 IF( ACCUM )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	245 $ CALL SLASET( 'ALL', KDU, KDU, ZERO, ONE, U, LDU )
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 * ==== Near-the-diagonal bulge chase. The following loop
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	248 * . performs the near-the-diagonal part of a small bulge
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	249 * . multi-shift QR sweep. Each 6*NBMPS-2 column diagonal
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	250 * . chunk extends from column INCOL to column NDCOL
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	251 * . (including both column INCOL and column NDCOL). The
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	252 * . following loop chases a 3*NBMPS column long chain of
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	253 * . NBMPS bulges 3*NBMPS-2 columns to the right. (INCOL
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	254 * . may be less than KTOP and and NDCOL may be greater than
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	255 * . KBOT indicating phantom columns from which to chase
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	256 * . bulges before they are actually introduced or to which
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	257 * . to chase bulges beyond column KBOT.) ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	258 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	259 DO 150 KRCOL = INCOL, MIN( INCOL+3*NBMPS-3, KBOT-2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	260 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	261 * ==== Bulges number MTOP to MBOT are active double implicit
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	262 * . shift bulges. There may or may not also be small
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	263 * . 2-by-2 bulge, if there is room. The inactive bulges
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	264 * . (if any) must wait until the active bulges have moved
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	265 * . down the diagonal to make room. The phantom matrix
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	266 * . paradigm described above helps keep track. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	267 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	268 MTOP = MAX( 1, ( ( KTOP-1 )-KRCOL+2 ) / 3+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	269 MBOT = MIN( NBMPS, ( KBOT-KRCOL ) / 3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	270 M22 = MBOT + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	271 BMP22 = ( MBOT.LT.NBMPS ) .AND. ( KRCOL+3*( M22-1 ) ).EQ.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	272 $ ( KBOT-2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	273 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	274 * ==== Generate reflections to chase the chain right
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	275 * . one column. (The minimum value of K is KTOP-1.) ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	276 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	277 DO 20 M = MTOP, MBOT
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	278 K = KRCOL + 3*( M-1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	279 IF( K.EQ.KTOP-1 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	280 CALL SLAQR1( 3, H( KTOP, KTOP ), LDH, SR( 2*M-1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	281 $ SI( 2M-1 ), SR( 2M ), SI( 2*M ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	282 $ V( 1, M ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	283 ALPHA = V( 1, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	284 CALL SLARFG( 3, ALPHA, V( 2, M ), 1, V( 1, M ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	285 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	286 BETA = H( K+1, K )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	287 V( 2, M ) = H( K+2, K )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	288 V( 3, M ) = H( K+3, K )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	289 CALL SLARFG( 3, BETA, V( 2, M ), 1, V( 1, M ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	290 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	291 * ==== A Bulge may collapse because of vigilant
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	292 * . deflation or destructive underflow. (The
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	293 * . initial bulge is always collapsed.) Use
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	294 * . the two-small-subdiagonals trick to try
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	295 * . to get it started again. If V(2,M).NE.0 and
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	296 * . V(3,M) = H(K+3,K+1) = H(K+3,K+2) = 0, then
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	297 * . this bulge is collapsing into a zero
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	298 * . subdiagonal. It will be restarted next
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	299 * . trip through the loop.)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	300 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	301 IF( V( 1, M ).NE.ZERO .AND.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	302 $ ( V( 3, M ).NE.ZERO .OR. ( H( K+3,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	303 $ K+1 ).EQ.ZERO .AND. H( K+3, K+2 ).EQ.ZERO ) ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	304 $ THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	305 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	306 * ==== Typical case: not collapsed (yet). ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	307 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	308 H( K+1, K ) = BETA
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	309 H( K+2, K ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	310 H( K+3, K ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	311 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	312 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	313 * ==== Atypical case: collapsed. Attempt to
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	314 * . reintroduce ignoring H(K+1,K). If the
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	315 * . fill resulting from the new reflector
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	316 * . is too large, then abandon it.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	317 * . Otherwise, use the new one. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	318 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	319 CALL SLAQR1( 3, H( K+1, K+1 ), LDH, SR( 2*M-1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	320 $ SI( 2M-1 ), SR( 2M ), SI( 2*M ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	321 $ VT )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	322 SCL = ABS( VT( 1 ) ) + ABS( VT( 2 ) ) +
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	323 $ ABS( VT( 3 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	324 IF( SCL.NE.ZERO ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	325 VT( 1 ) = VT( 1 ) / SCL
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	326 VT( 2 ) = VT( 2 ) / SCL
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	327 VT( 3 ) = VT( 3 ) / SCL
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 * ==== The following is the traditional and
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	331 * . conservative two-small-subdiagonals
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	332 * . test. ====
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 IF( ABS( H( K+1, K ) )*( ABS( VT( 2 ) )+
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	335 $ ABS( VT( 3 ) ) ).GT.ULPABS( VT( 1 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	336 $ ( ABS( H( K, K ) )+ABS( H( K+1,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	337 $ K+1 ) )+ABS( H( K+2, K+2 ) ) ) ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	338 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	339 * ==== Starting a new bulge here would
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	340 * . create non-negligible fill. If
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	341 * . the old reflector is diagonal (only
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	342 * . possible with underflows), then
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	343 * . change it to I. Otherwise, use
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	344 * . it with trepidation. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	345 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	346 IF( V( 2, M ).EQ.ZERO .AND. V( 3, M ).EQ.ZERO )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	347 $ THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	348 V( 1, M ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	349 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	350 H( K+1, K ) = BETA
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	351 H( K+2, K ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	352 H( K+3, K ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	353 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	354 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	355 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	356 * ==== Stating a new bulge here would
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	357 * . create only negligible fill.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	358 * . Replace the old reflector with
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	359 * . the new one. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	360 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	361 ALPHA = VT( 1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	362 CALL SLARFG( 3, ALPHA, VT( 2 ), 1, VT( 1 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	363 REFSUM = H( K+1, K ) + H( K+2, K )*VT( 2 ) +
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	364 $ H( K+3, K )*VT( 3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	365 H( K+1, K ) = H( K+1, K ) - VT( 1 )*REFSUM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	366 H( K+2, K ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	367 H( K+3, K ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	368 V( 1, M ) = VT( 1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	369 V( 2, M ) = VT( 2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	370 V( 3, M ) = VT( 3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	371 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	372 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	373 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	374 20 CONTINUE
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 * ==== Generate a 2-by-2 reflection, if needed. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	377 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	378 K = KRCOL + 3*( M22-1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	379 IF( BMP22 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	380 IF( K.EQ.KTOP-1 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	381 CALL SLAQR1( 2, H( K+1, K+1 ), LDH, SR( 2*M22-1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	382 $ SI( 2M22-1 ), SR( 2M22 ), SI( 2*M22 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	383 $ V( 1, M22 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	384 BETA = V( 1, M22 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	385 CALL SLARFG( 2, BETA, V( 2, M22 ), 1, V( 1, M22 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	386 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	387 BETA = H( K+1, K )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	388 V( 2, M22 ) = H( K+2, K )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	389 CALL SLARFG( 2, BETA, V( 2, M22 ), 1, V( 1, M22 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	390 H( K+1, K ) = BETA
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	391 H( K+2, K ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	392 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	393 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	394 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	395 * ==== Initialize V(1,M22) here to avoid possible undefined
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	396 * . variable problems later. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	397 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	398 V( 1, M22 ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	399 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	400 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	401 * ==== Multiply H by reflections from the left ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	402 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	403 IF( ACCUM ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	404 JBOT = MIN( NDCOL, KBOT )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	405 ELSE IF( WANTT ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	406 JBOT = N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	407 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	408 JBOT = KBOT
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	409 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	410 DO 40 J = MAX( KTOP, KRCOL ), JBOT
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	411 MEND = MIN( MBOT, ( J-KRCOL+2 ) / 3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	412 DO 30 M = MTOP, MEND
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	413 K = KRCOL + 3*( M-1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	414 REFSUM = V( 1, M )( H( K+1, J )+V( 2, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	415 $ H( K+2, J )+V( 3, M )*H( K+3, J ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	416 H( K+1, J ) = H( K+1, J ) - REFSUM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	417 H( K+2, J ) = H( K+2, J ) - REFSUM*V( 2, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	418 H( K+3, J ) = H( K+3, J ) - REFSUM*V( 3, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	419 30 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	420 40 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	421 IF( BMP22 ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	422 K = KRCOL + 3*( M22-1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	423 DO 50 J = MAX( K+1, KTOP ), JBOT
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	424 REFSUM = V( 1, M22 )( H( K+1, J )+V( 2, M22 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	425 $ H( K+2, J ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	426 H( K+1, J ) = H( K+1, J ) - REFSUM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	427 H( K+2, J ) = H( K+2, J ) - REFSUM*V( 2, M22 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	428 50 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	429 END IF
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 * ==== Multiply H by reflections from the right.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	432 * . Delay filling in the last row until the
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	433 * . vigilant deflation check is complete. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	434 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	435 IF( ACCUM ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	436 JTOP = MAX( KTOP, INCOL )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	437 ELSE IF( WANTT ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	438 JTOP = 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	439 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	440 JTOP = KTOP
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	441 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	442 DO 90 M = MTOP, MBOT
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	443 IF( V( 1, M ).NE.ZERO ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	444 K = KRCOL + 3*( M-1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	445 DO 60 J = JTOP, MIN( KBOT, K+3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	446 REFSUM = V( 1, M )( H( J, K+1 )+V( 2, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	447 $ H( J, K+2 )+V( 3, M )*H( J, K+3 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	448 H( J, K+1 ) = H( J, K+1 ) - REFSUM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	449 H( J, K+2 ) = H( J, K+2 ) - REFSUM*V( 2, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	450 H( J, K+3 ) = H( J, K+3 ) - REFSUM*V( 3, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	451 60 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	452 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	453 IF( ACCUM ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	454 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	455 * ==== Accumulate U. (If necessary, update Z later
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	456 * . with with an efficient matrix-matrix
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	457 * . multiply.) ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	458 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	459 KMS = K - INCOL
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	460 DO 70 J = MAX( 1, KTOP-INCOL ), KDU
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	461 REFSUM = V( 1, M )( U( J, KMS+1 )+V( 2, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	462 $ U( J, KMS+2 )+V( 3, M )*U( J, KMS+3 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	463 U( J, KMS+1 ) = U( J, KMS+1 ) - REFSUM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	464 U( J, KMS+2 ) = U( J, KMS+2 ) - REFSUM*V( 2, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	465 U( J, KMS+3 ) = U( J, KMS+3 ) - REFSUM*V( 3, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	466 70 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	467 ELSE IF( WANTZ ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	468 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	469 * ==== U is not accumulated, so update Z
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	470 * . now by multiplying by reflections
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	471 * . from the right. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	472 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	473 DO 80 J = ILOZ, IHIZ
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	474 REFSUM = V( 1, M )( Z( J, K+1 )+V( 2, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	475 $ Z( J, K+2 )+V( 3, M )*Z( J, K+3 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	476 Z( J, K+1 ) = Z( J, K+1 ) - REFSUM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	477 Z( J, K+2 ) = Z( J, K+2 ) - REFSUM*V( 2, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	478 Z( J, K+3 ) = Z( J, K+3 ) - REFSUM*V( 3, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	479 80 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	480 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	481 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	482 90 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	483 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	484 * ==== Special case: 2-by-2 reflection (if needed) ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	485 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	486 K = KRCOL + 3*( M22-1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	487 IF( BMP22 .AND. ( V( 1, M22 ).NE.ZERO ) ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	488 DO 100 J = JTOP, MIN( KBOT, K+3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	489 REFSUM = V( 1, M22 )( H( J, K+1 )+V( 2, M22 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	490 $ H( J, K+2 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	491 H( J, K+1 ) = H( J, K+1 ) - REFSUM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	492 H( J, K+2 ) = H( J, K+2 ) - REFSUM*V( 2, M22 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	493 100 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	494 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	495 IF( ACCUM ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	496 KMS = K - INCOL
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	497 DO 110 J = MAX( 1, KTOP-INCOL ), KDU
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	498 REFSUM = V( 1, M22 )( U( J, KMS+1 )+V( 2, M22 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	499 $ U( J, KMS+2 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	500 U( J, KMS+1 ) = U( J, KMS+1 ) - REFSUM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	501 U( J, KMS+2 ) = U( J, KMS+2 ) - REFSUM*V( 2, M22 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	502 110 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	503 ELSE IF( WANTZ ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	504 DO 120 J = ILOZ, IHIZ
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	505 REFSUM = V( 1, M22 )( Z( J, K+1 )+V( 2, M22 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	506 $ Z( J, K+2 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	507 Z( J, K+1 ) = Z( J, K+1 ) - REFSUM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	508 Z( J, K+2 ) = Z( J, K+2 ) - REFSUM*V( 2, M22 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	509 120 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	510 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	511 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	512 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	513 * ==== Vigilant deflation check ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	514 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	515 MSTART = MTOP
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	516 IF( KRCOL+3*( MSTART-1 ).LT.KTOP )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	517 $ MSTART = MSTART + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	518 MEND = MBOT
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	519 IF( BMP22 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	520 $ MEND = MEND + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	521 IF( KRCOL.EQ.KBOT-2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	522 $ MEND = MEND + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	523 DO 130 M = MSTART, MEND
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	524 K = MIN( KBOT-1, KRCOL+3*( M-1 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	525 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	526 * ==== The following convergence test requires that
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	527 * . the tradition small-compared-to-nearby-diagonals
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	528 * . criterion and the Ahues & Tisseur (LAWN 122, 1997)
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	529 * . criteria both be satisfied. The latter improves
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	530 * . accuracy in some examples. Falling back on an
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	531 * . alternate convergence criterion when TST1 or TST2
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	532 * . is zero (as done here) is traditional but probably
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	533 * . unnecessary. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	534 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	535 IF( H( K+1, K ).NE.ZERO ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	536 TST1 = ABS( H( K, K ) ) + ABS( H( K+1, K+1 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	537 IF( TST1.EQ.ZERO ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	538 IF( K.GE.KTOP+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	539 $ TST1 = TST1 + ABS( H( K, K-1 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	540 IF( K.GE.KTOP+2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	541 $ TST1 = TST1 + ABS( H( K, K-2 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	542 IF( K.GE.KTOP+3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	543 $ TST1 = TST1 + ABS( H( K, K-3 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	544 IF( K.LE.KBOT-2 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	545 $ TST1 = TST1 + ABS( H( K+2, K+1 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	546 IF( K.LE.KBOT-3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	547 $ TST1 = TST1 + ABS( H( K+3, K+1 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	548 IF( K.LE.KBOT-4 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	549 $ TST1 = TST1 + ABS( H( K+4, K+1 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	550 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	551 IF( ABS( H( K+1, K ) ).LE.MAX( SMLNUM, ULP*TST1 ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	552 $ THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	553 H12 = MAX( ABS( H( K+1, K ) ), ABS( H( K, K+1 ) ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	554 H21 = MIN( ABS( H( K+1, K ) ), ABS( H( K, K+1 ) ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	555 H11 = MAX( ABS( H( K+1, K+1 ) ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	556 $ ABS( H( K, K )-H( K+1, K+1 ) ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	557 H22 = MIN( ABS( H( K+1, K+1 ) ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	558 $ ABS( H( K, K )-H( K+1, K+1 ) ) )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	559 SCL = H11 + H12
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	560 TST2 = H22*( H11 / SCL )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	561 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	562 IF( TST2.EQ.ZERO .OR. H21*( H12 / SCL ).LE.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	563 $ MAX( SMLNUM, ULP*TST2 ) )H( K+1, K ) = ZERO
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	564 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	565 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	566 130 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	567 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	568 * ==== Fill in the last row of each bulge. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	569 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	570 MEND = MIN( NBMPS, ( KBOT-KRCOL-1 ) / 3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	571 DO 140 M = MTOP, MEND
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	572 K = KRCOL + 3*( M-1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	573 REFSUM = V( 1, M )V( 3, M )H( K+4, K+3 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	574 H( K+4, K+1 ) = -REFSUM
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	575 H( K+4, K+2 ) = -REFSUM*V( 2, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	576 H( K+4, K+3 ) = H( K+4, K+3 ) - REFSUM*V( 3, M )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	577 140 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	578 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	579 * ==== End of near-the-diagonal bulge chase. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	580 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	581 150 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	582 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	583 * ==== Use U (if accumulated) to update far-from-diagonal
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	584 * . entries in H. If required, use U to update Z as
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	585 * . well. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	586 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	587 IF( ACCUM ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	588 IF( WANTT ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	589 JTOP = 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	590 JBOT = N
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	591 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	592 JTOP = KTOP
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	593 JBOT = KBOT
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	594 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	595 IF( ( .NOT.BLK22 ) .OR. ( INCOL.LT.KTOP ) .OR.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	596 $ ( NDCOL.GT.KBOT ) .OR. ( NS.LE.2 ) ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	597 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	598 * ==== Updates not exploiting the 2-by-2 block
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	599 * . structure of U. K1 and NU keep track of
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	600 * . the location and size of U in the special
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	601 * . cases of introducing bulges and chasing
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	602 * . bulges off the bottom. In these special
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	603 * . cases and in case the number of shifts
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	604 * . is NS = 2, there is no 2-by-2 block
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	605 * . structure to exploit. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	606 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	607 K1 = MAX( 1, KTOP-INCOL )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	608 NU = ( KDU-MAX( 0, NDCOL-KBOT ) ) - K1 + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	609 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	610 * ==== Horizontal Multiply ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	611 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	612 DO 160 JCOL = MIN( NDCOL, KBOT ) + 1, JBOT, NH
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	613 JLEN = MIN( NH, JBOT-JCOL+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	614 CALL SGEMM( 'C', 'N', NU, JLEN, NU, ONE, U( K1, K1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	615 $ LDU, H( INCOL+K1, JCOL ), LDH, ZERO, WH,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	616 $ LDWH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	617 CALL SLACPY( 'ALL', NU, JLEN, WH, LDWH,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	618 $ H( INCOL+K1, JCOL ), LDH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	619 160 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	620 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	621 * ==== Vertical multiply ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	622 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	623 DO 170 JROW = JTOP, MAX( KTOP, INCOL ) - 1, NV
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	624 JLEN = MIN( NV, MAX( KTOP, INCOL )-JROW )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	625 CALL SGEMM( 'N', 'N', JLEN, NU, NU, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	626 $ H( JROW, INCOL+K1 ), LDH, U( K1, K1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	627 $ LDU, ZERO, WV, LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	628 CALL SLACPY( 'ALL', JLEN, NU, WV, LDWV,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	629 $ H( JROW, INCOL+K1 ), LDH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	630 170 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	631 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	632 * ==== Z multiply (also vertical) ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	633 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	634 IF( WANTZ ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	635 DO 180 JROW = ILOZ, IHIZ, NV
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	636 JLEN = MIN( NV, IHIZ-JROW+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	637 CALL SGEMM( 'N', 'N', JLEN, NU, NU, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	638 $ Z( JROW, INCOL+K1 ), LDZ, U( K1, K1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	639 $ LDU, ZERO, WV, LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	640 CALL SLACPY( 'ALL', JLEN, NU, WV, LDWV,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	641 $ Z( JROW, INCOL+K1 ), LDZ )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	642 180 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	643 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	644 ELSE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	645 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	646 * ==== Updates exploiting U's 2-by-2 block structure.
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	647 * . (I2, I4, J2, J4 are the last rows and columns
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	648 * . of the blocks.) ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	649 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	650 I2 = ( KDU+1 ) / 2
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	651 I4 = KDU
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	652 J2 = I4 - I2
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	653 J4 = KDU
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	654 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	655 * ==== KZS and KNZ deal with the band of zeros
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	656 * . along the diagonal of one of the triangular
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	657 * . blocks. ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	658 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	659 KZS = ( J4-J2 ) - ( NS+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	660 KNZ = NS + 1
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	661 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	662 * ==== Horizontal multiply ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	663 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	664 DO 190 JCOL = MIN( NDCOL, KBOT ) + 1, JBOT, NH
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	665 JLEN = MIN( NH, JBOT-JCOL+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	666 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	667 * ==== Copy bottom of H to top+KZS of scratch ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	668 * (The first KZS rows get multiplied by zero.) ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	669 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	670 CALL SLACPY( 'ALL', KNZ, JLEN, H( INCOL+1+J2, JCOL ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	671 $ LDH, WH( KZS+1, 1 ), LDWH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	672 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	673 * ==== Multiply by U21' ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	674 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	675 CALL SLASET( 'ALL', KZS, JLEN, ZERO, ZERO, WH, LDWH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	676 CALL STRMM( 'L', 'U', 'C', 'N', KNZ, JLEN, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	677 $ U( J2+1, 1+KZS ), LDU, WH( KZS+1, 1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	678 $ LDWH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	679 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	680 * ==== Multiply top of H by U11' ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	681 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	682 CALL SGEMM( 'C', 'N', I2, JLEN, J2, ONE, U, LDU,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	683 $ H( INCOL+1, JCOL ), LDH, ONE, WH, LDWH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	684 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	685 * ==== Copy top of H bottom of WH ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	686 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	687 CALL SLACPY( 'ALL', J2, JLEN, H( INCOL+1, JCOL ), LDH,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	688 $ WH( I2+1, 1 ), LDWH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	689 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	690 * ==== Multiply by U21' ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	691 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	692 CALL STRMM( 'L', 'L', 'C', 'N', J2, JLEN, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	693 $ U( 1, I2+1 ), LDU, WH( I2+1, 1 ), LDWH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	694 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	695 * ==== Multiply by U22 ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	696 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	697 CALL SGEMM( 'C', 'N', I4-I2, JLEN, J4-J2, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	698 $ U( J2+1, I2+1 ), LDU,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	699 $ H( INCOL+1+J2, JCOL ), LDH, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	700 $ WH( I2+1, 1 ), LDWH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	701 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	702 * ==== Copy it back ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	703 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	704 CALL SLACPY( 'ALL', KDU, JLEN, WH, LDWH,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	705 $ H( INCOL+1, JCOL ), LDH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	706 190 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	707 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	708 * ==== Vertical multiply ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	709 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	710 DO 200 JROW = JTOP, MAX( INCOL, KTOP ) - 1, NV
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	711 JLEN = MIN( NV, MAX( INCOL, KTOP )-JROW )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	712 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	713 * ==== Copy right of H to scratch (the first KZS
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	714 * . columns get multiplied by zero) ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	715 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	716 CALL SLACPY( 'ALL', JLEN, KNZ, H( JROW, INCOL+1+J2 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	717 $ LDH, WV( 1, 1+KZS ), LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	718 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	719 * ==== Multiply by U21 ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	720 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	721 CALL SLASET( 'ALL', JLEN, KZS, ZERO, ZERO, WV, LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	722 CALL STRMM( 'R', 'U', 'N', 'N', JLEN, KNZ, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	723 $ U( J2+1, 1+KZS ), LDU, WV( 1, 1+KZS ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	724 $ LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	725 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	726 * ==== Multiply by U11 ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	727 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	728 CALL SGEMM( 'N', 'N', JLEN, I2, J2, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	729 $ H( JROW, INCOL+1 ), LDH, U, LDU, ONE, WV,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	730 $ LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	731 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	732 * ==== Copy left of H to right of scratch ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	733 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	734 CALL SLACPY( 'ALL', JLEN, J2, H( JROW, INCOL+1 ), LDH,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	735 $ WV( 1, 1+I2 ), LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	736 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	737 * ==== Multiply by U21 ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	738 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	739 CALL STRMM( 'R', 'L', 'N', 'N', JLEN, I4-I2, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	740 $ U( 1, I2+1 ), LDU, WV( 1, 1+I2 ), LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	741 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	742 * ==== Multiply by U22 ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	743 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	744 CALL SGEMM( 'N', 'N', JLEN, I4-I2, J4-J2, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	745 $ H( JROW, INCOL+1+J2 ), LDH,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	746 $ U( J2+1, I2+1 ), LDU, ONE, WV( 1, 1+I2 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	747 $ LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	748 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	749 * ==== Copy it back ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	750 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	751 CALL SLACPY( 'ALL', JLEN, KDU, WV, LDWV,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	752 $ H( JROW, INCOL+1 ), LDH )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	753 200 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	754 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	755 * ==== Multiply Z (also vertical) ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	756 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	757 IF( WANTZ ) THEN
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	758 DO 210 JROW = ILOZ, IHIZ, NV
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	759 JLEN = MIN( NV, IHIZ-JROW+1 )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	760 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	761 * ==== Copy right of Z to left of scratch (first
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	762 * . KZS columns get multiplied by zero) ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	763 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	764 CALL SLACPY( 'ALL', JLEN, KNZ,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	765 $ Z( JROW, INCOL+1+J2 ), LDZ,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	766 $ WV( 1, 1+KZS ), LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	767 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	768 * ==== Multiply by U12 ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	769 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	770 CALL SLASET( 'ALL', JLEN, KZS, ZERO, ZERO, WV,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	771 $ LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	772 CALL STRMM( 'R', 'U', 'N', 'N', JLEN, KNZ, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	773 $ U( J2+1, 1+KZS ), LDU, WV( 1, 1+KZS ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	774 $ LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	775 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	776 * ==== Multiply by U11 ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	777 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	778 CALL SGEMM( 'N', 'N', JLEN, I2, J2, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	779 $ Z( JROW, INCOL+1 ), LDZ, U, LDU, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	780 $ WV, LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	781 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	782 * ==== Copy left of Z to right of scratch ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	783 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	784 CALL SLACPY( 'ALL', JLEN, J2, Z( JROW, INCOL+1 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	785 $ LDZ, WV( 1, 1+I2 ), LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	786 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	787 * ==== Multiply by U21 ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	788 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	789 CALL STRMM( 'R', 'L', 'N', 'N', JLEN, I4-I2, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	790 $ U( 1, I2+1 ), LDU, WV( 1, 1+I2 ),
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	791 $ LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	792 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	793 * ==== Multiply by U22 ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	794 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	795 CALL SGEMM( 'N', 'N', JLEN, I4-I2, J4-J2, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	796 $ Z( JROW, INCOL+1+J2 ), LDZ,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	797 $ U( J2+1, I2+1 ), LDU, ONE,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	798 $ WV( 1, 1+I2 ), LDWV )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	799 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	800 * ==== Copy the result back to Z ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	801 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	802 CALL SLACPY( 'ALL', JLEN, KDU, WV, LDWV,
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	803 $ Z( JROW, INCOL+1 ), LDZ )
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	804 210 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	805 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	806 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	807 END IF
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	808 220 CONTINUE
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	809 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	810 * ==== End of SLAQR5 ====
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	811 *
82be108cc558 First attempt at single precision tyeps David Bateman <dbateman@free.fr> parents: diff changeset	812 END

Mercurial > hg > octave-nkf

annotate libcruft/lapack/slaqr5.f @ 9409:ed392e13c117