dgesdd: Handle norm nan value

[jschueller]

Closes #469

[langou]

Thanks Julien. We probably want to propagate the fix to SGESDD, CGESDD, and ZGESDD then.

[jschueller]

ok, done

[thijssteel]

minor suggestion, return info = N+1 instead of 1, that lets users distinguish non-convergence from NAN entries. Otherwise, looks good.

[langou]

Right, same as @thijssteel. I like the direct return if there is a NaN in the input matrix A. I am not sure about INFO=1 though.
I was planning to write in the comments:
INFO = 1: NaN present in input matrix A; or
INFO > 0: DBDSDC did not converge, updating process failed.
But the overlap of INFO > 0 and INFO =1 was not great. (Though I could live with it.)

What about INFO=-4 ? That would be for "invalid input A".

That would be consistent with INFO = -2 for example when the input M is invalid. Or INFO = -3 when the input N is invalid, etc.

I do not think we return INFO = -4 (invalid input) for an invalid array yet, but I think we can start on this.

So: INFO = -4 ?

[jschueller]

yes!

[codecov]

Codecov Report

Merging #471 (22f2ee7) into master (adbc021) will decrease coverage by 0.00%.
The diff coverage is 33.33%.

@@            Coverage Diff             @@
##           master     #471      +/-   ##
==========================================
- Coverage   83.23%   83.23%   -0.01%     
==========================================
  Files        1820     1820              
  Lines      170835   170847      +12     
==========================================
+ Hits       142195   142199       +4     
- Misses      28640    28648       +8     

Impacted Files	Coverage Δ
SRC/cgesdd.f	95.70% <33.33%> (-0.30%)	⬇️
SRC/dgesdd.f	94.84% <33.33%> (-0.42%)	⬇️
SRC/sgesdd.f	94.84% <33.33%> (-0.42%)	⬇️
SRC/zgesdd.f	95.70% <33.33%> (-0.30%)	⬇️

---

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update adbc021...22f2ee7. Read the comment docs.

[langou]

So now, codecov is not happy because we do not have a test to check that INFO=-4 is triggered correctly. Oh well. So, during the tester, we should do a 4x4 matrix (for example), put a NaN in it, send it to xGESDD and make sure xGESDD returns INFO=-4. So, maybe we put this (adding a check for INFO=-4 to xGESDD) on our TODO list and accept this Pull Request from @jschueller. BTW: thanks @jschueller

[langou]

Of note, Jim Demmel made a nice presentation a year ago on "Correctness of Floating-Point Programs: Exception Handling and Reproducibility". The slides are https://people.eecs.berkeley.edu/~demmel/ExceptionHandling_Reproducibility_Feb20.pptx
Slides 18-21 are the most relevant to this discussion. And I think @jschueller patch is more or less what is suggested on slide 21

[martin-frbg]

This has now caused a regression in NumPy by way of OpenBLAS - numpy/numpy#18914 so there appear to be some codebases who rely on the previous behaviour of silently propagating the NaN

[jschueller]

it should probably added to the doc then