Hi could you please help me out with the QR decomposition code...I have modified it a bit so that I don't have to include the nrutil file... all the magnitude values are correct, however I get a sign error quite to often...I check for all the possible solutions, but with no success...could you guys please review my code and I have 2 examples where errors are prominent(Note: The R doesn't always result in error...But when I am running QR as a part of large matrix the error is quite prominent)...Also the error is not only associate with zero, when I consider a larger matrix all the most of the matrix elements result in error..Attached is the code and C and matlab example's discrepencies...thanks for all your help in advanced

...I have 2 examples where errors are prominent

Firstly, QR decomposition is useful for solving systems of linear equations if the coefficient matrix is non-singular. The matrix you posted is singular.

Due to roundoff error, the matrix is, apparently, not reported as singular since, in qrdcmp.h, tests for values equal to zero fail with very small (but not exactly zero) values.

Secondly, I don't see any "prominent" differences between the Matlab results and the NR results in the examples you posted. (See Footnote.)

If you your definition of "prominent" errors includes displayed values of 0.0000 not being equal to -0.0000, or if you want to see values closer to the calculated values, just print everything with scientific notation. For your 5x5 example, here's what I get from NR3 qrdcmp:

R matrix of the decomposition:
-2.924038e+01 -3.112134e+01 -3.300230e+01 -3.488326e+01 -3.676422e+01
0.000000e+00 -1.209127e+00 -2.418254e+00 -3.627381e+00 -4.836508e+00
0.000000e+00 0.000000e+00 -1.054998e-15 -2.818357e-15 -3.157377e-15
0.000000e+00 0.000000e+00 0.000000e+00 3.400170e-15 -1.144897e-18
0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 -9.872438e-19

The very small numbers on the last three rows of this rank-2 matrix represent round-off error from the true zero values, and your fixed point printout shows them as 0.0000 or -0.0000.

Maybe you can experiment with "improving" the decomposition algorithm (so that singular matrices can be detected) for some matrices by changing floating point comparisons in qrdcmp from

"something == 0.0"

to

"fabs(something) < some_small_floating_point_number

The problem with that is that it might take a non-singular matrix with an unfortunately large condition number and report that it is singular. I mean, in general, how would you select the "some_small_floating_point_number" to use in the comparisons?

Bottom line:

You can multiply Q times R and see whether the residuals are sufficiently small. This may be adequate to give you some confidence in your decomposition code. I'm not sure whether this is has any real value in trying to solve a system of equations whose coefficient matrix is singular.

Better yet (probably) for practical analysis of the matrix itself:

Run the matrix through svdcmp and you might get more useful information.


Regards,

Dave

Footnote:
If you would like to see more comparisons, please paste the numbers into your post. Don't use screen shots. In addition to taking more bandwidth, they make it impossible to copy/paste the data to example files in case anyone wants try to reproduce your problems.