[CONTACT]

[ABOUT]

[POLICY]

[ADVERTISE]

SUBROUTINE N,NN,P,QTR,

Found at: ftp.icm.edu.pl:70/packages/netlib/port/dq7rad.f

      SUBROUTINE DQ7RAD(N, NN, P, QTR, QTRSET, RMAT, W, Y)
C
C  ***  ADD ROWS W TO QR FACTORIZATION WITH R MATRIX RMAT AND
C  ***  Q**T * RESIDUAL = QTR.  Y = NEW COMPONENTS OF RESIDUAL
C  ***  CORRESPONDING TO W.  QTR, Y REFERENCED ONLY IF QTRSET = .TRUE.
C
      LOGICAL QTRSET
      INTEGER N, NN, P
      DOUBLE PRECISION QTR(P), RMAT(1), W(NN,P), Y(N)
C     DIMENSION RMAT(P*(P+1)/2)
C/+
      DOUBLE PRECISION DSQRT
C/
      DOUBLE PRECISION DD7TPR, DR7MDC, DV2NRM
      EXTERNAL DD7TPR, DR7MDC,DV2AXY, DV7SCL, DV2NRM
C
C  ***  LOCAL VARIABLES  ***
C
      INTEGER I, II, IJ, IP1, J, K, NK
      DOUBLE PRECISION ARI, QRI, RI, S, T, WI
      DOUBLE PRECISION BIG, BIGRT, ONE, TINY, TINYRT, ZERO
C/7
      SAVE BIGRT, TINY, TINYRT
C/
      DATA BIG/-1.D+0/, BIGRT/-1.D+0/, ONE/1.D+0/, TINY/0.D+0/,
     1     TINYRT/0.D+0/, ZERO/0.D+0/
C
C------------------------------ BODY -----------------------------------
C
      IF (TINY .GT. ZERO) GO TO 10
         TINY = DR7MDC(1)
         BIG = DR7MDC(6)
         IF (TINY*BIG .LT. ONE) TINY = ONE / BIG
 10   K = 1
      NK = N
      II = 0
      DO 180 I = 1, P
         II = II + I
         IP1 = I + 1
         IJ = II + I
         IF (NK .LE. 1) T = DABS(W(K,I))
         IF (NK .GT. 1) T = DV2NRM(NK, W(K,I))
         IF (T .LT. TINY) GOTO  180
         RI = RMAT(II)
         IF (RI .NE. ZERO) GO TO 100
            IF (NK .GT. 1) GO TO 30
               IJ = II
               DO 20 J = I, P
                  RMAT(IJ) = W(K,J)
                  IJ = IJ + J
 20               CONTINUE
               IF (QTRSET) QTR(I) = Y(K)
               W(K,I) = ZERO
               GO TO 999
 30         WI = W(K,I)
            IF (BIGRT .GT. ZERO) GO TO 40
               BIGRT = DR7MDC(5)
               TINYRT = DR7MDC(2)
 40         IF (T .LE. TINYRT) GO TO 50
            IF (T .GE. BIGRT) GO TO 50
               IF (WI .LT. ZERO) T = -T
               WI = WI + T
               S = DSQRT(T * WI)
               GO TO 70
 50         S = DSQRT(T)
            IF (WI .LT. ZERO) GO TO 60
               WI = WI + T
               S = S * DSQRT(WI)
               GO TO 70
 60         T = -T
            WI = WI + T
            S = S * DSQRT(-WI)
 70         W(K,I) = WI
            CALL DV7SCL(NK, W(K,I), ONE/S, W(K,I))
            RMAT(II) = -T
            IF (.NOT. QTRSET) GO TO 80
            CALL DV2AXY(NK, Y(K), -DD7TPR(NK,Y(K),W(K,I)), W(K,I), Y(K))
            QTR(I) = Y(K)
 80         IF (IP1 .GT. P) GO TO 999
            DO 90 J = IP1, P
               CALL DV2AXY(NK, W(K,J), -DD7TPR(NK,W(K,J),W(K,I)),
     1                    W(K,I), W(K,J))
               RMAT(IJ) = W(K,J)
               IJ = IJ + J
 90            CONTINUE
            IF (NK .LE. 1) GO TO 999
            K = K + 1
            NK = NK - 1
            GO TO 180
C
 100     ARI = DABS(RI)
         IF (ARI .GT. T) GO TO 110
            T = T * DSQRT(ONE + (ARI/T)**2)
            GO TO 120
 110     T = ARI * DSQRT(ONE + (T/ARI)**2)
 120     IF (RI .LT. ZERO) T = -T
         RI = RI + T
         RMAT(II) = -T
         S = -RI / T
         IF (NK .LE. 1) GO TO 150
         CALL DV7SCL(NK, W(K,I), ONE/RI, W(K,I))
         IF (.NOT. QTRSET) GO TO 130
            QRI = QTR(I)
            T = S * ( QRI  +  DD7TPR(NK, Y(K), W(K,I)) )
            QTR(I) = QRI + T
 130     IF (IP1 .GT. P) GO TO 999
         IF (QTRSET) CALL DV2AXY(NK, Y(K), T, W(K,I), Y(K))
         DO 140 J = IP1, P
            RI = RMAT(IJ)
            T = S * ( RI  +  DD7TPR(NK, W(K,J), W(K,I)) )
            CALL DV2AXY(NK, W(K,J), T, W(K,I), W(K,J))
            RMAT(IJ) = RI + T
            IJ = IJ + J
 140        CONTINUE
         GO TO 180
C
 150     WI = W(K,I) / RI
         W(K,I) = WI
         IF (.NOT. QTRSET) GO TO 160
            QRI = QTR(I)
            T = S * ( QRI + Y(K)*WI )
            QTR(I) = QRI + T
 160     IF (IP1 .GT. P) GO TO 999
         IF (QTRSET) Y(K) = T*WI + Y(K)
         DO 170 J = IP1, P
            RI = RMAT(IJ)
            T = S * (RI + W(K,J)*WI)
            W(K,J) = W(K,J) + T*WI
            RMAT(IJ) = RI + T
            IJ = IJ + J
 170        CONTINUE
 180     CONTINUE
C
 999  RETURN
C  ***  LAST LINE OF DQ7RAD FOLLOWS  ***
      END

		
.


AD:

NEW PAGES:

[ODDNUGGET]

[GOPHER]