diff --git a/SRC/dlaed2.f b/SRC/dlaed2.f
index 6d1357df9..9605c7af1 100644
--- a/SRC/dlaed2.f
+++ b/SRC/dlaed2.f
@@ -73,7 +73,7 @@
 *>         On entry, D contains the eigenvalues of the two submatrices to
 *>         be combined.
 *>         On exit, D contains the trailing (N-K) updated eigenvalues
-*>         (those which were deflated) sorted into increasing order.
+*>         (those which were deflated) sorted into decreasing order.
 *> \endverbatim
 *>
 *> \param[in,out] Q
diff --git a/SRC/dlaed8.f b/SRC/dlaed8.f
index 367f8d9f5..94ac0a9df 100644
--- a/SRC/dlaed8.f
+++ b/SRC/dlaed8.f
@@ -83,7 +83,7 @@
 *>          D is DOUBLE PRECISION array, dimension (N)
 *>         On entry, the eigenvalues of the two submatrices to be
 *>         combined.  On exit, the trailing (N-K) updated eigenvalues
-*>         (those which were deflated) sorted into increasing order.
+*>         (those which were deflated) sorted into decreasing order.
 *> \endverbatim
 *>
 *> \param[in,out] Q
diff --git a/SRC/dlasd2.f b/SRC/dlasd2.f
index 72803bd74..a88b9b457 100644
--- a/SRC/dlasd2.f
+++ b/SRC/dlasd2.f
@@ -86,7 +86,7 @@
 *>         On entry D contains the singular values of the two submatrices
 *>         to be combined.  On exit D contains the trailing (N-K) updated
 *>         singular values (those which were deflated) sorted into
-*>         increasing order.
+*>         decreasing order.
 *> \endverbatim
 *>
 *> \param[out] Z
@@ -217,7 +217,7 @@
 *>          IDXQ is INTEGER array, dimension(N)
 *>         This contains the permutation which separately sorts the two
 *>         sub-problems in D into ascending order.  Note that entries in
-*>         the first hlaf of this permutation must first be moved one
+*>         the first half of this permutation must first be moved one
 *>         position backward; and entries in the second half
 *>         must first have NL+1 added to their values.
 *> \endverbatim
@@ -451,7 +451,7 @@ SUBROUTINE DLASD2( NL, NR, SQRE, K, D, Z, ALPHA, BETA, U, LDU,
 *
 *        Check if singular values are close enough to allow deflation.
 *
-         IF( ABS( D( J )-D( JPREV ) ).LE.TOL ) THEN
+         IF( ( D( J )-D( JPREV ) ).LE.TOL ) THEN
 *
 *           Deflation is possible.
 *
@@ -487,7 +487,14 @@ SUBROUTINE DLASD2( NL, NR, SQRE, K, D, Z, ALPHA, BETA, U, LDU,
             END IF
             COLTYP( JPREV ) = 4
             K2 = K2 - 1
-            IDXP( K2 ) = JPREV
+*
+*           Insert the deflated index in the correct position in IDXP.
+*           If J - JPREV is greater than 1, the indices in between
+*           must be shifted to preserve the correct output order.
+*
+            DO 105 JP = JPREV, J - 1
+               IDXP( K2 + J - 1 - JP ) = JP
+  105       CONTINUE
             JPREV = J
          ELSE
             K = K + 1
diff --git a/SRC/dlasd7.f b/SRC/dlasd7.f
index a4ea5ae38..522b57d8a 100644
--- a/SRC/dlasd7.f
+++ b/SRC/dlasd7.f
@@ -99,7 +99,7 @@
 *>         On entry D contains the singular values of the two submatrices
 *>         to be combined. On exit D contains the trailing (N-K) updated
 *>         singular values (those which were deflated) sorted into
-*>         increasing order.
+*>         decreasing order.
 *> \endverbatim
 *>
 *> \param[out] Z
@@ -453,7 +453,7 @@ SUBROUTINE DLASD7( ICOMPQ, NL, NR, SQRE, K, D, Z, ZW, VF, VFW,
 *
 *        Check if singular values are close enough to allow deflation.
 *
-         IF( ABS( D( J )-D( JPREV ) ).LE.TOL ) THEN
+         IF( ( D( J )-D( JPREV ) ).LE.TOL ) THEN
 *
 *           Deflation is possible.
 *
@@ -489,7 +489,14 @@ SUBROUTINE DLASD7( ICOMPQ, NL, NR, SQRE, K, D, Z, ZW, VF, VFW,
             CALL DROT( 1, VF( JPREV ), 1, VF( J ), 1, C, S )
             CALL DROT( 1, VL( JPREV ), 1, VL( J ), 1, C, S )
             K2 = K2 - 1
-            IDXP( K2 ) = JPREV
+*
+*           Insert the deflated index in the correct position in IDXP.
+*           If J - JPREV is greater than 1, the indices in between
+*           must be shifted to preserve the correct output order.
+*
+            DO 85 JP = JPREV, J - 1
+               IDXP( K2 + J - 1 - JP ) = JP
+   85       CONTINUE
             JPREV = J
          ELSE
             K = K + 1
diff --git a/SRC/slaed2.f b/SRC/slaed2.f
index 9150bea15..24c62d5cf 100644
--- a/SRC/slaed2.f
+++ b/SRC/slaed2.f
@@ -73,7 +73,7 @@
 *>         On entry, D contains the eigenvalues of the two submatrices to
 *>         be combined.
 *>         On exit, D contains the trailing (N-K) updated eigenvalues
-*>         (those which were deflated) sorted into increasing order.
+*>         (those which were deflated) sorted into decreasing order.
 *> \endverbatim
 *>
 *> \param[in,out] Q
diff --git a/SRC/slaed8.f b/SRC/slaed8.f
index 47684d387..c51ccc934 100644
--- a/SRC/slaed8.f
+++ b/SRC/slaed8.f
@@ -83,7 +83,7 @@
 *>          D is REAL array, dimension (N)
 *>         On entry, the eigenvalues of the two submatrices to be
 *>         combined.  On exit, the trailing (N-K) updated eigenvalues
-*>         (those which were deflated) sorted into increasing order.
+*>         (those which were deflated) sorted into decreasing order.
 *> \endverbatim
 *>
 *> \param[in,out] Q
diff --git a/SRC/slasd2.f b/SRC/slasd2.f
index dfb203af4..116d42b2d 100644
--- a/SRC/slasd2.f
+++ b/SRC/slasd2.f
@@ -86,7 +86,7 @@
 *>         On entry D contains the singular values of the two submatrices
 *>         to be combined.  On exit D contains the trailing (N-K) updated
 *>         singular values (those which were deflated) sorted into
-*>         increasing order.
+*>         decreasing order.
 *> \endverbatim
 *>
 *> \param[out] Z
@@ -217,7 +217,7 @@
 *>          IDXQ is INTEGER array, dimension (N)
 *>         This contains the permutation which separately sorts the two
 *>         sub-problems in D into ascending order.  Note that entries in
-*>         the first hlaf of this permutation must first be moved one
+*>         the first half of this permutation must first be moved one
 *>         position backward; and entries in the second half
 *>         must first have NL+1 added to their values.
 *> \endverbatim
@@ -451,7 +451,7 @@ SUBROUTINE SLASD2( NL, NR, SQRE, K, D, Z, ALPHA, BETA, U, LDU,
 *
 *        Check if singular values are close enough to allow deflation.
 *
-         IF( ABS( D( J )-D( JPREV ) ).LE.TOL ) THEN
+         IF( ( D( J )-D( JPREV ) ).LE.TOL ) THEN
 *
 *           Deflation is possible.
 *
@@ -487,7 +487,14 @@ SUBROUTINE SLASD2( NL, NR, SQRE, K, D, Z, ALPHA, BETA, U, LDU,
             END IF
             COLTYP( JPREV ) = 4
             K2 = K2 - 1
-            IDXP( K2 ) = JPREV
+*
+*           Insert the deflated index in the correct position in IDXP.
+*           If J - JPREV is greater than 1, the indices in between
+*           must be shifted to preserve the correct output order.
+*
+            DO 105 JP = JPREV, J - 1
+               IDXP( K2 + J - 1 - JP ) = JP
+  105       CONTINUE
             JPREV = J
          ELSE
             K = K + 1
diff --git a/SRC/slasd7.f b/SRC/slasd7.f
index 0bff53e76..4c850dc75 100644
--- a/SRC/slasd7.f
+++ b/SRC/slasd7.f
@@ -99,7 +99,7 @@
 *>         On entry D contains the singular values of the two submatrices
 *>         to be combined. On exit D contains the trailing (N-K) updated
 *>         singular values (those which were deflated) sorted into
-*>         increasing order.
+*>         decreasing order.
 *> \endverbatim
 *>
 *> \param[out] Z
@@ -453,7 +453,7 @@ SUBROUTINE SLASD7( ICOMPQ, NL, NR, SQRE, K, D, Z, ZW, VF, VFW,
 *
 *        Check if singular values are close enough to allow deflation.
 *
-         IF( ABS( D( J )-D( JPREV ) ).LE.TOL ) THEN
+         IF( ( D( J )-D( JPREV ) ).LE.TOL ) THEN
 *
 *           Deflation is possible.
 *
@@ -489,7 +489,14 @@ SUBROUTINE SLASD7( ICOMPQ, NL, NR, SQRE, K, D, Z, ZW, VF, VFW,
             CALL SROT( 1, VF( JPREV ), 1, VF( J ), 1, C, S )
             CALL SROT( 1, VL( JPREV ), 1, VL( J ), 1, C, S )
             K2 = K2 - 1
-            IDXP( K2 ) = JPREV
+*
+*           Insert the deflated index in the correct position in IDXP.
+*           If J - JPREV is greater than 1, the indices in between
+*           must be shifted to preserve the correct output order.
+*
+            DO 85 JP = JPREV, J - 1
+               IDXP( K2 + J - 1 - JP ) = JP
+   85       CONTINUE
             JPREV = J
          ELSE
             K = K + 1