Unexpected result calculating the determinant of a singular matrix (42S)

[Dave Britten]

When calculating the determinant of the matrix [[-2,1,3][1,2,1][3,1,-2]] on my 42S, I would expect to get 0 as it is a singular matrix. But the 42S says it's 3.30000000001E-12. How does the 42S calculate determinants that would lead to that result? And how do I identify when it's giving me a suspicious result? I've done very little linear algebra, so there might be something simple and obvious going on here.

[Valentin Albillo]

.
Hi, Dave:

(10-21-2019 01:54 AM)Dave Britten Wrote:  When calculating the determinant of the matrix [[-2,1,3][1,2,1][3,1,-2]] on my 42S, I would expect to get 0 as it is a singular matrix. But the 42S says it's 3.30000000001E-12. How does the 42S calculate determinants that would lead to that result? And how do I identify when it's giving me a suspicious result? I've done very little linear algebra, so there might be something simple and obvious going on here.

The 42S uses LU-decomposition to compute determinants. This process involves divisions so inexact terms are produced and thus rounding errors do creep in and that's why you don't get an exact result sometimes, even if the matrix has all integer elements and it's as small as 2x2.

For an exact way to compute determinants download and have a look at my PDF paper:

Exact Determinants and Permanents

which includes a program and revealing examples.

Regards.
V.

[Dave Britten]

Okay, I kind of wondered if there were some tricks like that going on. I remember seeing mention of LU decompositions in the Advantage Module manual, so I'm assuming that gives a similar result. The 48G, on the other hand, does return a determinant of 0 for that particular matrix, as does the 3x3 linear system solver from the 32SII manual.

[Valentin Albillo]

(10-21-2019 02:30 AM)Dave Britten Wrote:  Okay, I kind of wondered if there were some tricks like that going on. I remember seeing mention of LU decompositions in the Advantage Module manual, so I'm assuming that gives a similar result. The 48G, on the other hand, does return a determinant of 0 for that particular matrix, as does the 3x3 linear system solver from the 32SII manual.

The 48G "cheats": it detects that all elements are integer and forces the result (computed using LU as well) to be the nearest integer, which succeeds sometimes and fails some others.

Try the 7x7 matrix I give as an example in the linked paper in the 48G and see if you get 1 as the determinant. You won't, the cheat doesn't work.

V.

[Thomas Okken]

The 48G gets exactly zero on [[-.2 .1 .3][.1 .2 .1][.3 .1 -.2]] as well. Wouldn't that defeat the cheat?

[Valentin Albillo]

.
Hi, Thomas:

(10-21-2019 03:11 AM)Thomas Okken Wrote:  The 48G gets exactly zero on [[-.2 .1 .3][.1 .2 .1][.3 .1 -.2]] as well. Wouldn't that defeat the cheat?

No. It checks that the values of all elements are integer, not their types. The floating-point constant 2. has the integer value 2.

Whether the cheat is activated or not depends on a system flag, have a look at Messages #9, 40 and 49 in this thread (PDF document) at my site:

https://albillo.hpcalc.org/threads/HP%20...chines.pdf

This is discussed in detail in other similar threads available there.

Best regards.
V.

[pinkman]

Read Thomas’question carefully

[Moggul]

The Free 42 on my phone gives 0. Different methodology or different precision?

[Thomas Okken]

(10-21-2019 08:41 AM)Moggul Wrote:  The Free 42 on my phone gives 0. Different methodology or different precision?

Free42 uses LU decomposition to calculate determinants as well. The implementation is based on the one from "Numerical Recipes in C," without the "TINY" fudge factor when encountering a zero pivot. If the algorithm were identical to the one in the HP-42S, I'd expect the same kinds of errors, just smaller because of the extra precision, but in actual fact, it returns exactly zero, just like the 48G.

(10-21-2019 04:23 AM)Valentin Albillo Wrote:  

(10-21-2019 03:11 AM)Thomas Okken Wrote:  The 48G gets exactly zero on [[-.2 .1 .3][.1 .2 .1][.3 .1 -.2]] as well. Wouldn't that defeat the cheat?

No. It checks that the values of all elements are integer, not their types. The floating-point constant 2. has the integer value 2.

I think you misread my post. I tried Dave's example, divided by 10. Not inexact numbers, but actual non-integer values.

[John Keith]

The HP 50 in approximate mode and with flag -54 (the "cheat" flag) set*, returns 0 for both Dave's and Thomas's matrices. I'm surprised that the 42s gives an inexact result, I assumed the internal code for matrix math was essentially the same as the 48/49 series.


* Note that setting the flag removes the "cheat". Also, according to the HP 50 AUR, the basis for setting tiny elements to 0 is that intermediate results are less than 1E-14. It does not say anything about integer values.

[ttw]

(10-21-2019 11:58 AM)John Keith Wrote:  The HP 50 in approximate mode and with flag -54 (the "cheat" flag) set*, returns 0 for both Dave's and Thomas's matrices. I'm surprised that the 42s gives an inexact result, I assumed the internal code for matrix math was essentially the same as the 48/49 series.


* Note that setting the flag removes the "cheat". Also, according to the HP 50 AUR, the basis for setting tiny elements to 0 is that intermediate results are less than 1E-14. It does not say anything about integer values.

Approximately true. However, I've been working on some number theory stuff (which I'll post if I get it in shape for public use) that use really big (greater than 2^64) integers. Some seemingly integer stuff ends up being floating point. Conversions of big numbers out of binary comes to mine. I've got most things to work.

Reducing some things with big integer multiple of things like Sqrt(2) or (Sqrt(5)-1)/2 and the like need careful handling. Generally FLOOR and CEIL work well. FXND can be a problem as I found some case (I can't reproduce it) where would convert a number to a floating point.

That's for theory. In practice, I can just keep the numerators and denominators of stuff separate and use really close rational approximations for the irrational numbers in the final step.

[Dave Britten]

(10-21-2019 11:58 AM)John Keith Wrote:  The HP 50 in approximate mode and with flag -54 (the "cheat" flag) set*, returns 0 for both Dave's and Thomas's matrices. I'm surprised that the 42s gives an inexact result, I assumed the internal code for matrix math was essentially the same as the 48/49 series.

It looks like the "cheat" was added with the 48G. I just tried calculating the determinant with my 48SX, and I get the same 3.30000000001E-12 as the 42S. I don't have a 28S/C handy, but I would expect they do the same as the 42S and 48SX.

[Thomas Okken]

But the 48G returns exactly zero regardless of whether flag -54 is set or clear...

(Free42 Binary does not return zero, so the fact that Free42 Decimal does is apparently just a coincidence, not the result of a better algorithm. There is no "tiny element is zero" cheat in play either way.)

(The HP-15C returns 2e-9 (or rather, "15C Scientific Calculator by Vicinno," on iOS, but that should be the same thing).)

[Claudio L.]

(10-21-2019 03:06 PM)Thomas Okken Wrote:  But the 48G returns exactly zero regardless of whether flag -54 is set or clear...

(Free42 Binary does not return zero, so the fact that Free42 Decimal does is apparently just a coincidence, not the result of a better algorithm. There is no "tiny element is zero" cheat in play either way.)

(The HP-15C returns 2e-9 (or rather, "15C Scientific Calculator by Vicinno," on iOS, but that should be the same thing).)

May or may not be a coincidence. There are division-free decomposition algorithms that preserve integers, so if an integer matrix is given, all operations remain integer up until the end. I don't know if this is the case for Free42 or the 48g, I know it is the case with newRPL.
My point is, you don't necessarily need a cheat to get a good result with integer matrices.

Here's one algorithm, for example (not the one I used, but same idea).

I think this one is the one used by newRPL.

[Thomas Okken]

(10-21-2019 04:46 PM)Claudio L. Wrote:  

(10-21-2019 03:06 PM)Thomas Okken Wrote:  (Free42 Binary does not return zero, so the fact that Free42 Decimal does is apparently just a coincidence, not the result of a better algorithm. There is no "tiny element is zero" cheat in play either way.)

May or may not be a coincidence. There are division-free decomposition algorithms that preserve integers, so if an integer matrix is given, all operations remain integer up until the end. I don't know if this is the case for Free42 or the 48g, I know it is the case with newRPL.

Free42 does not use a division-free algorithm, that's why I said getting zero for Dave's example was a coincidence. Free42 Decimal and Free42 Binary use the same algorithm, but one returns zero and the other does not. The only difference between the two is the floating-point system used.

[toml_12953]

(10-21-2019 01:54 AM)Dave Britten Wrote:  When calculating the determinant of the matrix [[-2,1,3][1,2,1][3,1,-2]] on my 42S, I would expect to get 0 as it is a singular matrix. But the 42S says it's 3.30000000001E-12. How does the 42S calculate determinants that would lead to that result? And how do I identify when it's giving me a suspicious result? I've done very little linear algebra, so there might be something simple and obvious going on here.

Prime gets 0 in Home screen

[John Keith]

(10-21-2019 12:48 PM)ttw Wrote:  However, I've been working on some number theory stuff (which I'll post if I get it in shape for public use) that use really big (greater than 2^64) integers. Some seemingly integer stuff ends up being floating point. Conversions of big numbers out of binary comes to mine. I've got most things to work.

Reducing some things with big integer multiple of things like Sqrt(2) or (Sqrt(5)-1)/2 and the like need careful handling. Generally FLOOR and CEIL work well. FXND can be a problem as I found some case (I can't reproduce it) where would convert a number to a floating point.

Conversion of binary numbers to reals will lose precision if the value is greater than 10^12. To convert large binary numbers to exact integers, you can try
->STR 3. OVER SIZE 1. - SUB OBJ-> in decimal mode.

I can't see how FXND would return reals as long as you are in exact mode (check flags -3 and -105). If you want floating-point numbers with more than 12 digits you will have to use LongFloat.

I am also interested in number theory and I would like to see what you come up with.

Apologies to others for taking this thread off-topic.

[ttw]

I can't reproduce the FXND today. I will post the following when I get a "useful" version of convert fraction to partial quotients (not to hard, still debating on input style: two integers or fraction). At times I need statistics for the results; I'll probably just add some post-processing, things like sum of PCs or max or alternating sum and difference.

I also have some programs that convert a list of partial quotients into the quadratic irrational which that as the repeated part. Going the other to is fun. I still have a bit of work to do as this one is slow if using fractions but very complex if dealing with numerator and denominator separately. Some care is needed as the integer part (first PC) has to be handled separately, then there's a non-repeating part. Getting it in general is a bit tedious.

I may just post a "good" version with arbitrary choices and let anyone who wants to use these modify them to suit.

I've been using integers which are larger than 2^64 which can cause some problems. However keeping everything in integers is pretty helpful.

I did find out a few funny things: IDIV2 is very slow compared to separate parts (even IQUOT and IREMAINDER), the matrix form for continued fractions is really slow too.

[Valentin Albillo]

(10-21-2019 11:37 AM)Thomas Okken Wrote:  

(10-21-2019 04:23 AM)Valentin Albillo Wrote:  No. It checks that the values of all elements are integer, not their types. The floating-point constant 2. has the integer value 2.

I think you misread my post. I tried Dave's example, divided by 10. Not inexact numbers, but actual non-integer values.

Yes, absolutely. I read your post very late at night (actually almost dawn) and in the smallish screen of a tablet and without my reading glasses the decimal point was absolutely invisible to me. I simply assumed that your message said "1.", say, where it actually said ".1". My mistake, sorry.

Best regards.
V.

[Werner]

Hi everyone.

The 48GX (and up) does not check whether the matrix elements are integer - it determines the least significant digit in the input (say it is of the order 10^s) and with Flag -54 clear it will round the result to 10^(s*n), with n the order of the matrix.

(10-21-2019 02:36 AM)Valentin Albillo Wrote:  Try the 7x7 matrix I give as an example in the linked paper in the 48G and see if you get 1 as the determinant. You won't, the cheat doesn't work.

? of course the cheat works. With Flag -54 clear, the 48GX returns 1 exactly, with Flag -54 set it returns .999945522778. The condition number is about 10^11, and the 48GX works with 15 digits internally, so we get 15-11=4 correct digts.

Also, Valentin, the 42S uses a*b-c*d when calculating the determinant of a 2x2 system, as you can see when you calculate the determinant of

1 2
3 1

and

1 2 0
3 1 0
0 0 1

The former returns -5 exactly, the latter -5.00000000001

Best Regards,
Werner