Mercurial > octave-nkf

comparison libcruft/lapack/slaqr1.f @ 7789:82be108cc558

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
First attempt at single precision tyeps * * * corrections to qrupdate single precision routines * * * prefer demotion to single over promotion to double * * * Add single precision support to log2 function * * * Trivial PROJECT file update * * * Cache optimized hermitian/transpose methods * * * Add tests for tranpose/hermitian and ChangeLog entry for new transpose code
author David Bateman <dbateman@free.fr>
date Sun, 27 Apr 2008 22:34:17 +0200
parents
children
comparison
equal deleted inserted replaced
7788:45f5faba05a2 7789:82be108cc558
1 SUBROUTINE SLAQR1( N, H, LDH, SR1, SI1, SR2, SI2, V )
2 *
3 * -- LAPACK auxiliary routine (version 3.1) --
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 * November 2006
6 *
7 * .. Scalar Arguments ..
8 REAL SI1, SI2, SR1, SR2
9 INTEGER LDH, N
10 * ..
11 * .. Array Arguments ..
12 REAL H( LDH, * ), V( * )
13 * ..
14 *
15 * Given a 2-by-2 or 3-by-3 matrix H, SLAQR1 sets v to a
16 * scalar multiple of the first column of the product
17 *
18 * (*) K = (H - (sr1 + i*si1)*I)*(H - (sr2 + i*si2)*I)
19 *
20 * scaling to avoid overflows and most underflows. It
21 * is assumed that either
22 *
23 * 1) sr1 = sr2 and si1 = -si2
24 * or
25 * 2) si1 = si2 = 0.
26 *
27 * This is useful for starting double implicit shift bulges
28 * in the QR algorithm.
29 *
30 *
31 * N (input) integer
32 * Order of the matrix H. N must be either 2 or 3.
33 *
34 * H (input) REAL array of dimension (LDH,N)
35 * The 2-by-2 or 3-by-3 matrix H in (*).
36 *
37 * LDH (input) integer
38 * The leading dimension of H as declared in
39 * the calling procedure. LDH.GE.N
40 *
41 * SR1 (input) REAL
42 * SI1 The shifts in (*).
43 * SR2
44 * SI2
45 *
46 * V (output) REAL array of dimension N
47 * A scalar multiple of the first column of the
48 * matrix K in (*).
49 *
50 * ================================================================
51 * Based on contributions by
52 * Karen Braman and Ralph Byers, Department of Mathematics,
53 * University of Kansas, USA
54 *
55 * ================================================================
56 *
57 * .. Parameters ..
58 REAL ZERO
59 PARAMETER ( ZERO = 0.0e0 )
60 * ..
61 * .. Local Scalars ..
62 REAL H21S, H31S, S
63 * ..
64 * .. Intrinsic Functions ..
65 INTRINSIC ABS
66 * ..
67 * .. Executable Statements ..
68 IF( N.EQ.2 ) THEN
69 S = ABS( H( 1, 1 )-SR2 ) + ABS( SI2 ) + ABS( H( 2, 1 ) )
70 IF( S.EQ.ZERO ) THEN
71 V( 1 ) = ZERO
72 V( 2 ) = ZERO
73 ELSE
74 H21S = H( 2, 1 ) / S
75 V( 1 ) = H21S*H( 1, 2 ) + ( H( 1, 1 )-SR1 )*
76 $ ( ( H( 1, 1 )-SR2 ) / S ) - SI1*( SI2 / S )
77 V( 2 ) = H21S*( H( 1, 1 )+H( 2, 2 )-SR1-SR2 )
78 END IF
79 ELSE
80 S = ABS( H( 1, 1 )-SR2 ) + ABS( SI2 ) + ABS( H( 2, 1 ) ) +
81 $ ABS( H( 3, 1 ) )
82 IF( S.EQ.ZERO ) THEN
83 V( 1 ) = ZERO
84 V( 2 ) = ZERO
85 V( 3 ) = ZERO
86 ELSE
87 H21S = H( 2, 1 ) / S
88 H31S = H( 3, 1 ) / S
89 V( 1 ) = ( H( 1, 1 )-SR1 )*( ( H( 1, 1 )-SR2 ) / S ) -
90 $ SI1*( SI2 / S ) + H( 1, 2 )*H21S + H( 1, 3 )*H31S
91 V( 2 ) = H21S*( H( 1, 1 )+H( 2, 2 )-SR1-SR2 ) +
92 $ H( 2, 3 )*H31S
93 V( 3 ) = H31S*( H( 1, 1 )+H( 3, 3 )-SR1-SR2 ) +
94 $ H21S*H( 3, 2 )
95 END IF
96 END IF
97 END