xLASQ1/xLASQ2 return INFO = 0 with NaN output on a NaN bidiagonal (DBDSQR without vectors, DBDSDC('N'), DGESVD('N','N'))

[rmlarsen]

xLASQ1/xLASQ2 (the dqds path) return INFO = 0 with NaN in the output when the input bidiagonal contains a NaN, so xBDSQR without vectors, xBDSDC with COMPQ = 'N' and xGESVD with JOBU = JOBVT = 'N' report success on a NaN matrix. Every sibling path reports failure on the same input: DBDSQR with vectors, DSTEQR and DSTERF return INFO = N-1 (non-convergence) and DGESDD returns INFO = -4 (its own NaN check from #469).

Reference-LAPACK master f96546f, gfortran 13.3, -O2, reference BLAS. N = 30, D = 1, 2, ..., 30, E = 0.5, one entry replaced by a NaN:

Routine	NaN in D(1)	NaN in D(15)	NaN in D(30)
DLASQ1	INFO = 0, 30 NaN	INFO = 0, 16 NaN	XERBLA from DLASCL (#1387)
DBDSQR('U', N, 0, 0, 0, ...)	INFO = 0, 30 NaN	INFO = 0, 16 NaN	XERBLA from DLASCL (#1387)
DBDSQR('U', N, N, 0, 0, ...) (with VT)	INFO = 29	INFO = 29	INFO = 29
DSTEQR('N', ...)	INFO = 29	INFO = 29	INFO = 29
DSTERF	INFO = 29	INFO = 29	INFO = 29

Through the drivers, with the same bidiagonal stored as a dense 30x30 matrix and A(1,1) = NaN: DGESVD('N', 'N', ...) returns INFO = 0 with all 30 singular values NaN, DGESVD('A', 'A', ...) returns INFO = 29, DGESDD('N', ...) returns INFO = -4. DBDSDC('U', 'N', ...) returns INFO = 0 for a NaN in any of D(1..29) (it calls DLASDQ, which calls DBDSQR without vectors). A NaN in E behaves the same way (E(1): INFO = 0, 2 NaN; E(N-1): INFO = 0, 30 NaN). SLASQ1 is identical.

Mechanism. The deflation tests in xLASQ3 are written as "keep iterating while Z(...) .GT. TOL2*(...)":

      IF( Z( NN-5 ).GT.TOL2*( SIGMA+Z( NN-3 ) ) .AND.
     $    Z( NN-2*PP-4 ).GT.TOL2*Z( NN-7 ) )
     $   GO TO 30
   20 CONTINUE
      Z( 4*N0-3 ) = Z( 4*N0+PP-3 ) + SIGMA
      N0 = N0 - 1
      GO TO 10

Every comparison with a NaN operand is false, so a NaN q or e is deflated at label 20 (or 40) as a converged singular value, the NaN never reaches a test that could set INFO, and the recurrence continues on the rest of the array. The other routines fail their convergence tests for the same reason and hit the iteration cap instead, which is what produces INFO = N-1.

Reproducer (link against the reference library):

      program lasq1_nan
      implicit none
      integer, parameter :: n = 30
      double precision :: d(n), e(n-1), work(4*n), one, nan
      integer :: info, i
      one = 1d0
      nan = sqrt( -one )
      do i = 1, n
         d(i) = dble( i )
      end do
      e = 0.5d0
      d(1) = nan
      call dlasq1( n, d, e, work, info )
      print *, 'info =', info, ' NaN count =', count( d /= d )
      end program

prints info = 0 NaN count = 30.

Suggested fix. xLASQ1 already computes the scaling norm from D and E; a NaN scan there (DISNAN( xLANST( 'M', N, D, E ) ), since xLANST propagates NaN) returning INFO = 1 would cover every position and make xBDSQR report INFO = N-1 through its existing "try QR after dqds failed" path, matching the with-vectors result. PR #1387 adds that check for the D(N) position only (where the NaN reaches SIGMX through MAX and would otherwise abort in DLASCL); widening it to the full norm is a one-line change to that PR if that is preferred over a separate fix. Related: #469 (DGESDD), #1382 (xLALSD/xBDSDC with vectors).