Hi,

I have a question regard the bidiagonalization portion of the SVD algorithm, specifially the following code:

for(l=k;l>=1;l--)
{
nm=l-1;
if((float)(fabs(rv1[l])+anorm)==anorm)
{
flag=0;
break;
}
if((float)(fabs(w[nm])+anorm)==anorm)
break;
}

Here, l ranges from k:1, so nm ranges from k-1:0. Also, w is a vector with indices from 1:n. So, it would seem that in the second if statement, we would be asking for w[0] which doesn't exists. Am I correct here? This seems like an error to me. I would appreciate any comments.

thanks,
wes

Hi Wes,

The comment at this point in the text says
"Note that rv1[1] is always zero."
This means that the first if statement in the segment you gave will always be true when the index is 1, so the break is triggered. The line you were worried about is never reached in this case.

Saul Teukolsky

Great, thanks.

wes

I'd like to note that I got a bug at exactly this point. To be precise, it wasn't NR's code; I had translated it from C to Java, and converted float to double. Tests showed that
rv1[1] printed as 0.0, but
Math.abs(rv1[1]) + anorm == anorm
was false. I fixed the problem by changing the if line to
if (l == 1 || [the == test]) {
flag = false;
break;
}

The answer to my post seems to be Teukolsky's post in the Official Bug Reports forum, 08-05-2005.