41CL Quiz: Determinant of 30x30 anti-Identity matrix

[Valentin Albillo]

(05-19-2018 12:33 AM)Gene Wrote:  I should have said for these matrices as specified here. Some with only the 0 and 1 elements which **ought** to be an integer answer are getting round off with 0 and 1's.

Can we really have a pivot that isn't a 0 or 1 in this case?

But of course. For the "anti-identity" matrices, the first that gives a non-integer determinant is the 5x5 case (12-digit HP-71B/Math ROM, 15-digit internally.)

Using a short routine I wrote, for general (0,1) square matrices I've readily found 5x5 and 4x4 cases with non-integer determinants, some of them quite peculiar.

Ditto for (1,3) 2x2 matrices, and so on and so forth ...

Regards.
V.

[Ángel Martin]

Glad to see this thread growing a life of its own... FWIW I'm knee-deep wading in the Advantage/CCD MCODE investigating a way to adapt the Matrices/MSYS/DET code to use the CL Y-registers instead... conceptually quite simple but the code structures make it rather tricky. Stay tuned.

[Valentin Albillo]

.
Self-quote to add numerical examples to what I said in my previous post:

(05-19-2018 04:50 AM)Valentin Albillo Wrote:  For the "anti-identity" matrices, the first that gives a non-integer determinant is the 5x5 case (12-digit HP-71B/Math ROM, 15-digit internally.)

Anti-identity:

>MAT DISP A;

      0 1 1 1 1
      1 0 1 1 1
      1 1 0 1 1
      1 1 1 0 1
      1 1 1 1 0

>DET(A)

      3.99999999999

Up to 30x30, the maximum error is found for the 10x10 case:

      Det = -9.00000000003

Quote:I've readily found 5x5 and 4x4 cases with non-integer determinants, some of them quite peculiar.

5x5 case:

>MAT DISP A;

      1 0 1 1 1
      1 1 0 1 0
      0 1 0 1 1
      1 1 1 0 1
      0 1 1 1 0

>DET(A)

      5.00000000001

4x4 case:

This one is most peculiar:

>MAT DISP A;

      0 0 0 1
      0 0 1 0
      1 1 1 1
      1 1 1 0

>DET(A)

      1.E-24 (!!)

3x3 case:

>MAT DISP A;

      1 1 1
      0 0 1
      1 1 0

>DET(A)

      -.000000000001

2x2 case: (1-3 elements)

>MAT DISP A;

      1 2
      3 1

>DET(A)

      -5.00000000001

And this deserves a full facepalm or two:

>MAT DISP A;

      3 1
      3 1

>DET(A)

      .000000000003

Of course the exact determinant is 0 because both rows are identical. The fact that it does not detect it but instead goes on a spree using a full-drawn algorithm with pivot(s) and such for so ultra-simple a case it truly astounding.

Best regards and have a nice weekend.
V.
.

[Gene]

Thank you, Valentin. Food for me to chew on.

Ever trying to learn.

[Werner]

(05-19-2018 11:53 PM)Valentin Albillo Wrote:  >MAT DISP A;

      3 1
      3 1

>DET(A)

      .000000000003

Of course the exact determinant is 0 because both rows are identical. The fact that it does not detect it but instead goes on a spree using a full-drawn algorithm with pivot(s) and such for so ultra-simple a case it truly astounding.

Best regards and have a nice weekend.
V.
.

The venerable 42S solves 2x2 cases separately, and so does not fall into this trap.
(but Free42 does ;-))
That, however, does not solve the problem as the equivalent

3 1 0
3 1 0
0 0 1

yields the determinant 3.e-12 as well. Not much to be done about that, save for the 48G's (and successor's) approach to return an integer determinant for a matrix with integer coefficients (it's actually even better than that but that would require a longer explanation ;-)

Cheers, Werner

[Valentin Albillo]

.
Hi, Werner:

(05-20-2018 09:21 AM)Werner Wrote:  The venerable 42S solves 2x2 cases separately, and so does not fall into this trap. (but Free42 does ;-)) That, however, does not solve the problem as the equivalent

3 1 0
3 1 0
0 0 1

yields the determinant 3.e-12 as well. Not much to be done about that, save for the 48G's (and successor's) approach to return an integer determinant for a matrix with integer coefficients.

I think that something very simple could and should have been done, namely to use the expanded formula for 2x2 and 3x3 cases, which entails very, very few multiplications and additions and no divisions whatsoever. This would be both faster and exact as compared to using the general reduction algorithm on such simple cases.

HP did whatever needed to be done to ensure that 2^3 gives exactly 8 but didn't care that det(3,1,3,1) isn't 0 as it should.

Quote:(it's actually even better than that but that would require a longer explanation ;-)

As we say here, "prisa no tenemos y tarde no es", i.e. "It's not too late and we're not in a hurry" so go ahead with that longer explanation, please, I'm all ears ... :-D

Regards.
V.
.

[J-F Garnier]

(05-20-2018 09:21 AM)Werner Wrote:  ...
for the 48G's (and successor's) approach to return an integer determinant for a matrix with integer coefficients ...

For me, this clearly illustrates the move from a calculator-for-enginners to a calculator-for-students (math oriented).

For an engineer, it doesn't matter so much that the determinant is exactly 0 or some 1e-12. Moreover, engineers usually deal with real, floating point numbers, rarely with pure integers.

J-F

[Werner]

To illustrate what the 48G (and successors) do to calculate the determinant, take
the Hilbert matrix of order 5, multiplied by 2520 to make all elements integers.

[[ 2520 1260 840 630 504 ]
[ 1260 840 630 504 420 ]
[ 840 630 504 420 360 ]
[ 630 504 420 360 315 ]
[ 504 420 360 315 280 ]]

The adjustment of the determinant is controlled by flag -54, Tiny element.
When that flag is clear, the determinant is adjusted, when set it is not.

So, with the above matrix on the stack, DET will show

381024.000008 -54 Set
381024 -54 Clear

Now divide the matrix by 10, introducing a few fractional numbers.
The result is

3.81024000008 -54 Set
3.81024 -54 Clear

Divide it by 10^10 again (so the original matrix is divided by 10^11)

3.81024000008e-50 -54 Set
3.81024e-50 -54 Clear

It is clear that the algorithm does not just check for integers.
What it does is it determines the scaling factor to guarantee that the least
significant digit of the computed determinant, divided by this
scaling factor, has a nonnegative exponent, and therefore the
ratio can be rounded to integer to possibly achieve greater
accuracy. (These are not my words, but Paul McClellan's ;-)

This scaling factor is 10^(N*s), where s is the exponent of the least significant digit
of all elements of the NxN matrix.
Some provisions are being made:
- flag -54 set, then the scaling factor is set to 0, meaning no rounding is to be performed
- When all elements of the matrix are zero, the scaling factor is set to 1
- for n>79, the scaling factor is set to zero as even in extended precision we risk overflow

Cheers, Werner

[Werner]

(05-20-2018 03:39 PM)J-F Garnier Wrote:  For me, this clearly illustrates the move from a calculator-for-enginners to a calculator-for-students (math oriented).

For an engineer, it doesn't matter so much that the determinant is exactly 0 or some 1e-12. Moreover, engineers usually deal with real, floating point numbers, rarely with pure integers.

J-F

I've always wondered if the algorithm as I described above could deliver a wrongly rounded result, in which case I would have to agree with you.
Otherwise, all it does is delivering a more accurate result, and I have no problem whatsoever with that. For those who do, just set Flag -54.

Cheers, Werner

[J-F Garnier]

Thanks for your detailed explanations about the determinant adjustment on the 48G series.

So I changed my mind, I recognize that it is not just a trick to artificially round the determinant results in integer matrix problems, but is generally an improvement of the accuracy, in a way similar to the much older 2^3 accuracy fix mentioned by Valentin.

You quoted Paul McClellan, do you have any references to articles or forum posts where he may have described more on the subject?

Quote:I've always wondered if the algorithm as I described above could deliver a wrongly rounded result, in which case I would have to agree with you.

Yes, and the fact that the adjustment can be disabled with flag -54 may suggest that, in some situations, it may be better not to adjust the determinant.

J-F

[toml_12953]

(05-20-2018 03:39 PM)J-F Garnier Wrote:  For an engineer, it doesn't matter so much that the determinant is exactly 0 or some 1e-12. Moreover, engineers usually deal with real, floating point numbers, rarely with pure integers.

J-F

Good point. If we can get to the moon using slide rules with about three significant figures, twelve sig figs should be enough to do just about anything!

[ijabbott]

I suppose just returning the determinant of the anti-identity matrix "by inspection" would be considered cheating? E.g the determinant of the NxN anti-identity matrix is N-1 for odd N, or 1-N for even N.

[Gerson W. Barbosa]

(05-22-2018 11:17 PM)ijabbott Wrote:  I suppose just returning the determinant of the anti-identity matrix "by inspection" would be considered cheating? E.g the determinant of the NxN anti-identity matrix is N-1 for odd N, or 1-N for even N.

I would not consider that cheating, although I suggested so in post #13. BTW, there is a typo there: the RPL program returns -29, not 29 (I expected someone would notice it, but apparently no one did).

[brickviking]

(05-22-2018 11:40 PM)Gerson W. Barbosa Wrote:  

(05-22-2018 11:17 PM)ijabbott Wrote:  I suppose just returning the determinant of the anti-identity matrix "by inspection" would be considered cheating? E.g the determinant of the NxN anti-identity matrix is N-1 for odd N, or 1-N for even N.

I would not consider that cheating, although I suggested so in post #13. BTW, there is a typo there: the RPL program returns -29, not 29 (I expected someone would notice it, but apparently no one did).

I noticed it (the typo), but failed to notify anyone else of this fact. And I did note the values swapping signs on increment. I thought that was a bit peculiar, but given I know nothing much about matrix operations, I chalked it up to my inexperience.

And it seems that x49gp is faster than the real 50g by quite some margin, indicating about 3.1 seconds to cough up a -29 instead of the 16.5-ish seconds I'd been seeing. Clockspeed seems to affect this time on the 50g, but not by as much as I'd hoped. I can't raise the clockspeed past 101 MHz before the calculator starts either behaving erratically, or just plain locking up.

(Post 231)

[Gerson W. Barbosa]

(05-23-2018 01:30 AM)brickviking Wrote:  

(05-22-2018 11:40 PM)Gerson W. Barbosa Wrote:  BTW, there is a typo there: the RPL program returns -29, not 29 (I expected someone would notice it, but apparently no one did).

I noticed it, but failed to notify anyone else of this fact. And I did note the values swapping signs on increment.

I meant noticing the typo.
Yes, that property allows for a fast and compact program when the only goal is finding the determinant.

[Ángel Martin]

(05-19-2018 05:14 AM)Ángel Martin Wrote:  Glad to see this thread growing a life of its own... FWIW I'm knee-deep wading in the Advantage/CCD MCODE investigating a way to adapt the Matrices/MSYS/DET code to use the CL Y-registers instead... conceptually quite simple but the code structures make it rather tricky. Stay tuned.

Not easy but after a couple of days reverse-engineering the original (and largely convoluted!) MCODE progress has been made: behold the Y-Advantage Module, which uses the CL Y-registers in full functionality. Just use the control character "Y" (instead of the original "R") and the matrix will be created/used from the Expanded registers area, from YR000 to YR1023.

Testing is holding up well so far, including the major operations (MDET, MINV, MSYS).
The "Y" control character replaces the same code used before with "R". I don;t know if both can be enabled simultaneously, so further investigation is in order.

To be continued.

[Claudio L.]

This thread took an interesting turn.
For those interested, here's a couple of links that were very useful to me:

Fraction-Free Methods for Determinants

And (more advanced), very similar concepts but applied to factorization methods:

Roundoff-Error-Free Algorithm for Solving linear systems...

In both cases the full text is available for free to download, just point and click.

[pier4r]

Thanks for sharing!

[Ángel Martin]

(05-23-2018 11:31 AM)Ángel Martin Wrote:  .... To be continued.

The Y-Advantage module is finally ready for beta testing. I'll spare you the details but it's been tricky to modify the original code in all needed points to take advantage of the CL Y-Registers, but I *believe* I figured them all out...

BTW the benchmark with the 30x30 anti-identity matrix is now 16 seconds to calculate its determinant - down from the 11 minutes using the FOCAL program. And the result is -29.000000000 exact to 9 decimal places.

Not bad! ;-)

[rprosperi]

(05-31-2018 07:40 AM)Ángel Martin Wrote:  

(05-23-2018 11:31 AM)Ángel Martin Wrote:  .... To be continued.

The Y-Advantage module is finally ready for beta testing. I'll spare you the details but it's been tricky to modify the original code in all needed points to take advantage of the CL Y-Registers, but I *believe* I figured them all out...

BTW the benchmark with the 30x30 anti-identity matrix is now 16 seconds to calculate its determinant - down from the 11 minutes using the FOCAL program. And the result is -29.000000000 exact to 9 decimal places.

Not bad! ;-)

Amazing results Angel! It's truly impressive what you are able to accomplish once you put your mind to it. Well done!

Which manual provides the best description of using the Y-Registers?