LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dlartgs.f
Go to the documentation of this file.
1*> \brief \b DLARTGS generates a plane rotation designed to introduce a bulge in implicit QR iteration for the bidiagonal SVD problem.
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download DLARTGS + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlartgs.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlartgs.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlartgs.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE DLARTGS( X, Y, SIGMA, CS, SN )
22*
23* .. Scalar Arguments ..
24* DOUBLE PRECISION CS, SIGMA, SN, X, Y
25* ..
26*
27*
28*> \par Purpose:
29* =============
30*>
31*> \verbatim
32*>
33*> DLARTGS generates a plane rotation designed to introduce a bulge in
34*> Golub-Reinsch-style implicit QR iteration for the bidiagonal SVD
35*> problem. X and Y are the top-row entries, and SIGMA is the shift.
36*> The computed CS and SN define a plane rotation satisfying
37*>
38*> [ CS SN ] . [ X^2 - SIGMA ] = [ R ],
39*> [ -SN CS ] [ X * Y ] [ 0 ]
40*>
41*> with R nonnegative. If X^2 - SIGMA and X * Y are 0, then the
42*> rotation is by PI/2.
43*> \endverbatim
44*
45* Arguments:
46* ==========
47*
48*> \param[in] X
49*> \verbatim
50*> X is DOUBLE PRECISION
51*> The (1,1) entry of an upper bidiagonal matrix.
52*> \endverbatim
53*>
54*> \param[in] Y
55*> \verbatim
56*> Y is DOUBLE PRECISION
57*> The (1,2) entry of an upper bidiagonal matrix.
58*> \endverbatim
59*>
60*> \param[in] SIGMA
61*> \verbatim
62*> SIGMA is DOUBLE PRECISION
63*> The shift.
64*> \endverbatim
65*>
66*> \param[out] CS
67*> \verbatim
68*> CS is DOUBLE PRECISION
69*> The cosine of the rotation.
70*> \endverbatim
71*>
72*> \param[out] SN
73*> \verbatim
74*> SN is DOUBLE PRECISION
75*> The sine of the rotation.
76*> \endverbatim
77*
78* Authors:
79* ========
80*
81*> \author Univ. of Tennessee
82*> \author Univ. of California Berkeley
83*> \author Univ. of Colorado Denver
84*> \author NAG Ltd.
85*
86*> \ingroup lartgs
87*
88* =====================================================================
89 SUBROUTINE dlartgs( X, Y, SIGMA, CS, SN )
90*
91* -- LAPACK computational routine --
92* -- LAPACK is a software package provided by Univ. of Tennessee, --
93* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
94*
95* .. Scalar Arguments ..
96 DOUBLE PRECISION CS, SIGMA, SN, X, Y
97* ..
98*
99* ===================================================================
100*
101* .. Parameters ..
102 DOUBLE PRECISION NEGONE, ONE, ZERO
103 parameter( negone = -1.0d0, one = 1.0d0, zero = 0.0d0 )
104* ..
105* .. Local Scalars ..
106 DOUBLE PRECISION R, S, THRESH, W, Z
107* ..
108* .. External Subroutines ..
109 EXTERNAL dlartgp
110* ..
111* .. External Functions ..
112 DOUBLE PRECISION DLAMCH
113 EXTERNAL dlamch
114* .. Executable Statements ..
115*
116 thresh = dlamch('E')
117*
118* Compute the first column of B**T*B - SIGMA^2*I, up to a scale
119* factor.
120*
121 IF( (sigma .EQ. zero .AND. abs(x) .LT. thresh) .OR.
122 $ (abs(x) .EQ. sigma .AND. y .EQ. zero) ) THEN
123 z = zero
124 w = zero
125 ELSE IF( sigma .EQ. zero ) THEN
126 IF( x .GE. zero ) THEN
127 z = x
128 w = y
129 ELSE
130 z = -x
131 w = -y
132 END IF
133 ELSE IF( abs(x) .LT. thresh ) THEN
134 z = -sigma*sigma
135 w = zero
136 ELSE
137 IF( x .GE. zero ) THEN
138 s = one
139 ELSE
140 s = negone
141 END IF
142 z = s * (abs(x)-sigma) * (s+sigma/x)
143 w = s * y
144 END IF
145*
146* Generate the rotation.
147* CALL DLARTGP( Z, W, CS, SN, R ) might seem more natural;
148* reordering the arguments ensures that if Z = 0 then the rotation
149* is by PI/2.
150*
151 CALL dlartgp( w, z, sn, cs, r )
152*
153 RETURN
154*
155* End DLARTGS
156*
157 END
158
subroutine dlartgp(f, g, cs, sn, r)
DLARTGP generates a plane rotation so that the diagonal is nonnegative.
Definition dlartgp.f:95
subroutine dlartgs(x, y, sigma, cs, sn)
DLARTGS generates a plane rotation designed to introduce a bulge in implicit QR iteration for the bid...
Definition dlartgs.f:90