1
2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 |
SUBROUTINE DROTM(N,DX,INCX,DY,INCY,DPARAM)
* .. Scalar Arguments .. INTEGER INCX,INCY,N * .. * .. Array Arguments .. DOUBLE PRECISION DPARAM(5),DX(*),DY(*) * .. * * Purpose * ======= * * APPLY THE MODIFIED GIVENS TRANSFORMATION, H, TO THE 2 BY N MATRIX * * (DX**T) , WHERE **T INDICATES TRANSPOSE. THE ELEMENTS OF DX ARE IN * (DY**T) * * DX(LX+I*INCX), I = 0 TO N-1, WHERE LX = 1 IF INCX .GE. 0, ELSE * LX = (-INCX)*N, AND SIMILARLY FOR SY USING LY AND INCY. * WITH DPARAM(1)=DFLAG, H HAS ONE OF THE FOLLOWING FORMS.. * * DFLAG=-1.D0 DFLAG=0.D0 DFLAG=1.D0 DFLAG=-2.D0 * * (DH11 DH12) (1.D0 DH12) (DH11 1.D0) (1.D0 0.D0) * H=( ) ( ) ( ) ( ) * (DH21 DH22), (DH21 1.D0), (-1.D0 DH22), (0.D0 1.D0). * SEE DROTMG FOR A DESCRIPTION OF DATA STORAGE IN DPARAM. * * Arguments * ========= * * N (input) INTEGER * number of elements in input vector(s) * * DX (input/output) DOUBLE PRECISION array, dimension N * double precision vector with N elements * * INCX (input) INTEGER * storage spacing between elements of DX * * DY (input/output) DOUBLE PRECISION array, dimension N * double precision vector with N elements * * INCY (input) INTEGER * storage spacing between elements of DY * * DPARAM (input/output) DOUBLE PRECISION array, dimension 5 * DPARAM(1)=DFLAG * DPARAM(2)=DH11 * DPARAM(3)=DH21 * DPARAM(4)=DH12 * DPARAM(5)=DH22 * * ===================================================================== * * .. Local Scalars .. DOUBLE PRECISION DFLAG,DH11,DH12,DH21,DH22,TWO,W,Z,ZERO INTEGER I,KX,KY,NSTEPS * .. * .. Data statements .. DATA ZERO,TWO/0.D0,2.D0/ * .. * DFLAG = DPARAM(1) IF (N.LE.0 .OR. (DFLAG+TWO.EQ.ZERO)) RETURN IF (INCX.EQ.INCY.AND.INCX.GT.0) THEN * NSTEPS = N*INCX IF (DFLAG.LT.ZERO) THEN DH11 = DPARAM(2) DH12 = DPARAM(4) DH21 = DPARAM(3) DH22 = DPARAM(5) DO I = 1,NSTEPS,INCX W = DX(I) Z = DY(I) DX(I) = W*DH11 + Z*DH12 DY(I) = W*DH21 + Z*DH22 END DO ELSE IF (DFLAG.EQ.ZERO) THEN DH12 = DPARAM(4) DH21 = DPARAM(3) DO I = 1,NSTEPS,INCX W = DX(I) Z = DY(I) DX(I) = W + Z*DH12 DY(I) = W*DH21 + Z END DO ELSE DH11 = DPARAM(2) DH22 = DPARAM(5) DO I = 1,NSTEPS,INCX W = DX(I) Z = DY(I) DX(I) = W*DH11 + Z DY(I) = -W + DH22*Z END DO END IF ELSE KX = 1 KY = 1 IF (INCX.LT.0) KX = 1 + (1-N)*INCX IF (INCY.LT.0) KY = 1 + (1-N)*INCY * IF (DFLAG.LT.ZERO) THEN DH11 = DPARAM(2) DH12 = DPARAM(4) DH21 = DPARAM(3) DH22 = DPARAM(5) DO I = 1,N W = DX(KX) Z = DY(KY) DX(KX) = W*DH11 + Z*DH12 DY(KY) = W*DH21 + Z*DH22 KX = KX + INCX KY = KY + INCY END DO ELSE IF (DFLAG.EQ.ZERO) THEN DH12 = DPARAM(4) DH21 = DPARAM(3) DO I = 1,N W = DX(KX) Z = DY(KY) DX(KX) = W + Z*DH12 DY(KY) = W*DH21 + Z KX = KX + INCX KY = KY + INCY END DO ELSE DH11 = DPARAM(2) DH22 = DPARAM(5) DO I = 1,N W = DX(KX) Z = DY(KY) DX(KX) = W*DH11 + Z DY(KY) = -W + DH22*Z KX = KX + INCX KY = KY + INCY END DO END IF END IF RETURN END |