Skip to content

Commit 5b0bc5e

Browse files
xORBDB5/xUNBDB5: ensure xORBDB6 input has unit norm
Call sites may run into problems when this subroutine computes nonzero vectors that are very small in norm.
1 parent 7eb4057 commit 5b0bc5e

File tree

4 files changed

+32
-4
lines changed

4 files changed

+32
-4
lines changed

SRC/cunbdb5.f

Lines changed: 8 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -179,7 +179,7 @@ SUBROUTINE CUNBDB5( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
179179
REAL EPS, NORM, SCL, SSQ
180180
* ..
181181
* .. External Subroutines ..
182-
EXTERNAL CLASSQ, CUNBDB6, XERBLA
182+
EXTERNAL CLASSQ, CUNBDB6, CSCAL, XERBLA
183183
* ..
184184
* .. External Functions ..
185185
REAL SLAMCH, SCNRM2
@@ -227,6 +227,13 @@ SUBROUTINE CUNBDB5( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
227227
NORM = SCL * SQRT( SSQ )
228228
*
229229
IF( NORM .GT. N * EPS ) THEN
230+
* Scale vector to unit norm to avoid problems in the caller code.
231+
* Computing the reciprocal is undesirable but
232+
* * xLASCL cannot be used because of the vector increments and
233+
* * the round-off error has a negligible impact on
234+
* orthogonalization.
235+
CALL CSCAL( M1, ONE / NORM, X1, INCX1 )
236+
CALL CSCAL( M2, ONE / NORM, X2, INCX2 )
230237
CALL CUNBDB6( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
231238
$ LDQ2, WORK, LWORK, CHILDINFO )
232239
*

SRC/dorbdb5.f

Lines changed: 8 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -179,7 +179,7 @@ SUBROUTINE DORBDB5( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
179179
DOUBLE PRECISION EPS, NORM, SCL, SSQ
180180
* ..
181181
* .. External Subroutines ..
182-
EXTERNAL DLASSQ, DORBDB6, XERBLA
182+
EXTERNAL DLASSQ, DORBDB6, DSCAL, XERBLA
183183
* ..
184184
* .. External Functions ..
185185
DOUBLE PRECISION DLAMCH, DNRM2
@@ -227,6 +227,13 @@ SUBROUTINE DORBDB5( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
227227
NORM = SCL * SQRT( SSQ )
228228
*
229229
IF( NORM .GT. N * EPS ) THEN
230+
* Scale vector to unit norm to avoid problems in the caller code.
231+
* Computing the reciprocal is undesirable but
232+
* * xLASCL cannot be used because of the vector increments and
233+
* * the round-off error has a negligible impact on
234+
* orthogonalization.
235+
CALL DSCAL( M1, ONE / NORM, X1, INCX1 )
236+
CALL DSCAL( M2, ONE / NORM, X2, INCX2 )
230237
CALL DORBDB6( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
231238
$ LDQ2, WORK, LWORK, CHILDINFO )
232239
*

SRC/sorbdb5.f

Lines changed: 8 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -179,7 +179,7 @@ SUBROUTINE SORBDB5( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
179179
REAL EPS, NORM, SCL, SSQ
180180
* ..
181181
* .. External Subroutines ..
182-
EXTERNAL SLASSQ, SORBDB6, XERBLA
182+
EXTERNAL SLASSQ, SORBDB6, SSCAL, XERBLA
183183
* ..
184184
* .. External Functions ..
185185
REAL SLAMCH, SNRM2
@@ -227,6 +227,13 @@ SUBROUTINE SORBDB5( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
227227
NORM = SCL * SQRT( SSQ )
228228
*
229229
IF( NORM .GT. N * EPS ) THEN
230+
* Scale vector to unit norm to avoid problems in the caller code.
231+
* Computing the reciprocal is undesirable but
232+
* * xLASCL cannot be used because of the vector increments and
233+
* * the round-off error has a negligible impact on
234+
* orthogonalization.
235+
CALL SSCAL( M1, ONE / NORM, X1, INCX1 )
236+
CALL SSCAL( M2, ONE / NORM, X2, INCX2 )
230237
CALL SORBDB6( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
231238
$ LDQ2, WORK, LWORK, CHILDINFO )
232239
*

SRC/zunbdb5.f

Lines changed: 8 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -179,7 +179,7 @@ SUBROUTINE ZUNBDB5( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
179179
DOUBLE PRECISION EPS, NORM, SCL, SSQ
180180
* ..
181181
* .. External Subroutines ..
182-
EXTERNAL ZLASSQ, ZUNBDB6, XERBLA
182+
EXTERNAL ZLASSQ, ZUNBDB6, ZSCAL, XERBLA
183183
* ..
184184
* .. External Functions ..
185185
DOUBLE PRECISION DLAMCH, DZNRM2
@@ -227,6 +227,13 @@ SUBROUTINE ZUNBDB5( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
227227
NORM = SCL * SQRT( SSQ )
228228
*
229229
IF( NORM .GT. N * EPS ) THEN
230+
* Scale vector to unit norm to avoid problems in the caller code.
231+
* Computing the reciprocal is undesirable but
232+
* * xLASCL cannot be used because of the vector increments and
233+
* * the round-off error has a negligible impact on
234+
* orthogonalization.
235+
CALL ZSCAL( M1, ONE / NORM, X1, INCX1 )
236+
CALL ZSCAL( M2, ONE / NORM, X2, INCX2 )
230237
CALL ZUNBDB6( M1, M2, N, X1, INCX1, X2, INCX2, Q1, LDQ1, Q2,
231238
$ LDQ2, WORK, LWORK, CHILDINFO )
232239
*

0 commit comments

Comments
 (0)