1*> \brief \b DROTMG
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*  Definition:
9*  ===========
10*
11*       SUBROUTINE DROTMG(DD1,DD2,DX1,DY1,DPARAM)
12*
13*       .. Scalar Arguments ..
14*       DOUBLE PRECISION DD1,DD2,DX1,DY1
15*       ..
16*       .. Array Arguments ..
17*       DOUBLE PRECISION DPARAM(5)
18*       ..
19*
20*
21*> \par Purpose:
22*  =============
23*>
24*> \verbatim
25*>
26*>    CONSTRUCT THE MODIFIED GIVENS TRANSFORMATION MATRIX H WHICH ZEROS
27*>    THE SECOND COMPONENT OF THE 2-VECTOR  (DSQRT(DD1)*DX1,DSQRT(DD2)*>    DY2)**T.
28*>    WITH DPARAM(1)=DFLAG, H HAS ONE OF THE FOLLOWING FORMS..
29*>
30*>    DFLAG=-1.D0     DFLAG=0.D0        DFLAG=1.D0     DFLAG=-2.D0
31*>
32*>      (DH11  DH12)    (1.D0  DH12)    (DH11  1.D0)    (1.D0  0.D0)
33*>    H=(          )    (          )    (          )    (          )
34*>      (DH21  DH22),   (DH21  1.D0),   (-1.D0 DH22),   (0.D0  1.D0).
35*>    LOCATIONS 2-4 OF DPARAM CONTAIN DH11, DH21, DH12, AND DH22
36*>    RESPECTIVELY. (VALUES OF 1.D0, -1.D0, OR 0.D0 IMPLIED BY THE
37*>    VALUE OF DPARAM(1) ARE NOT STORED IN DPARAM.)
38*>
39*>    THE VALUES OF GAMSQ AND RGAMSQ SET IN THE DATA STATEMENT MAY BE
40*>    INEXACT.  THIS IS OK AS THEY ARE ONLY USED FOR TESTING THE SIZE
41*>    OF DD1 AND DD2.  ALL ACTUAL SCALING OF DATA IS DONE USING GAM.
42*>
43*> \endverbatim
44*
45*  Arguments:
46*  ==========
47*
48*> \param[in,out] DD1
49*> \verbatim
50*>          DD1 is DOUBLE PRECISION
51*> \endverbatim
52*>
53*> \param[in,out] DD2
54*> \verbatim
55*>          DD2 is DOUBLE PRECISION
56*> \endverbatim
57*>
58*> \param[in,out] DX1
59*> \verbatim
60*>          DX1 is DOUBLE PRECISION
61*> \endverbatim
62*>
63*> \param[in] DY1
64*> \verbatim
65*>          DY1 is DOUBLE PRECISION
66*> \endverbatim
67*>
68*> \param[out] DPARAM
69*> \verbatim
70*>          DPARAM is DOUBLE PRECISION array, dimension (5)
71*>     DPARAM(1)=DFLAG
72*>     DPARAM(2)=DH11
73*>     DPARAM(3)=DH21
74*>     DPARAM(4)=DH12
75*>     DPARAM(5)=DH22
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 double_blas_level1
87*
88*  =====================================================================
89      SUBROUTINE DROTMG(DD1,DD2,DX1,DY1,DPARAM)
90*
91*  -- Reference BLAS level1 routine --
92*  -- Reference BLAS 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 DD1,DD2,DX1,DY1
97*     ..
98*     .. Array Arguments ..
99      DOUBLE PRECISION DPARAM(5)
100*     ..
101*
102*  =====================================================================
103*
104*     .. Local Scalars ..
105      DOUBLE PRECISION DFLAG,DH11,DH12,DH21,DH22,DP1,DP2,DQ1,DQ2,DTEMP,
106     $                 DU,GAM,GAMSQ,ONE,RGAMSQ,TWO,ZERO
107*     ..
108*     .. Intrinsic Functions ..
109      INTRINSIC DABS
110*     ..
111*     .. Data statements ..
112*
113      DATA ZERO,ONE,TWO/0.D0,1.D0,2.D0/
114      DATA GAM,GAMSQ,RGAMSQ/4096.D0,16777216.D0,5.9604645D-8/
115*     ..
116
117      IF (DD1.LT.ZERO) THEN
118*        GO ZERO-H-D-AND-DX1..
119         DFLAG = -ONE
120         DH11 = ZERO
121         DH12 = ZERO
122         DH21 = ZERO
123         DH22 = ZERO
124*
125         DD1 = ZERO
126         DD2 = ZERO
127         DX1 = ZERO
128      ELSE
129*        CASE-DD1-NONNEGATIVE
130         DP2 = DD2*DY1
131         IF (DP2.EQ.ZERO) THEN
132            DFLAG = -TWO
133            DPARAM(1) = DFLAG
134            RETURN
135         END IF
136*        REGULAR-CASE..
137         DP1 = DD1*DX1
138         DQ2 = DP2*DY1
139         DQ1 = DP1*DX1
140*
141         IF (DABS(DQ1).GT.DABS(DQ2)) THEN
142            DH21 = -DY1/DX1
143            DH12 = DP2/DP1
144*
145            DU = ONE - DH12*DH21
146*
147           IF (DU.GT.ZERO) THEN
148             DFLAG = ZERO
149             DD1 = DD1/DU
150             DD2 = DD2/DU
151             DX1 = DX1*DU
152           ELSE
153*            This code path if here for safety. We do not expect this
154*            condition to ever hold except in edge cases with rounding
155*            errors. See DOI: 10.1145/355841.355847
156             DFLAG = -ONE
157             DH11 = ZERO
158             DH12 = ZERO
159             DH21 = ZERO
160             DH22 = ZERO
161*
162             DD1 = ZERO
163             DD2 = ZERO
164             DX1 = ZERO
165           END IF
166         ELSE
167
168            IF (DQ2.LT.ZERO) THEN
169*              GO ZERO-H-D-AND-DX1..
170               DFLAG = -ONE
171               DH11 = ZERO
172               DH12 = ZERO
173               DH21 = ZERO
174               DH22 = ZERO
175*
176               DD1 = ZERO
177               DD2 = ZERO
178               DX1 = ZERO
179            ELSE
180               DFLAG = ONE
181               DH11 = DP1/DP2
182               DH22 = DX1/DY1
183               DU = ONE + DH11*DH22
184               DTEMP = DD2/DU
185               DD2 = DD1/DU
186               DD1 = DTEMP
187               DX1 = DY1*DU
188            END IF
189         END IF
190
191*     PROCEDURE..SCALE-CHECK
192         IF (DD1.NE.ZERO) THEN
193            DO WHILE ((DD1.LE.RGAMSQ) .OR. (DD1.GE.GAMSQ))
194               IF (DFLAG.EQ.ZERO) THEN
195                  DH11 = ONE
196                  DH22 = ONE
197                  DFLAG = -ONE
198               ELSE
199                  DH21 = -ONE
200                  DH12 = ONE
201                  DFLAG = -ONE
202               END IF
203               IF (DD1.LE.RGAMSQ) THEN
204                  DD1 = DD1*GAM**2
205                  DX1 = DX1/GAM
206                  DH11 = DH11/GAM
207                  DH12 = DH12/GAM
208               ELSE
209                  DD1 = DD1/GAM**2
210                  DX1 = DX1*GAM
211                  DH11 = DH11*GAM
212                  DH12 = DH12*GAM
213               END IF
214            ENDDO
215         END IF
216
217         IF (DD2.NE.ZERO) THEN
218            DO WHILE ( (DABS(DD2).LE.RGAMSQ) .OR. (DABS(DD2).GE.GAMSQ) )
219               IF (DFLAG.EQ.ZERO) THEN
220                  DH11 = ONE
221                  DH22 = ONE
222                  DFLAG = -ONE
223               ELSE
224                  DH21 = -ONE
225                  DH12 = ONE
226                  DFLAG = -ONE
227               END IF
228               IF (DABS(DD2).LE.RGAMSQ) THEN
229                  DD2 = DD2*GAM**2
230                  DH21 = DH21/GAM
231                  DH22 = DH22/GAM
232               ELSE
233                  DD2 = DD2/GAM**2
234                  DH21 = DH21*GAM
235                  DH22 = DH22*GAM
236               END IF
237            END DO
238         END IF
239
240      END IF
241
242      IF (DFLAG.LT.ZERO) THEN
243         DPARAM(2) = DH11
244         DPARAM(3) = DH21
245         DPARAM(4) = DH12
246         DPARAM(5) = DH22
247      ELSE IF (DFLAG.EQ.ZERO) THEN
248         DPARAM(3) = DH21
249         DPARAM(4) = DH12
250      ELSE
251         DPARAM(2) = DH11
252         DPARAM(5) = DH22
253      END IF
254
255      DPARAM(1) = DFLAG
256      RETURN
257*
258*     End of DROTMG
259*
260      END
261