SUBROUTINE SSMALLSIGMA(FLAG,METHOD,M,N,A,B,TAUA,TAUB,TAU,
     $     WORK2, INDRV3, INDRV4, INDRV5, INDRV6, INDRV7, RTMIN, RTMAX)

      IMPLICIT NONE
      
      REAL A(*), B(*), WORK2(*)
      INTEGER N, M, J, METHOD, FLAG
      INTEGER INDRV3, INDRV4, INDRV5, INDRV6, INDRV7
      REAL TMP1, TMP2, TMP3, TMP4, TMP5
      REAL TAUA(4), TAUB, TAU
      REAL C1, S1, T, RTMIN, RTMAX
      
      REAL ONE, ZERO, HALF, TWO, CONST
      PARAMETER (ONE = 1.0E0, ZERO = 0.0E0, HALF = 0.5E0, TWO = 2.0E0)
      PARAMETER (CONST = 0.75E0)
      EXTERNAL SFMA0
      REAL SFMA0

      METHOD = 1
      IF (TAUB .LE. TAU) GO TO 160
*         
      TAUB = TAU
*     
      TMP4 = A(M)-TAU
      TMP5 = A(M)+TAU
      IF (TMP4 .LE. ZERO) THEN
         GO TO 160
      ELSE
         TMP3 = SQRT(TMP4)*SQRT(TMP5)
      ENDIF
      DO J = M, N-2
         CALL SLARTG2(TMP3,B(J),C1,S1,T,RTMIN,RTMAX)
         TMP4 = SFMA0(C1,A(J+1),-TAU)
         TMP5 = SFMA0(C1,A(J+1),TAU)
         IF (TMP4 .LE. ZERO) THEN
            GO TO 160
         ELSE
            TMP3 = SQRT(TMP4)*SQRT(TMP5)
         ENDIF
      ENDDO
      CALL SLARTG2(TMP3,B(N-1),C1,S1,T,RTMIN,RTMAX)
      TMP4 = SFMA0(C1,A(N),-TAU)
      IF (TMP4 .LT. ZERO) THEN
         TMP2 = C1*A(N)
         IF (TMP2 .LT. TAU) THEN
            IF (FLAG .EQ. 0) THEN
               TAUA(METHOD+1) = TMP2
            ELSE
               TMP3 = CONST*TAU
               IF (TMP3 .GE. TAU) TMP3 = HALF*TAU
               TAUA(METHOD+1) = MAX(TMP2,TMP3)
            ENDIF
         ELSE
            GO TO 160
         ENDIF
      ELSE
         TAUA(METHOD+1) = TAU
      ENDIF
      IF (TAUA(METHOD+1) .GT. TAUA(METHOD)) THEN
         METHOD = METHOD+1
      ENDIF
      IF (METHOD .GE. 2) RETURN
*
 160  WORK2(INDRV7+M:INDRV7+N) = ONE/A(M:N) 

      WORK2(INDRV5+N) = WORK2(INDRV7+N)
      TMP3 = WORK2(INDRV5+N)
      DO J = N-1,M,-1
         WORK2(INDRV3+J) = B(J)*WORK2(INDRV7+J)
         WORK2(INDRV5+J) = WORK2(INDRV7+J)+
     $        WORK2(INDRV3+J)*WORK2(INDRV5+J+1)
         TMP3 = MAX(TMP3,WORK2(INDRV5+J))
      ENDDO
      TMP3 = ONE/TMP3
      WORK2(INDRV5+M) = (WORK2(INDRV5+M)*TMP3)*WORK2(INDRV7+M)
      TMP1 = WORK2(INDRV5+M)
      DO J = M+1,N
         WORK2(INDRV4+J) = B(J-1)*WORK2(INDRV7+J)
         WORK2(INDRV5+J) = (WORK2(INDRV5+J)*TMP3)*WORK2(INDRV7+J)+
     $        WORK2(INDRV4+J)*WORK2(INDRV5+J-1)
         TMP1 = MAX(TMP1,WORK2(INDRV5+J))
      ENDDO
      IF (TMP1 .GT. ZERO) THEN
         TAUA(METHOD+1) = MIN(SQRT(TMP3/TMP1),TAUB)
         IF (TAUA(METHOD+1) .GT. TAUA(METHOD)) THEN
            METHOD = METHOD+1
         ENDIF
      ENDIF
      
      TMP1 = ONE/TMP1
      WORK2(INDRV5+N) = WORK2(INDRV5+N)*TMP1
      WORK2(INDRV6+N) = WORK2(INDRV5+N)*WORK2(INDRV7+N)
      TMP3 = WORK2(INDRV6+N)
      DO J = N-1,M,-1
         WORK2(INDRV5+J) = WORK2(INDRV5+J)*TMP1
         WORK2(INDRV6+J) = WORK2(INDRV5+J)*WORK2(INDRV7+J)+
     $        WORK2(INDRV3+J)*WORK2(INDRV6+J+1)
         TMP3 = MAX(TMP3,WORK2(INDRV6+J))
      ENDDO
      TMP3 = ONE/TMP3
      WORK2(INDRV6+M) = (WORK2(INDRV6+M)*TMP3)*WORK2(INDRV7+M)
      TMP1 = WORK2(INDRV6+M)/WORK2(INDRV5+M)
      DO J = M+1,N
         WORK2(INDRV6+J) = (WORK2(INDRV6+J)*TMP3)*WORK2(INDRV7+J)
     $        +WORK2(INDRV4+J)*WORK2(INDRV6+J-1)
         TMP1 = MAX(TMP1,WORK2(INDRV6+J)/WORK2(INDRV5+J))
      ENDDO
      IF (TMP1 .GT. ZERO) THEN
         TAUA(METHOD+1) = MIN(SQRT(TMP3/TMP1),TAUB)
         IF (TAUA(METHOD+1) .GT. TAUA(METHOD)) THEN
            METHOD = METHOD+1
         ENDIF
      ENDIF
      
      TMP1 = ONE/TMP1
      WORK2(INDRV6+N) = WORK2(INDRV6+N)*TMP1
      WORK2(INDRV5+N) = WORK2(INDRV6+N)*WORK2(INDRV7+N)
      TMP3 = WORK2(INDRV5+N)
      DO J = N-1,M,-1
         WORK2(INDRV6+J) = WORK2(INDRV6+J)*TMP1
         WORK2(INDRV5+J) = WORK2(INDRV6+J)*WORK2(INDRV7+J)+
     $        WORK2(INDRV3+J)*WORK2(INDRV5+J+1)
         TMP3 = MAX(TMP3,WORK2(INDRV5+J))
      ENDDO
      TMP3 = ONE/TMP3
      TMP2 = (WORK2(INDRV5+M)*TMP3)*WORK2(INDRV7+M)
      TMP1 = TMP2/WORK2(INDRV6+M)
      DO J = M+1,N
         TMP2 = (WORK2(INDRV5+J)*TMP3)*WORK2(INDRV7+J)
     $        +WORK2(INDRV4+J)*TMP2
         TMP1 = MAX(TMP1,TMP2/WORK2(INDRV6+J))
      ENDDO
      IF (TMP1 .GT. ZERO) THEN
         TAUA(METHOD+1) = MIN(SQRT(TMP3/TMP1),TAUB)
         IF (TAUA(METHOD+1) .GT. TAUA(METHOD)) THEN
            METHOD = METHOD+1
         ENDIF
      ENDIF

      IF (METHOD .GE. 2) RETURN
      
 170  WORK2(INDRV5+N) = ONE/A(N)
      TMP3 = WORK2(INDRV5+N)
      DO J = N-1,M,-1
         WORK2(INDRV3+J) = B(J)/A(J)
         WORK2(INDRV5+J) = ONE/A(J)+
     $        WORK2(INDRV3+J)*WORK2(INDRV5+J+1)
         TMP3 = MAX(TMP3,WORK2(INDRV5+J))
      ENDDO
      WORK2(INDRV5+M) = (WORK2(INDRV5+M)/TMP3)/A(M)
      TMP1 = WORK2(INDRV5+M)
      DO J = M+1,N
         WORK2(INDRV4+J) = B(J-1)/A(J)
         WORK2(INDRV5+J) = (WORK2(INDRV5+J)/TMP3)/A(J)+
     $        WORK2(INDRV4+J)*WORK2(INDRV5+J-1)
         TMP1 = MAX(TMP1,WORK2(INDRV5+J))
      ENDDO
      IF (TMP3 .GT. ZERO .AND. TMP1 .GT. ZERO) THEN
         TAUA(METHOD+1) = MIN(ONE/SQRT(TMP1)/SQRT(TMP3),TAUB)
         IF (TAUA(METHOD+1) .GT. TAUA(METHOD)) THEN
            METHOD = METHOD+1
         ENDIF
      ENDIF
      
      WORK2(INDRV5+N) = WORK2(INDRV5+N)/TMP1
      WORK2(INDRV6+N) = WORK2(INDRV5+N)/A(N)
      TMP3 = WORK2(INDRV6+N)
      DO J = N-1,M,-1
         WORK2(INDRV5+J) = WORK2(INDRV5+J)/TMP1
         WORK2(INDRV6+J) = WORK2(INDRV5+J)/A(J)+
     $        WORK2(INDRV3+J)*WORK2(INDRV6+J+1)
         TMP3 = MAX(TMP3,WORK2(INDRV6+J))
      ENDDO
      WORK2(INDRV6+M) = (WORK2(INDRV6+M)/TMP3)/A(M)
      TMP1 = WORK2(INDRV6+M)/WORK2(INDRV5+M)
      DO J = M+1,N
         WORK2(INDRV6+J) = (WORK2(INDRV6+J)/TMP3)/A(J)+
     $        WORK2(INDRV4+J)*WORK2(INDRV6+J-1)
         TMP1 = MAX(TMP1,WORK2(INDRV6+J)/WORK2(INDRV5+J))
      ENDDO
      IF (TMP3 .GT. ZERO .AND. TMP1 .GT. ZERO) THEN
         TAUA(METHOD+1) = MIN(ONE/SQRT(TMP1)/SQRT(TMP3),TAUB)
         IF (TAUA(METHOD+1) .GT. TAUA(METHOD)) THEN
            METHOD = METHOD+1
         ENDIF
      ENDIF

      WORK2(INDRV6+N) = WORK2(INDRV6+N)/TMP1
      WORK2(INDRV5+N) = WORK2(INDRV6+N)/A(N)
      TMP3 = WORK2(INDRV5+N)
      DO J = N-1,M,-1
         WORK2(INDRV6+J) = WORK2(INDRV6+J)/TMP1
         WORK2(INDRV5+J) = WORK2(INDRV6+J)/A(J)+
     $        WORK2(INDRV3+J)*WORK2(INDRV5+J+1)
         TMP3 = MAX(TMP3,WORK2(INDRV5+J))
      ENDDO
      TMP2 = (WORK2(INDRV5+M)/TMP3)/A(M)
      TMP1 = TMP2/WORK2(INDRV6+M)
      DO J = M+1,N
         TMP2 = (WORK2(INDRV5+J)/TMP3)/A(J)+
     $        WORK2(INDRV4+J)*TMP2
         TMP1 = MAX(TMP1,TMP2/WORK2(INDRV6+J))
      ENDDO
      IF (TMP3 .GT. ZERO .AND. TMP1 .GT. ZERO) THEN
         TAUA(METHOD+1) = MIN(ONE/SQRT(TMP1)/SQRT(TMP3),TAUB)
         IF (TAUA(METHOD+1) .GT. TAUA(METHOD)) THEN
            METHOD = METHOD+1
         ENDIF
      ENDIF

      IF (METHOD .GE. 2) RETURN

      TAUA(METHOD+1) = A(N)-HALF*B(N-1)
      IF (TAUA(METHOD+1) .LE. ZERO) RETURN
      DO J=N-1,M+1,-1
         TMP2 = A(J)-HALF*(B(J)+B(J-1))
         IF (TMP2 .LE. ZERO) RETURN
         TAUA(METHOD+1)=MIN(TAUA(METHOD+1),TMP2)
      ENDDO
      TMP2=A(M)-HALF*B(M)
      IF (TMP2 .LE. ZERO) RETURN
      TAUA(METHOD+1)=MIN(TAUA(METHOD+1),TMP2)
      IF (TAUA(METHOD+1) .GT. TAUA(METHOD)) THEN
         METHOD = METHOD+1
      ENDIF
      
      RETURN
      END SUBROUTINE SSMALLSIGMA

      SUBROUTINE SOQDQ(N, A, B, WORK2, INDRV5, INDRV6, SIGMA1,
     $     EPS, RTMIN, RTMAX, FLAG)

      IMPLICIT NONE
      INTEGER N, INDRV5, INDRV6, J, FLAG
      REAL A(*), B(*), WORK2(*)
      REAL TMP1, C1, S1
      REAL EPS, RTMIN, RTMAX, SIGMA1
      REAL ONE, ZERO, HALF, TWO
      PARAMETER (ONE = 1.0E0, ZERO = 0.0E0, HALF = 0.5E0, TWO = 2.0E0)

      FLAG = 0
      
      TMP1 = A(1)
      IF (TMP1 .LE. SIGMA1) THEN
         TMP1 = ZERO
         FLAG = 1
      ENDIF
      DO J = 1, N-1
         IF (B(J) .LE. EPS*TMP1) B(J) = ZERO
         IF (TMP1 .EQ. ZERO .AND. B(J) .EQ. ZERO) THEN
            C1 = ZERO
            S1 = ONE
            WORK2(INDRV5+J) = ZERO
         ELSE
            CALL SLARTG2(TMP1,B(J),C1,S1,WORK2(INDRV5+J),RTMIN,RTMAX)
         ENDIF
         WORK2(INDRV6+J) = S1*A(J+1)
         TMP1 = C1*A(J+1)
         IF (TMP1 .LE. SIGMA1) THEN
            TMP1 = ZERO
            FLAG = 1
         ENDIF
      ENDDO
      WORK2(INDRV5+N) = TMP1

      END SUBROUTINE SOQDQ
      
      SUBROUTINE SOQDSQ(L,N0,M0,A,B,SU,LDSU,WORK2,WORK3,
     $     IWORK,IWORK2,INFO)

      IMPLICIT NONE
      
      INTEGER N0, N1, M0, INFO, LDSU, FLAG, OLDN, IWORK(N0),
     $     IWORK2(M0)
      REAL A(*), B(*), WORK2(*), WORK3(2*N0,M0)
      REAL SU(LDSU, *)
      INTEGER N, I, I0, J1, J, L
      INTEGER INDRV3, INDRV4, INDRV5, INDRV6, INDRV7, METHOD
      REAL TMP1, TMP2, TMP3, TMP4, TMP5
      REAL TAU, TAUA(4), TAUB
      REAL SIGMA, SIGMA1, DESIG, SIGMA2, T, DESIG0
      REAL C1, S1, DMIN1
      REAL EPS, TOL, SAFMIN, SAFMAX, RTMIN, RTMAX
      
      REAL ONE, ZERO, HALF, TWO, CONST
      PARAMETER (ONE = 1.0E0, ZERO = 0.0E0, HALF = 0.5E0, TWO = 2.0E0)
      PARAMETER (CONST = 0.75E0)
      REAL TEN
      PARAMETER (TEN = 10.0E0)
*     
      EXTERNAL SLAMCH
      REAL SLAMCH
      EXTERNAL SFMA0
      REAL SFMA0
      LOGICAL LSAME
      EXTERNAL LSAME
*
!$omp parallel
!$omp single
      TAUA(1) = ZERO
      
      INDRV3 = 0
      INDRV4 = INDRV3+N0
      INDRV5 = INDRV4+N0
      INDRV6 = INDRV5+N0
      INDRV7 = INDRV6+N0
*     
      EPS = SLAMCH( 'Precision' )
      SAFMIN = SLAMCH( 'Safe minimum' )
      RTMIN = SQRT( SAFMIN )
      SAFMAX = ONE/SAFMIN
      RTMAX = SQRT( SAFMAX/2.0E0 )
      TOL = TEN*EPS
*     
      DO I = 1, N0-1
         IF (A(I) .LT. ZERO) THEN
            A(I) = -A(I)
            B(I) = -B(I)
!            CALL SSCAL(N0,-ONE,SU(1,I),1)
            IWORK(I) = -1
         ELSE
            IWORK(I) = 1
         ENDIF
         IF (B(I) .LT. ZERO) THEN
            B(I) = -B(I)
            A(I+1) = -A(I+1)
         ENDIF
      ENDDO
      IF (A(N0) .LT. ZERO) THEN
         A(N0) = -A(N0)
!         CALL SSCAL(N0,-ONE,SU(1,N0),1)
         IWORK(N0) = -1
      ELSE
         IWORK(N0) = 1
      ENDIF
*
      OLDN = -1
      N = N0
*
 3000 SIGMA = ZERO
      SIGMA1 = SAFMIN
      DESIG = ZERO
      SIGMA2 = TOL*SIGMA
      
      DO I = 1, M0+1
*
 15      IF (N .LE. N0-L) GO TO 4000
*         
         IF (N .EQ. 1) THEN
            
            CALL SLARTG7(SIGMA,DESIG,A(N),A(N),DESIG0,RTMIN)

!$omp task firstprivate(N) private(N1,J1,J,I0,TMP1,TMP2,TMP3,TMP4,C1,S1)
            N1 = N0+1-N
            DO J = 1,N0
               SU(J,N1)=ZERO
            ENDDO
            SU(N,N1)=ONE
            DO I0=I-1,1,-1
               J1 = IWORK2(I0)
               TMP2 = SU( J1, N1)
               DO J = J1-1,1,-1
                  TMP4 = WORK3(2*J-1,I0)
                  TMP3 = WORK3(2*J,I0)
                  TMP1 = SU( J, N1)
                  IF (SIGN(ONE,TMP4) .GT. ZERO) THEN
                     SU( J+1, N1) = TMP4*(TMP3*TMP2+TMP1)+TMP2
                     TMP2 = TMP4*(TMP3*TMP1-TMP2)+TMP1
                  ELSE
                     TMP2 = -TMP2
                     SU( J+1, N1) = TMP4*(TMP3*TMP1+TMP2)+TMP1
                     TMP2 = TMP4*(TMP3*TMP2-TMP1)+TMP2
                  ENDIF
               ENDDO
               SU( 1, N1) = TMP2
            ENDDO
            DO J = N0,1,-1
               IF (IWORK(J) .LT. 0) THEN
                  SU(J,N1) = - SU(J,N1)
               ENDIF
            ENDDO
!$omp end task
            
            GO TO 4000
            
         ENDIF
*
         IF (B(N-1) .LE. SIGMA2) THEN

            IF(A(N) .EQ. ZERO) THEN
               A(N) = SIGMA

!$omp task firstprivate(N) private(N1,J1,J,I0,TMP1,TMP2,TMP3,TMP4,C1,S1)
               N1 = N0+1-N
               DO J = 1,N0
                  SU(J,N1)=ZERO
               ENDDO
               SU(N,N1)=ONE
               DO I0=I-1,1,-1
                  J1 = IWORK2(I0)
                  TMP2 = SU( J1, N1)
                  DO J = J1-1,1,-1
                     TMP4 = WORK3(2*J-1,I0)
                     TMP3 = WORK3(2*J,I0)
                     TMP1 = SU( J, N1)
                     IF (SIGN(ONE,TMP4) .GT. ZERO) THEN
                        SU( J+1, N1) = TMP4*(TMP3*TMP2+TMP1)+TMP2
                        TMP2 = TMP4*(TMP3*TMP1-TMP2)+TMP1
                     ELSE
                        TMP2 = -TMP2
                        SU( J+1, N1) = TMP4*(TMP3*TMP1+TMP2)+TMP1
                        TMP2 = TMP4*(TMP3*TMP2-TMP1)+TMP2
                     ENDIF
                  ENDDO
                  SU( 1, N1) = TMP2
               ENDDO
               DO J = N0,1,-1
                  IF (IWORK(J) .LT. 0) THEN
                     SU(J,N1) = - SU(J,N1)
                  ENDIF
               ENDDO
!$omp end task
               
               N = N-1
               IF (N .LE. 0) THEN
                  GO TO 4000
               ENDIF

               GO TO 15
               
            ENDIF

            CALL SLARTG7(SIGMA,DESIG,A(N),T,DESIG0,RTMIN)
            IF(T .LE. SIGMA) THEN
               A(N) = T

!$omp task firstprivate(N) private(N1,J1,J,I0,TMP1,TMP2,TMP3,TMP4,C1,S1)
               N1 = N0+1-N
               DO J = 1,N0
                  SU(J,N1)=ZERO
               ENDDO
               SU(N,N1)=ONE
               DO I0=I-1,1,-1
                  J1 = IWORK2(I0)
                  TMP2 = SU( J1, N1)
                  DO J = J1-1,1,-1
                     TMP4 = WORK3(2*J-1,I0)
                     TMP3 = WORK3(2*J,I0)
                     TMP1 = SU( J, N1)
                     IF (SIGN(ONE,TMP4) .GT. ZERO) THEN
                        SU( J+1, N1) = TMP4*(TMP3*TMP2+TMP1)+TMP2
                        TMP2 = TMP4*(TMP3*TMP1-TMP2)+TMP1
                     ELSE
                        TMP2 = -TMP2
                        SU( J+1, N1) = TMP4*(TMP3*TMP1+TMP2)+TMP1
                        TMP2 = TMP4*(TMP3*TMP2-TMP1)+TMP2
                     ENDIF
                  ENDDO
                  SU( 1, N1) = TMP2
               ENDDO
               DO J = N0,1,-1
                  IF (IWORK(J) .LT. 0) THEN
                     SU(J,N1) = - SU(J,N1)
                  ENDIF
               ENDDO
!$omp end task
               
               N = N-1
               IF (N .LE. 0) THEN
                  GO TO 4000
               ENDIF

               GO TO 15
               
            ENDIF

            IF (N .EQ. 2) THEN
               TAU = A(N-1)
            ELSE
               TAUB = MINVAL(A(1:N-2))
               CALL SLAS2(A(N-2), B(N-2), A(N-1), TMP2, TMP3)
               TAU = MIN(TMP2,A(N-1))
               CALL SSMALLSIGMA(0,METHOD,1,N-1,A,B,TAUA,TAUB,TAU,
     $              WORK2,INDRV3, INDRV4, INDRV5, INDRV6, INDRV7,
     $              RTMIN, RTMAX)
               TAU = TAUA(METHOD)
            ENDIF
               
            IF(A(N) .LE. TAU) THEN
               A(N) = T

!$omp task firstprivate(N) private(N1,J1,J,I0,TMP1,TMP2,TMP3,TMP4,C1,S1)
               N1 = N0+1-N
               DO J = 1,N0
                  SU(J,N1)=ZERO
               ENDDO
               SU(N,N1)=ONE
               DO I0=I-1,1,-1
                  J1 = IWORK2(I0)
                  TMP2 = SU( J1, N1)
                  DO J = J1-1,1,-1
                     TMP4 = WORK3(2*J-1,I0)
                     TMP3 = WORK3(2*J,I0)
                     TMP1 = SU( J, N1)
                     IF (SIGN(ONE,TMP4) .GT. ZERO) THEN
                        SU( J+1, N1) = TMP4*(TMP3*TMP2+TMP1)+TMP2
                        TMP2 = TMP4*(TMP3*TMP1-TMP2)+TMP1
                     ELSE
                        TMP2 = -TMP2
                        SU( J+1, N1) = TMP4*(TMP3*TMP1+TMP2)+TMP1
                        TMP2 = TMP4*(TMP3*TMP2-TMP1)+TMP2
                     ENDIF
                  ENDDO
                  SU( 1, N1) = TMP2
               ENDDO
               DO J = N0,1,-1
                  IF (IWORK(J) .LT. 0) THEN
                     SU(J,N1) = - SU(J,N1)
                  ENDIF
               ENDDO
!$omp end task
                
               N = N-1
               IF (N .LE. 0) THEN
                  GO TO 4000
               ENDIF

               GO TO 15
               
            ENDIF
         ENDIF
*
         IF (I .EQ. M0+1) EXIT
*
         DO J = 1, N-1
            IF (B(J) .LE. SIGMA2) THEN
               B(J) = ZERO
            ENDIF
         ENDDO
*
         CALL SLAS2(A(N-1), B(N-1), A(N), TMP2, TMP3)
         TAU = MIN(TMP2,A(N))
         IF (TAU .LE. SIGMA1) THEN
            CALL SOQDQ(N, A, B, WORK2, INDRV5, INDRV6, SIGMA1,
     $           EPS, RTMIN, RTMAX, FLAG)
            IF (FLAG .EQ. 1) GO TO 360
         ENDIF

         IF (OLDN .EQ. N) THEN
            TAUB = DMIN1
         ELSE
            TAUB = MINVAL(A(1:N-1))
         ENDIF
         IF (TAUB .LE. SIGMA1) THEN
            CALL SOQDQ(N, A, B, WORK2, INDRV5, INDRV6, SIGMA1,
     $           EPS, RTMIN, RTMAX, FLAG)
            IF (FLAG .EQ. 1) GO TO 360
         ENDIF
         
         CALL SSMALLSIGMA(1,METHOD,1,N,A,B,TAUA,TAUB,TAU,
     $        WORK2, INDRV3, INDRV4, INDRV5, INDRV6, INDRV7,
     $        RTMIN, RTMAX)
         TAU = TAUA(METHOD)
*     
 125     IF (TAU .LE. ZERO) THEN
            CALL SOQDQ(N, A, B, WORK2, INDRV5, INDRV6, SIGMA1,
     $           EPS, RTMIN, RTMAX, FLAG)
            GO TO 360
         ENDIF
*     
         TMP4 = A(1)-TAU
         TMP5 = A(1)+TAU
         IF (TMP4 .LE. ZERO) THEN
            IF (TMP4 .EQ. ZERO .AND. B(1) .EQ. ZERO) THEN
               TMP3 = ZERO
            ELSE
               TMP2 = A(1)
               IF (TMP2 .GE. TAU) TMP2 = ZERO
               TMP3 = CONST*TAU
               IF (TMP3 .GE. TAU) TMP3 = HALF*TAU
               TMP2 = MAX(TMP2,TMP3)
               IF (METHOD .GE. 2) THEN
                  IF (TAUA(METHOD-1) .GE. TMP2) THEN
                     METHOD = METHOD-1
                     TAU = TAUA(METHOD)
                  ELSE
                     TAU = TMP2
                  ENDIF
               ELSE
                  TAU = TMP2
               ENDIF
               GO TO 125
            ENDIF
         ELSE
            TMP3 = SQRT(TMP4)*SQRT(TMP5)
         ENDIF
         
         DO J = 1, N-2
            CALL SLARTG2(TMP3,B(J),C1,S1,WORK2(INDRV5+J),RTMIN,RTMAX)
            WORK2(INDRV6+J) = S1*A(J+1)
            
            TMP4 = SFMA0(C1,A(J+1),-TAU)
            TMP5 = SFMA0(C1,A(J+1),TAU)
            IF (TMP4 .LE. ZERO) THEN
               IF (TMP4 .EQ. ZERO .AND. B(J+1) .EQ. ZERO) THEN
                  TMP3 = ZERO
               ELSE
                  TMP2 = C1*A(J+1)
                  IF (TMP2 .GE. TAU) TMP2 = ZERO
                  TMP3 = CONST*TAU
                  IF (TMP3 .GE. TAU) TMP3 = HALF*TAU
                  TMP2 = MAX(TMP2,TMP3)
                  IF (METHOD .GE. 2) THEN
                     IF (TAUA(METHOD-1) .GE. TMP2) THEN
                        METHOD = METHOD-1
                        TAU = TAUA(METHOD)
                     ELSE
                        TAU = TMP2
                     ENDIF
                  ELSE
                     TAU = TMP2
                  ENDIF
                  GO TO 125
               ENDIF
            ELSE
               TMP3 = SQRT(TMP4)*SQRT(TMP5)
            ENDIF
         ENDDO
         CALL SLARTG2(TMP3,B(N-1),C1,S1,WORK2(INDRV5+N-1),RTMIN,RTMAX)
         WORK2(INDRV6+N-1) = S1*A(N)
         
         TMP4 = SFMA0(C1,A(N),-TAU)
         TMP5 = SFMA0(C1,A(N),TAU)
         IF (TMP4 .LT. ZERO) THEN
            TMP2 = C1*A(N)
            IF (TMP2 .LT. TAU) THEN
               TMP3 = CONST*TAU
               IF (TMP3 .GE. TAU) TMP3 = HALF*TAU
               TMP2 = MAX(TMP2,TMP3)
               IF (METHOD .GE. 2) THEN
                  IF (TAUA(METHOD-1) .GE. TMP2) THEN
                     METHOD = METHOD-1
                     TAU = TAUA(METHOD)
                  ELSE
                     TAU = TMP2
                  ENDIF
               ELSE
                  TAU = TMP2
               ENDIF
               GO TO 125
            ENDIF
            TMP3 = ZERO
         ELSE
            TMP3 = SQRT(TMP4)*SQRT(TMP5)
         ENDIF
         
         IF (TMP3 .GE. A(N)) THEN
            TMP2 = TAU/A(N)
            IF (C1 .GE. TMP2) THEN
               TMP3 = MIN(A(N)*SQRT(C1-TMP2)*SQRT(C1+TMP2),A(N))
            ELSE
               TMP3 = A(N)
            ENDIF
         ENDIF
         WORK2(INDRV5+N) = TMP3
*
         CALL SLARTG7(SIGMA,DESIG,TAU,T,DESIG0,RTMIN)
*     
         SIGMA = T
         SIGMA1 = MAX(SQRT(EPS)*SIGMA,SAFMIN)
         DESIG = DESIG0
         SIGMA2 = TOL*SIGMA
*     
 360     OLDN = N

         TMP1 = WORK2(INDRV5+1)
         IF (TMP1 .LE. SAFMIN) TMP1 = ZERO
         DMIN1 = TMP1
         DO J = 1, N-2
            IF (WORK2(INDRV6+J) .LE. EPS*TMP1) WORK2(INDRV6+J) = ZERO
            CALL SLARTG2(TMP1,WORK2(INDRV6+J),C1,S1,A(J),RTMIN,RTMAX)
            IF (C1 .GE. S1) THEN
               WORK3(2*J-1,I)=SIGN(S1,ONE)
               WORK3(2*J,I)=-S1/(ONE+C1)
            ELSE
               WORK3(2*J-1,I)=SIGN(C1,-ONE)
               WORK3(2*J,I)=C1/(ONE+S1)
            ENDIF
            B(J) = S1*WORK2(INDRV5+J+1)
            TMP1 = C1*WORK2(INDRV5+J+1)
            IF (TMP1 .LE. SAFMIN) TMP1 = ZERO
            DMIN1 = MIN(DMIN1,TMP1)
         ENDDO
         IF (WORK2(INDRV6+N-1) .LE. EPS*TMP1) WORK2(INDRV6+N-1) = ZERO
         CALL SLARTG2(TMP1,WORK2(INDRV6+N-1),C1,S1,A(N-1),RTMIN,RTMAX)
         IF (C1 .GE. S1) THEN
            WORK3(2*N-3,I)=SIGN(S1,ONE)
            WORK3(2*N-2,I)=-S1/(ONE+C1)
         ELSE
            WORK3(2*N-3,I)=SIGN(C1,-ONE)
            WORK3(2*N-2,I)=C1/(ONE+S1)
         ENDIF
         B(N-1) = S1*WORK2(INDRV5+N)
         TMP1 = C1*WORK2(INDRV5+N)
         IF (TMP1 .LE. SAFMIN) TMP1 = ZERO
         A(N) = TMP1
*
 400     IWORK2(I)=N
      ENDDO
      INFO = 2
      GO TO 4100
*     
 4000 INFO = 0
 4100 CONTINUE
!$omp taskwait
!$omp end single
!$omp end parallel
      RETURN
*     
      END SUBROUTINE SOQDSQ