Add failing RREF Case

[Aweptimum]

As laid out in #472, this adds a case where the calculated rref is different from the expected one

[markrogoyski]

Hi @Aweptimum,

Thank you for opening a ticket and adding the test case. I'll take a look.

Mark

[Aweptimum]

It seems the issue is that when working with the two 0.0006114 entries as pivots, the row operations on the two lower-right entries of -3.738869*10^-7 result in 4.356*10^-17. Since it's below the margin of error, both of the lower diagonal entries are zeroed out. It would give the correct result otherwise.

[markrogoyski]

Have you tried manually setting setError() to something appropriate and then computing it?

[Aweptimum]

Yes, if the error is set to 10^-20 it passes.

[markrogoyski]

I'm inclined to think this isn't an issue, and just a natural result of numerical computations that involve divisions and comparisons. As soon as you divide numbers, you never really know what they are anymore due to floating-point precision errors. How do you determine if a very small number is intended, or zero? The library uses an arbitrary error threshold that works for most scenarios, but the option is there to override that when you know you are working with really small numbers.

How about we add this test case as a special case setting the error threshold and then call this one closed? Thoughts?

[Aweptimum]

I think there's something that can be done to improve the stability because the scaled up matrix passes. It's like the error needs to account for the relative magnitudes of the numbers in the matrix. Trouble is I'm having a hard time finding anything on improving the stability of gaussian elimination. But it's hard as a user to come up against a problem like this when the original given matrix only contains numbers on the order of 10¹ to 10³

[Beakerboy]

If all the elements are rational numbers, the RationalNumber class could be used. It would eliminate floating point errors. You can load an ObjectMatrix with RationalNumber objects, and the rref math should still work.

edit: Nevermind, I thought the method was abstracted to work on math Objects as well a built-in numeric types. If the RowDivide() and rowAdd() methods are implemented in ObjectMatrix, and the gaussianEliminiation method was abstracted more, it would work.

[Aweptimum]

After staring at it harder last week, I realized that the rref might be right because it came from a result of A - λI. The λ in question is in the eigenvalue array twice so it should have 2 solutions. Which means there should be 2 zero rows, right?

So yeah you can probably just close this. There's an issue elsewhere