The algorithm as it stands is prone to overflow, especially for large matrices.

For myself, I get overflow in this section of the accumuulation of the right-hand transformations:

for (k = l; k < n; ++k)
{
v[k][j] += s * v[k][i];
}

s: 9.3684343267160629e+305
v[k][i]: 444.83582005615625
v[k][j]: 1.#INF000000000000

Any suggestions for handling this? Can I do some sort of pre-processing normalisation on my data that I can undo again afterwards to get the right answer?

I have:

m: 2230
n: 651

I am calling the routine from svdfit, by the way, which zeroises small eigenvalues before calling svbksb, although that is not really relevant as overflow happens before control reaches that point in the code.

PS: I've found a way to reduce m and n significantly (now m = 1270 and n = 459) but I still get overflow.

Surprise, surprise... I'd made some mistakes in transcribing the code. One mistake was "for (k = 1 - 1; k < n; k++)" instead of "for (k = l - 1; k < n; k++)". The code is now working for me.

To the authors: please, please consider avoiding using "l" as a variable name. You acknowledge yourself in the FAQ that they are easily confused (especially in a font like Courier), so make things difficult for the user by using it? (Not everyone downloads code from the website.) Thanks.