SUBROUTINE DMDLVS(N0,A,B,WORK,INFO)
  
      IMPLICIT NONE
  
      DOUBLE PRECISION A(*), B(*), WORK(*)
      INTEGER N0,INFO,IINFO,MAXITER
      INTEGER N, M, I, J, K
      INTEGER INDRV1,INDRV2,INDRV3,INDRV4,INDRV5,INDRV6,INDRV7,INDRV8
      DOUBLE PRECISION TMP1, TMP2, TMP3, TAU, TAU1, TAU2, TAU3, S, T
      DOUBLE PRECISION TMP4, TMP5, TMP6, TMP7, TMP8, TMP9
      DOUBLE PRECISION EPS,SAFMIN,SCALE,SIGMX,TOL,TOL2
      DOUBLE PRECISION SIGMA,SIGMA1,DESIG,DESIG0
      DOUBLE PRECISION SCALE1
      DOUBLE PRECISION QMIN,EMAX
      DOUBLE PRECISION DN2,DN1,DN0

      DOUBLE PRECISION ONE, ZERO, HALF, TWO, HUNDRD, CONST
      PARAMETER ( ONE=1.0D0, ZERO=0.0D0, HALF=0.5D0 )
      PARAMETER ( TWO=2.0D0, HUNDRD=100.0D0, CONST=0.75D0 )

      EXTERNAL DLAMCH,DMDLVS2
      DOUBLE PRECISION DLAMCH
      EXTERNAL DISNAN
      LOGICAL DISNAN

      IF (N0 .EQ. 1) THEN
         A(1)=ABS(A(1))
         INFO = 0
         RETURN
      ENDIF

      IF (N0 .EQ. 2) THEN
         CALL DLAS2(A(1),B(1),A(2),A(2),A(1))
         INFO = 0
         RETURN
      ENDIF

      INDRV1 = 0
      INDRV2 = INDRV1+N0
      INDRV3 = INDRV2+N0
      INDRV4 = INDRV3+N0
      INDRV5 = INDRV4+N0
      INDRV6 = INDRV5+N0
      INDRV7 = INDRV6+N0
      INDRV8 = INDRV7+N0
  
      EPS = DLAMCH( 'Precision' )
      TOL = HUNDRD*EPS
      TOL2 = TOL**2
      SAFMIN = DLAMCH( 'Safe minimum' )
      SCALE = SQRT( HALF / SAFMIN )
      SCALE1 = ONE/SCALE
  
      WORK(INDRV1+1:INDRV1+N0)=ABS(A(1:N0))
      WORK(INDRV2+1:INDRV2+N0-1)=ABS(B(1:N0-1))

      SIGMX = MAX(MAXVAL(WORK(INDRV1+1:INDRV1+N0)),
     $     MAXVAL(WORK(INDRV2+1:INDRV2+N0-1)))

      IF (SIGMX .GT. ZERO) THEN
         CALL DLASCL( 'G', 0, 0, SIGMX, SCALE, N0, 1, 
     $        WORK(INDRV1+1), N0, IINFO )
         CALL DLASCL( 'G', 0, 0, SIGMX, SCALE, N0-1, 1, 
     $        WORK(INDRV2+1), N0-1, IINFO )
      ENDIF

      WORK(INDRV1+1:INDRV1+N0) = WORK(INDRV1+1:INDRV1+N0)**2
      WORK(INDRV2+1:INDRV2+N0-1) = WORK(INDRV2+1:INDRV2+N0-1)**2
  
      M = 1
      N = N0
*     
      IF (WORK(INDRV1+M) .LT. WORK(INDRV1+N)) THEN
         K=(N-M+1)/2+M-1
         DO J=M,K
            TMP1=WORK(INDRV1+J)
            WORK(INDRV1+J)=WORK(INDRV1+N+M-J)
            WORK(INDRV1+N+M-J)=TMP1
            TMP1=WORK(INDRV2+J)
            WORK(INDRV2+J)=WORK(INDRV2+N-1+M-J)
            WORK(INDRV2+N-1+M-J)=TMP1
         ENDDO
      END IF

      TMP1 = WORK(INDRV1+M)
      DO J = M, N-3
         IF (WORK(INDRV2+J) .LE. TOL2*TMP1) THEN
            WORK(INDRV2+J) = -ZERO
            WORK(INDRV4+J) = -ZERO
            M = J+1
            TMP1 = WORK(INDRV1+J+1)
         ELSE
            TMP1 = WORK(INDRV1+J+1)*(TMP1/(TMP1+WORK(INDRV2+J)))
         ENDIF
      ENDDO
*     
      IF (WORK(INDRV2+N-2) .LE. TOL2*TMP1) THEN
         WORK(INDRV2+N-2) = ZERO
         TMP1 = WORK(INDRV1+N-1)
      ELSE
         TMP1 = WORK(INDRV1+N-1)*(TMP1/(TMP1+WORK(INDRV2+N-2)))
      ENDIF
*     
      IF (WORK(INDRV2+N-1) .LE. TOL2*TMP1) THEN
         WORK(INDRV2+N-1) = ZERO
      ENDIF
*
      WORK(INDRV2+N) = ZERO
      WORK(INDRV4+N) = ZERO

 3000 SIGMA = -WORK(INDRV2+N)
      SIGMA1 = TOL2*SIGMA
      DESIG = -WORK(INDRV4+N)

      MAXITER=30*(N-M+1)
      DO I=1, MAXITER, 2
     
 15      IF (N-M+1 .EQ. 1) THEN
            
            CALL DMDLVS2(SIGMA,DESIG,WORK(INDRV1+N),T,DESIG0)
            WORK(INDRV1+N)=T
            
            GO TO 700
         ENDIF
     
         IF (N-M+1 .EQ. 2) THEN
            
            IF( WORK(INDRV1+N).GT.WORK(INDRV1+N-1) ) THEN
               S = WORK(INDRV1+N)
               WORK(INDRV1+N) = WORK(INDRV1+N-1)
               WORK(INDRV1+N-1) = S
            END IF
            T = HALF*( ( WORK(INDRV1+N-1)-WORK(INDRV1+N) )
     $           +WORK(INDRV2+N-1) )
            IF(WORK(INDRV2+N-1) .GT. WORK(INDRV1+N)*TOL2
     $           .AND. T .NE. ZERO) THEN
               S = WORK(INDRV1+N)*( WORK(INDRV2+N-1) / T )
               IF( S.LE.T ) THEN
                  S = WORK(INDRV1+N)*( WORK(INDRV2+N-1) / 
     $                 ( T*( ONE+SQRT( ONE+S / T ) ) ) )
               ELSE
                  S = WORK(INDRV1+N)*( WORK(INDRV2+N-1) / 
     $                 ( T+SQRT( T )*SQRT( T+S ) ) )
               END IF
               T = WORK(INDRV1+N-1) + ( S+WORK(INDRV2+N-1) )
               WORK(INDRV1+N) = WORK(INDRV1+N)
     $              *( WORK(INDRV1+N-1) / T )
               WORK(INDRV1+N-1) = T
            END IF
            
            CALL DMDLVS2(SIGMA,DESIG,WORK(INDRV1+N-1),T,DESIG0)
            WORK(INDRV1+N-1)=T
            CALL DMDLVS2(SIGMA,DESIG,WORK(INDRV1+N),T,DESIG0)
            WORK(INDRV1+N)=T
            
            GO TO 700
         ENDIF

         IF (WORK(INDRV2+N-2) .LE. SIGMA1) THEN
            
            IF( WORK(INDRV1+N).GT.WORK(INDRV1+N-1) ) THEN
               S = WORK(INDRV1+N)
               WORK(INDRV1+N) = WORK(INDRV1+N-1)
               WORK(INDRV1+N-1) = S
            END IF
            T = HALF*( ( WORK(INDRV1+N-1)-WORK(INDRV1+N) )
     $           +WORK(INDRV2+N-1) )
            IF(WORK(INDRV2+N-1) .GT. WORK(INDRV1+N)*TOL2
     $           .AND. T .NE. ZERO) THEN
               S = WORK(INDRV1+N)*( WORK(INDRV2+N-1) / T )
               IF( S.LE.T ) THEN
                  S = WORK(INDRV1+N)*( WORK(INDRV2+N-1) / 
     $                 ( T*( ONE+SQRT( ONE+S / T ) ) ) )
               ELSE
                  S = WORK(INDRV1+N)*( WORK(INDRV2+N-1) / 
     $                 ( T+SQRT( T )*SQRT( T+S ) ) )
               END IF
               T = WORK(INDRV1+N-1) + ( S+WORK(INDRV2+N-1) )
               WORK(INDRV1+N) = WORK(INDRV1+N)
     $              *( WORK(INDRV1+N-1) / T )
               WORK(INDRV1+N-1) = T
            END IF
            
            CALL DMDLVS2(SIGMA,DESIG,WORK(INDRV1+N-1),T,DESIG0)
            WORK(INDRV1+N-1)=T
            CALL DMDLVS2(SIGMA,DESIG,WORK(INDRV1+N),T,DESIG0)
            WORK(INDRV1+N)=T
            
            N = N - 2
            GO TO 15
         ENDIF

         IF (WORK(INDRV2+N-1) .LE. SIGMA1) THEN
            
            CALL DMDLVS2(SIGMA,DESIG,WORK(INDRV1+N),T,DESIG0)
            WORK(INDRV1+N)=T
            
            N = N - 1
            GO TO 15
         ENDIF
     
         IF (WORK(INDRV1+M) .LT. WORK(INDRV1+N)) THEN
            K=(N-M+1)/2+M-1
            DO J=M,K
               TMP1=WORK(INDRV1+J)
               WORK(INDRV1+J)=WORK(INDRV1+N+M-J)
               WORK(INDRV1+N+M-J)=TMP1
               TMP1=WORK(INDRV2+J)
               WORK(INDRV2+J)=WORK(INDRV2+N-1+M-J)
               WORK(INDRV2+N-1+M-J)=TMP1
            ENDDO
         END IF

         DN2 = WORK(INDRV1+M)
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         DO J=M, N-1
            TMP2 = WORK(INDRV2+J)/DN2
            TMP3 = ONE+TMP2
            WORK(INDRV3+J) = DN2*TMP3
            DN2 = WORK(INDRV1+J+1)/TMP3
            WORK(INDRV4+J) = TMP2*DN2
            IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         ENDDO
         WORK(INDRV3+N) = DN2
         
         IF (DISNAN(WORK(INDRV4+N-1))) GO TO 360
         
         GO TO 330

 360     DN2 = WORK(INDRV1+M)
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         DO J=M, N-3
            WORK(INDRV3+J) = DN2+WORK(INDRV2+J)
            TMP3 = WORK(INDRV1+J+1)/WORK(INDRV3+J)
            WORK(INDRV4+J) = WORK(INDRV2+J)*TMP3
            DN2 = DN2*TMP3
            IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         ENDDO

         WORK(INDRV3+N-2) = DN2+WORK(INDRV2+N-2)
         IF (SAFMIN*WORK(INDRV3+N-2) .LT. WORK(INDRV1+N-1) .AND.
     $        SAFMIN*WORK(INDRV1+N-1) .LT. WORK(INDRV3+N-2)) THEN
            TMP2=WORK(INDRV1+N-1)/WORK(INDRV3+N-2)
            WORK(INDRV4+N-2)=WORK(INDRV2+N-2)*TMP2
            DN2=DN2*TMP2
         ELSE
            WORK(INDRV4+N-2) = WORK(INDRV1+N-1)
     $           *(WORK(INDRV2+N-2)/WORK(INDRV3+N-2))
            DN2 = WORK(INDRV1+N-1)*(DN2/WORK(INDRV3+N-2))
         ENDIF
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO

         WORK(INDRV3+N-1) = DN2+WORK(INDRV2+N-1)
         IF (SAFMIN*WORK(INDRV3+N-1) .LT. WORK(INDRV1+N) .AND.
     $        SAFMIN*WORK(INDRV1+N) .LT. WORK(INDRV3+N-1)) THEN
            TMP2=WORK(INDRV1+N)/WORK(INDRV3+N-1)
            WORK(INDRV4+N-1)=WORK(INDRV2+N-1)*TMP2
            DN2=DN2*TMP2
         ELSE
            WORK(INDRV4+N-1) = WORK(INDRV1+N)
     $           *(WORK(INDRV2+N-1)/WORK(INDRV3+N-1))
            DN2 = WORK(INDRV1+N)*(DN2/WORK(INDRV3+N-1))
         ENDIF
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO

         WORK(INDRV3+N)=DN2

         IF (DISNAN(DN2)) GO TO 370

         GO TO 330
      
 370     DN2 = WORK(INDRV1+M)
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         DO J=M, N-3
            WORK(INDRV3+J) = DN2+WORK(INDRV2+J)
            IF (SAFMIN*WORK(INDRV3+J) .LT. WORK(INDRV1+J+1) .AND.  
     $           SAFMIN*WORK(INDRV1+J+1) .LT. WORK(INDRV3+J)) THEN
               TMP2=WORK(INDRV1+J+1)/WORK(INDRV3+J)
               WORK(INDRV4+J)=WORK(INDRV2+J)*TMP2
               DN2=DN2*TMP2
            ELSE
               WORK(INDRV4+J) = WORK(INDRV1+J+1)
     $              *(WORK(INDRV2+J)/WORK(INDRV3+J))
               DN2 = WORK(INDRV1+J+1)*(DN2/WORK(INDRV3+J))
            ENDIF
            IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         ENDDO

         WORK(INDRV3+N-2) = DN2+WORK(INDRV2+N-2)
         IF (SAFMIN*WORK(INDRV3+N-2) .LT. WORK(INDRV1+N-1) .AND.  
     $        SAFMIN*WORK(INDRV1+N-1) .LT. WORK(INDRV3+N-2)) THEN
            TMP2=WORK(INDRV1+N-1)/WORK(INDRV3+N-2)
            WORK(INDRV4+N-2)=WORK(INDRV2+N-2)*TMP2
            DN2=DN2*TMP2
         ELSE
            WORK(INDRV4+N-2) = WORK(INDRV1+N-1)
     $           *(WORK(INDRV2+N-2)/WORK(INDRV3+N-2))
            DN2 = WORK(INDRV1+N-1)*(DN2/WORK(INDRV3+N-2))
         ENDIF
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO

         WORK(INDRV3+N-1) = DN2+WORK(INDRV2+N-1)
         IF (SAFMIN*WORK(INDRV3+N-1) .LT. WORK(INDRV1+N) .AND.  
     $        SAFMIN*WORK(INDRV1+N) .LT. WORK(INDRV3+N-1)) THEN
            TMP2=WORK(INDRV1+N)/WORK(INDRV3+N-1)
            WORK(INDRV4+N-1)=WORK(INDRV2+N-1)*TMP2
            DN2=DN2*TMP2
         ELSE
            WORK(INDRV4+N-1) = WORK(INDRV1+N)
     $           *(WORK(INDRV2+N-1)/WORK(INDRV3+N-1))
            DN2 = WORK(INDRV1+N)*(DN2/WORK(INDRV3+N-1))
         ENDIF
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO

         WORK(INDRV3+N)=DN2

 330     TAU = WORK(INDRV3+N)
         TMP3 = WORK(INDRV3+N-1)
         IF( TAU.GT.TMP3 ) THEN
            S = TAU
            TAU = TMP3
            TMP3 = S
         END IF
         T = HALF*( ( TMP3-TAU )
     $        +WORK(INDRV4+N-1) )
         IF(WORK(INDRV4+N-1) .GT. TAU*TOL2
     $        .AND. T .NE. ZERO) THEN
            S = TAU*( WORK(INDRV4+N-1) / T )
            IF( S.LE.T ) THEN
               S = TAU*( WORK(INDRV4+N-1) /
     $              ( T*( ONE+SQRT( ONE+S / T ) ) ) )
            ELSE
               S = TAU*( WORK(INDRV4+N-1) /
     $              ( T+SQRT( T )*SQRT( T+S ) ) )
            END IF
            T = TMP3 + ( S+WORK(INDRV4+N-1) )
            TAU = TAU
     $           *( TMP3 / T )
            TMP3 = T
         END IF
         TAU = MIN(TAU,WORK(INDRV3+N))
         IF (TAU .EQ. ZERO) GO TO 1410

         TAU2 = MINVAL(WORK(INDRV3+M:INDRV3+N-1))
         IF (TAU2 .LE. TAU) THEN
            IF (TAU2 .EQ. ZERO) GO TO 1410
            TAU1 = TAU2
            GO TO 140
         ELSE
            TAU1 = TAU
         ENDIF
         
         DN2 = WORK(INDRV3+M)-TAU
         IF (DN2 .LE. ZERO) GO TO 140
         DO J=M, N-3
            DN2 = DN2*(WORK(INDRV3+J+1)/(DN2+WORK(INDRV4+J)))-TAU
            IF (DN2 .LE. ZERO) GO TO 140
         ENDDO
         
         TMP1 = DN2+WORK(INDRV4+N-2)
         IF (SAFMIN*TMP1 .LT. WORK(INDRV3+N-1) .AND.
     $        SAFMIN*WORK(INDRV3+N-1) .LT. TMP1) THEN
            DN1 = DN2*(WORK(INDRV3+N-1)/TMP1)-TAU
         ELSE
            DN1 = WORK(INDRV3+N-1)*(DN2/TMP1)-TAU
         ENDIF
         IF(DN1 .LE. ZERO) GO TO 140
         
         TMP1 = DN1+WORK(INDRV4+N-1)
         IF (SAFMIN*TMP1 .LT. WORK(INDRV3+N) .AND.
     $        SAFMIN*WORK(INDRV3+N) .LT. TMP1) THEN
            DN0 = DN1*(WORK(INDRV3+N)/TMP1)-TAU
         ELSE
            DN0 = WORK(INDRV3+N)*(DN1/TMP1)-TAU
         ENDIF
         IF(DN0 .LT. ZERO .AND. DN0+TAU .EQ. TAU) DN0=ZERO
         IF(DN0 .LT. ZERO) THEN
            TAU2 = MAX(DN0+TAU,ZERO)
            IF (TAU2 .EQ. TAU) TAU2 = ZERO
            TAU3 = CONST*TAU
            IF (TAU3 .EQ. TAU) TAU3 = HALF*TAU
            TAU = MAX(TAU2,TAU3)
         ENDIF

         IF (TAU .GT. ZERO) GO TO 125
         
 140     TMP7=SCALE/WORK(INDRV3+M)
         TMP1=TMP7
         TMP5=TMP7**2
         TMP3=TMP5
         TMP2=TMP3
         QMIN=TMP3
         DO J = M+1, N
            TMP4=SCALE/WORK(INDRV3+J)
            TMP8=(WORK(INDRV4+J-1)*SCALE1)*TMP4
            TMP7=TMP8*TMP7+TMP4
            TMP1=TMP1+TMP7
            TMP6=TMP7**2
            TMP3=TMP8*(TMP3+TMP5)+TMP6
            TMP2=TMP2+TMP3
            QMIN=MAX(QMIN,TMP3)
            TMP5=TMP6
         ENDDO
         
         TAU = MAX(ZERO,SCALE/TMP1,SCALE/SQRT(TMP2))
         TMP8 = DBLE(N-M+1)
         TMP9 = SCALE*(TMP8/
     $        (TMP1+SQRT(TMP8-ONE)*SQRT(TMP8*TMP2-TMP1**2)))
         IF (TMP9 .EQ. TMP9) THEN
            TAU = MAX(TAU,TMP9)
         ENDIF

         TAU2 = TAU1
         TMP8 = CEILING((TMP1/TMP2)*TMP1)
         IF (TMP8 .GT. ONE) THEN
            TMP9 = SCALE*(TMP8/
     $           (TMP1+SQRT(TMP8*TMP2-TMP1**2)/SQRT(TMP8-ONE)))
            IF (TMP9 .EQ. TMP9) THEN
               TAU2 = MIN(TMP9,SCALE/SQRT(QMIN),TAU2)
            ELSE
               TAU2 = MIN(SCALE/SQRT(QMIN),TAU2)
            ENDIF
         ELSE
            TAU2 = MIN(SCALE/SQRT(QMIN),TAU2)
         ENDIF
         
         TAU = MIN(TAU,TAU2)
         IF (TAU2 .LT. TWO*TAU) GO TO 125
         
 150     WORK(INDRV7+M:INDRV7+N) = ONE/SQRT(WORK(INDRV3+M:INDRV3+N))
         WORK(INDRV8+M:INDRV8+N-1) = SQRT(WORK(INDRV4+M:INDRV4+N-1))
         
         WORK(INDRV5+N) = WORK(INDRV7+N)
         TMP3 = WORK(INDRV5+N)
         DO J = N-1,M,-1
            WORK(INDRV1+J) = WORK(INDRV8+J)*WORK(INDRV7+J)
            WORK(INDRV5+J) = WORK(INDRV7+J)+
     $           WORK(INDRV1+J)*WORK(INDRV5+J+1)
            TMP3 = MAX(TMP3,WORK(INDRV5+J))
         ENDDO
         TMP3 = ONE/TMP3
         WORK(INDRV5+M) = (WORK(INDRV5+M)*TMP3)*WORK(INDRV7+M)
         TMP1 = WORK(INDRV5+M)
         DO J = M+1,N
            WORK(INDRV2+J) = WORK(INDRV8+J-1)*WORK(INDRV7+J)
            WORK(INDRV5+J) = (WORK(INDRV5+J)*TMP3)*WORK(INDRV7+J)+
     $           WORK(INDRV2+J)*WORK(INDRV5+J-1)
            TMP1 = MAX(TMP1,WORK(INDRV5+J))
         ENDDO
         TAU = ZERO
         IF (TMP1 .GT. ZERO) THEN
            TAU = MAX(TAU,TMP3/TMP1)
         ENDIF
         
         TMP1 = ONE/TMP1
         WORK(INDRV5+N) = WORK(INDRV5+N)*TMP1
         WORK(INDRV6+N) = WORK(INDRV5+N)*WORK(INDRV7+N)
         TMP3 = WORK(INDRV6+N)
         DO J = N-1,M,-1
            WORK(INDRV5+J) = WORK(INDRV5+J)*TMP1
            WORK(INDRV6+J) = WORK(INDRV5+J)*WORK(INDRV7+J)+
     $           WORK(INDRV1+J)*WORK(INDRV6+J+1)
            TMP3 = MAX(TMP3,WORK(INDRV6+J))
         ENDDO
         TMP3 = ONE/TMP3
         TMP2 = (WORK(INDRV6+M)*TMP3)*WORK(INDRV7+M)
         TMP1 = TMP2/WORK(INDRV5+M)
         DO J = M+1,N
            TMP2 = (WORK(INDRV6+J)*TMP3)*WORK(INDRV7+J)+
     $           WORK(INDRV2+J)*TMP2
            TMP1 = MAX(TMP1,TMP2/WORK(INDRV5+J))
         ENDDO
         IF (TMP1 .GT. ZERO) THEN
            TAU = MAX(TAU,TMP3/TMP1)
         ENDIF
         
         TAU = MIN(TAU,TAU1)
         IF (TAU .GT. ZERO) GO TO 125
         
         TAU=SQRT(WORK(INDRV3+N))-HALF*WORK(INDRV8+N-1)
         IF (TAU .LE. ZERO) GO TO 1410
         DO J=N-1,M+1,-1
            TMP2=SQRT(WORK(INDRV3+J))
     $           -HALF*(WORK(INDRV8+J)+WORK(INDRV8+J-1))
            IF (TMP2 .LE. ZERO) GO TO 1410
            TAU=MIN(TAU,TMP2)
         ENDDO
         TMP2=SQRT(WORK(INDRV3+M))-HALF*WORK(INDRV8+M)
         IF (TMP2 .LE. ZERO) GO TO 1410
         TAU=MIN(TAU,TMP2)**2
         
 125     IF (TAU .EQ. ZERO) GO TO 1410

         EMAX = ZERO
         DN2 = ONE
         WORK(INDRV5+M) = WORK(INDRV3+M)-TAU*DN2
         QMIN = WORK(INDRV5+M)
         DO J=M, N-3
            TMP1=WORK(INDRV4+J)/WORK(INDRV5+J)
            WORK(INDRV6+J) = WORK(INDRV3+J)*TMP1
            EMAX = MAX(EMAX,WORK(INDRV6+J))
            DN2=DN2*TMP1+ONE
            WORK(INDRV5+J+1) = WORK(INDRV3+J+1)-TAU*DN2
            QMIN=MIN(QMIN,WORK(INDRV5+J+1))
         ENDDO
         IF (QMIN .LE. ZERO .OR. EMAX*SAFMIN .GT. ONE) THEN
            TAU3 = CONST*TAU
            IF (TAU3 .EQ. TAU) TAU3 = HALF*TAU
            TAU = TAU3
            GO TO 125
         ENDIF

         TMP1=WORK(INDRV4+N-2)/WORK(INDRV5+N-2)
         WORK(INDRV6+N-2) = WORK(INDRV3+N-2)*TMP1
         IF (WORK(INDRV6+N-2)*SAFMIN .GT. ONE) THEN
            TAU3 = CONST*TAU
            IF (TAU3 .EQ. TAU) TAU3 = HALF*TAU
            TAU = TAU3
            GO TO 125
         ENDIF
         DN2=DN2*TMP1+ONE
         WORK(INDRV5+N-1) = WORK(INDRV3+N-1)-TAU*DN2
         IF (WORK(INDRV5+N-1) .LE. ZERO) THEN
            TAU2 = MAX(WORK(INDRV5+N-1)+TAU,ZERO)
            IF (TAU2 .EQ. TAU) TAU2 = ZERO
            TAU3 = CONST*TAU
            IF (TAU3 .EQ. TAU) TAU3 = HALF*TAU
            TAU = MAX(TAU2,TAU3)
            GO TO 125
         ENDIF

         TMP1=WORK(INDRV4+N-1)/WORK(INDRV5+N-1)
         WORK(INDRV6+N-1) = WORK(INDRV3+N-1)*TMP1
         IF (WORK(INDRV6+N-1)*SAFMIN .GT. ONE) THEN
            TAU3 = CONST*TAU
            IF (TAU3 .EQ. TAU) TAU3 = HALF*TAU
            TAU = TAU3
            GO TO 125
         ENDIF
         DN2=DN2*TMP1+ONE
         WORK(INDRV5+N) = WORK(INDRV3+N)-TAU*DN2
         IF (WORK(INDRV5+N) .LT. ZERO .AND.
     $        TAU*DN2+WORK(INDRV5+N) .EQ. TAU*DN2)
     $        WORK(INDRV5+N)=ZERO
         IF (WORK(INDRV5+N) .LT. ZERO) THEN
            TAU2 = MAX(WORK(INDRV5+N)+TAU,ZERO)
            IF (TAU2 .EQ. TAU) TAU2 = ZERO
            TAU3 = CONST*TAU
            IF (TAU3 .EQ. TAU) TAU3 = HALF*TAU
            TAU = MAX(TAU2,TAU3)
            GO TO 125
         ENDIF
         
         IF (DISNAN(WORK(INDRV5+N))) THEN
            TAU3 = CONST*TAU
            IF (TAU3 .EQ. TAU) TAU3 = HALF*TAU
            TAU = TAU3
            GO TO 125
         ENDIF

         CALL DMDLVS2(SIGMA,DESIG,TAU,T,DESIG0)
         SIGMA = T
         SIGMA1 = TOL2*SIGMA
         DESIG = DESIG0

         CALL DCOPY(N-M+1,WORK(INDRV5+M),1,WORK(INDRV3+M),1)
         CALL DCOPY(N-M,WORK(INDRV6+M),1,WORK(INDRV4+M),1)

 1410    DN2 = WORK(INDRV3+M)
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         DO J=M, N-1
            TMP2 = WORK(INDRV4+J)/DN2
            TMP3 = ONE+TMP2
            WORK(INDRV1+J) = DN2*TMP3
            DN2 = WORK(INDRV3+J+1)/TMP3
            WORK(INDRV2+J) = TMP2*DN2
            IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         ENDDO
         WORK(INDRV1+N) = DN2
         
         IF (DISNAN(WORK(INDRV2+N-1))) GO TO 2360
         
         GO TO 2380

 2360    DN2 = WORK(INDRV3+M)
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         DO J=M, N-3
            WORK(INDRV1+J) = DN2+WORK(INDRV4+J)
            TMP3 = WORK(INDRV3+J+1)/WORK(INDRV1+J)
            WORK(INDRV2+J) = WORK(INDRV4+J)*TMP3
            DN2 = DN2*TMP3
            IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         ENDDO

         WORK(INDRV1+N-2) = DN2+WORK(INDRV4+N-2)
         IF (WORK(INDRV1+N-2) .EQ. ZERO) THEN
            WORK(INDRV2+N-2) = ZERO
            DN2 = WORK(INDRV3+N-1)
         ELSE IF (SAFMIN*WORK(INDRV1+N-2) .LT. WORK(INDRV3+N-1) .AND.  
     $        SAFMIN*WORK(INDRV3+N-1) .LT. WORK(INDRV1+N-2)) THEN
            TMP2=WORK(INDRV3+N-1)/WORK(INDRV1+N-2)
            WORK(INDRV2+N-2)=WORK(INDRV4+N-2)*TMP2
            DN2=DN2*TMP2
         ELSE
            WORK(INDRV2+N-2) = WORK(INDRV3+N-1)
     $           *(WORK(INDRV4+N-2)/WORK(INDRV1+N-2))
            DN2 = WORK(INDRV3+N-1)*(DN2/WORK(INDRV1+N-2))
         ENDIF
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO

         WORK(INDRV1+N-1) = DN2+WORK(INDRV4+N-1)
         IF (WORK(INDRV1+N-1) .EQ. ZERO) THEN
            WORK(INDRV2+N-1) = ZERO
            DN2 = WORK(INDRV3+N)
         ELSE IF (SAFMIN*WORK(INDRV1+N-1) .LT. WORK(INDRV3+N) .AND.  
     $        SAFMIN*WORK(INDRV3+N) .LT. WORK(INDRV1+N-1)) THEN
            TMP2=WORK(INDRV3+N)/WORK(INDRV1+N-1)
            WORK(INDRV2+N-1)=WORK(INDRV4+N-1)*TMP2
            DN2=DN2*TMP2
         ELSE
            WORK(INDRV2+N-1) = WORK(INDRV3+N)
     $           *(WORK(INDRV4+N-1)/WORK(INDRV1+N-1))
            DN2 = WORK(INDRV3+N)*(DN2/WORK(INDRV1+N-1))
         ENDIF
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO

         WORK(INDRV1+N)=DN2

         IF (DISNAN(DN2)) GO TO 2370

         GO TO 2380
      
 2370    DN2 = WORK(INDRV3+M)
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         DO J=M, N-3
            WORK(INDRV1+J) = DN2+WORK(INDRV4+J)
            IF (WORK(INDRV1+J) .EQ. ZERO) THEN
               WORK(INDRV2+J) = ZERO
               DN2 = WORK(INDRV3+J+1)
            ELSE IF (SAFMIN*WORK(INDRV1+J) .LT. WORK(INDRV3+J+1) .AND.  
     $           SAFMIN*WORK(INDRV3+J+1) .LT. WORK(INDRV1+J)) THEN
               TMP2=WORK(INDRV3+J+1)/WORK(INDRV1+J)
               WORK(INDRV2+J)=WORK(INDRV4+J)*TMP2
               DN2=DN2*TMP2
            ELSE
               WORK(INDRV2+J) = WORK(INDRV3+J+1)
     $              *(WORK(INDRV4+J)/WORK(INDRV1+J))
               DN2 = WORK(INDRV3+J+1)*(DN2/WORK(INDRV1+J))
            ENDIF
            IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO
         ENDDO

         WORK(INDRV1+N-2) = DN2+WORK(INDRV4+N-2)
         IF (WORK(INDRV1+N-2) .EQ. ZERO) THEN
            WORK(INDRV2+N-2) = ZERO
            DN2 = WORK(INDRV3+N-1)
         ELSE IF (SAFMIN*WORK(INDRV1+N-2) .LT. WORK(INDRV3+N-1) .AND.  
     $        SAFMIN*WORK(INDRV3+N-1) .LT. WORK(INDRV1+N-2)) THEN
            TMP2=WORK(INDRV3+N-1)/WORK(INDRV1+N-2)
            WORK(INDRV2+N-2)=WORK(INDRV4+N-2)*TMP2
            DN2=DN2*TMP2
         ELSE
            WORK(INDRV2+N-2) = WORK(INDRV3+N-1)
     $           *(WORK(INDRV4+N-2)/WORK(INDRV1+N-2))
            DN2 = WORK(INDRV3+N-1)*(DN2/WORK(INDRV1+N-2))
         ENDIF
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO

         WORK(INDRV1+N-1) = DN2+WORK(INDRV4+N-1)
         IF (WORK(INDRV1+N-1) .EQ. ZERO) THEN
            WORK(INDRV2+N-1) = ZERO
            DN2 = WORK(INDRV3+N)
         ELSE IF (SAFMIN*WORK(INDRV1+N-1) .LT. WORK(INDRV3+N) .AND.  
     $        SAFMIN*WORK(INDRV3+N) .LT. WORK(INDRV1+N-1)) THEN
            TMP2=WORK(INDRV3+N)/WORK(INDRV1+N-1)
            WORK(INDRV2+N-1)=WORK(INDRV4+N-1)*TMP2
            DN2=DN2*TMP2
         ELSE
            WORK(INDRV2+N-1) = WORK(INDRV3+N)
     $           *(WORK(INDRV4+N-1)/WORK(INDRV1+N-1))
            DN2 = WORK(INDRV3+N)*(DN2/WORK(INDRV1+N-1))
         ENDIF
         IF (SIGMA .EQ. DN2+SIGMA) DN2 = ZERO

         WORK(INDRV1+N)=DN2

 2380    DO J= M, N-3
            IF (WORK(INDRV2+J) .LE. SIGMA1 .OR.
     $           WORK(INDRV4+J) .LE. TOL2*WORK(INDRV1+J)) THEN
               WORK(INDRV2+J)=-SIGMA
               WORK(INDRV4+J)=-DESIG
               M=J+1
            ENDIF
         ENDDO
         IF (WORK(INDRV4+N-2) .LE. TOL2*WORK(INDRV1+N-2)) THEN
            WORK(INDRV2+N-2)=ZERO
         ENDIF
         IF (WORK(INDRV4+N-1) .LE. TOL2*WORK(INDRV1+N-1)) THEN
            WORK(INDRV2+N-1)=ZERO
         ENDIF
         
      ENDDO
      INFO=2
      RETURN
  
 700  N=M-1
      IF (N .LE. 0) THEN
         GO TO 4000
      ENDIF
  
      DO J=M-2,1,-1
         IF (WORK(INDRV2+J) .LE. ZERO) THEN
            M=J+1
            GO TO 3000
         ENDIF
      ENDDO
      M=1
      GO TO 3000
  
 4000 A(1:N0)=SQRT(WORK(INDRV1+1:INDRV1+N0))
      CALL DLASRT( 'D', N0, A, IINFO )
      IF (SIGMX .GT. ZERO) THEN
         CALL DLASCL( 'G', 0, 0, SCALE, SIGMX, N0, 1, A, N0, IINFO )
      ENDIF
      RETURN
      END SUBROUTINE DMDLVS

      SUBROUTINE DMDLVS2(SIGMA,DESIG,TAU,SIGMA2,DESIG2)
      IMPLICIT NONE
      DOUBLE PRECISION SIGMA,DESIG,TAU,T,SIGMA2,DESIG2
      IF( TAU.LT.SIGMA ) THEN
         DESIG2 = DESIG + TAU
         T = SIGMA + DESIG2
         DESIG2 = DESIG2 - ( T-SIGMA )
      ELSE
         T = SIGMA + TAU
         DESIG2 = SIGMA - ( T-TAU ) + DESIG
      END IF
      SIGMA2 = T
      RETURN
      END SUBROUTINE DMDLVS2