text: *DECK DU12US
        SUBROUTINE DU12US (A, MDA, M, N, B, MDB, NB, MODE, KRANK, RNORM,
       +   H, W, IR, IC)
  C***BEGIN PROLOGUE  DU12US
  C***SUBSIDIARY
  C***PURPOSE  Subsidiary to DULSIA
  C***LIBRARY   SLATEC
  C***TYPE      DOUBLE PRECISION (U12US-S, DU12US-D)
  C***AUTHOR  (UNKNOWN)
  C***DESCRIPTION
  C
  C        Given the Householder LQ factorization of A, this
  C        subroutine solves the system AX=B. If the system
  C        is of reduced rank, this routine returns a solution
  C        according to the selected mode.
  C
  C       Note - If MODE.NE.2, W is never accessed.
  C
  C***SEE ALSO  DULSIA
  C***ROUTINES CALLED  DAXPY, DDOT, DNRM2, DSWAP
  C***REVISION HISTORY  (YYMMDD)
  C   810801  DATE WRITTEN
  C   890531  Changed all specific intrinsics to generic.  (WRB)
  C   890831  Modified array declarations.  (WRB)
  C   891214  Prologue converted to Version 4.0 format.  (BAB)
  C   900328  Added TYPE section.  (WRB)
  C***END PROLOGUE  DU12US
        IMPLICIT DOUBLE PRECISION (A-H,O-Z)
        DOUBLE PRECISION DDOT,DNRM2
        DIMENSION A(MDA,*),B(MDB,*),RNORM(*),H(*),W(*)
        INTEGER IC(*),IR(*)
  C***FIRST EXECUTABLE STATEMENT  DU12US
        K=KRANK
        KP1=K+1
  C
  C        RANK=0
  C
        IF(K.GT.0) GO TO 410
        DO 404 JB=1,NB
        RNORM(JB)=DNRM2(M,B(1,JB),1)
    404 CONTINUE
        DO 406 JB=1,NB
        DO 406 I=1,N
        B(I,JB)=0.0D0
    406 CONTINUE
        RETURN
  C
  C     REORDER B TO REFLECT ROW INTERCHANGES
  C
    410 CONTINUE
        I=0
    412 I=I+1
        IF(I.EQ.M) GO TO 418
        J=IR(I)
        IF(J.EQ.I) GO TO 412
        IF(J.LT.0) GO TO 412
        IR(I)=-IR(I)
        DO 413 JB=1,NB
        RNORM(JB)=B(I,JB)
    413 CONTINUE
        IJ=I
    414 DO 415 JB=1,NB
        B(IJ,JB)=B(J,JB)
    415 CONTINUE
        IJ=J
        J=IR(IJ)
        IR(IJ)=-IR(IJ)
        IF(J.NE.I) GO TO 414
        DO 416 JB=1,NB
        B(IJ,JB)=RNORM(JB)
    416 CONTINUE
        GO TO 412
    418 CONTINUE
        DO 420 I=1,M
        IR(I)=ABS(IR(I))
    420 CONTINUE
  C
  C     IF A IS OF REDUCED RANK AND MODE=2,
  C     APPLY HOUSEHOLDER TRANSFORMATIONS TO B
  C
        IF(MODE.LT.2 .OR. K.EQ.M) GO TO 440
        MMK=M-K
        DO 430 JB=1,NB
        DO 425 J=1,K
        I=KP1-J
        TT=-DDOT(MMK,A(KP1,I),1,B(KP1,JB),1)/W(I)
        TT=TT-B(I,JB)
        CALL DAXPY(MMK,TT,A(KP1,I),1,B(KP1,JB),1)
        B(I,JB)=B(I,JB)+TT*W(I)
    425 CONTINUE
    430 CONTINUE
  C
  C     FIND NORMS OF RESIDUAL VECTOR(S)..(BEFORE OVERWRITE B)
  C
    440 DO 442 JB=1,NB
        RNORM(JB)=DNRM2((M-K),B(KP1,JB),1)
    442 CONTINUE
  C
  C     BACK SOLVE LOWER TRIANGULAR L
  C
        DO 450 JB=1,NB
        DO 448 I=1,K
        B(I,JB)=B(I,JB)/A(I,I)
        IF(I.EQ.K) GO TO 450
        IP1=I+1
        CALL DAXPY(K-I,-B(I,JB),A(IP1,I),1,B(IP1,JB),1)
    448 CONTINUE
    450 CONTINUE
  C
  C
  C      TRUNCATED SOLUTION
  C
        IF(K.EQ.N) GO TO 462
        DO 460 JB=1,NB
        DO 460 I=KP1,N
        B(I,JB)=0.0D0
    460 CONTINUE
  C
  C     APPLY HOUSEHOLDER TRANSFORMATIONS TO B
  C
    462 DO 470 I=1,K
        J=KP1-I
        TT=A(J,J)
        A(J,J)=H(J)
        DO 465 JB=1,NB
        BB=-DDOT(N-J+1,A(J,J),MDA,B(J,JB),1)/H(J)
        CALL DAXPY(N-J+1,BB,A(J,J),MDA,B(J,JB),1)
    465 CONTINUE
        A(J,J)=TT
    470 CONTINUE
  C
  C
  C     REORDER B TO REFLECT COLUMN INTERCHANGES
  C
        I=0
    482 I=I+1
        IF(I.EQ.N) GO TO 488
        J=IC(I)
        IF(J.EQ.I) GO TO 482
        IF(J.LT.0) GO TO 482
        IC(I)=-IC(I)
    484 CALL DSWAP(NB,B(J,1),MDB,B(I,1),MDB)
        IJ=IC(J)
        IC(J)=-IC(J)
        J=IJ
        IF(J.EQ.I) GO TO 482
        GO TO 484
    488 CONTINUE
        DO 490 I=1,N
        IC(I)=ABS(IC(I))
    490 CONTINUE
  C
  C        SOLUTION VECTORS ARE IN FIRST N ROWS OF B(,)
  C
        RETURN
        END
  