Yes the problem comes with a zero V.
You are correct that the case m=k needs to be protected:
- Code: Select all
If (k==m) then
lastv = k
In the case
DIRECT = 'F' and STOREV = 'C'
V is lower trapezoidal and the upper triangular part of V may be unset.
Further, the '1's on the diagonal are not stored and may also be unset.
Consider, in ILADLR, the pass up the third column, with the leading
lower triangular part of V contains zeros up to V(4,3);
V(3,3) is then accessed, and if this is zero then V(2,3) is referenced.
What we are really interested of course is lastv>k so we needn't be looking
up as far as the leading lower triangle part anyway.