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

[Gene]

The CL Expanded Memory ROM found on TOS. QRG is there.

[rprosperi]

(05-31-2018 03:42 PM)Gene Wrote:  The CL Expanded Memory ROM found on TOS. QRG is there.

Thx. I have this!! I forgot the Y-Registers in the CL are in fact the Expanded memory register set. You know, the half-life of unused knowledge is very short...

[Ángel Martin]

(05-31-2018 05:39 PM)rprosperi Wrote:  

(05-31-2018 03:42 PM)Gene Wrote:  The CL Expanded Memory ROM found on TOS. QRG is there.

Thx. I have this!! I forgot the Y-Registers in the CL are in fact the Expanded memory register set. You know, the half-life of unused knowledge is very short...

Also available at Monte's "Other Docs" page:
http://systemyde.com/hp41/documents.html

I've finally decided to incorporate the enhancements in the SANDMATRIX module, which has multiple additions over the original Matrix function set. I'm actually amazed to see this working, believe me when i tell you it wasn't trivial...

So the SandMatrix now offers the following options to create a matrix, depending on its name in ALPHA:

1. 1st. character "R" => array will be located in the standard registers; R00 to the Rmax. defined by SIZE.
2. 1st. character "Y" => array will be located in the expanded CL Y-Registers; Y000 to Y1023
3. any other 1st. character => array will be located in X-Memory

The new option is obviously #2.

All functions in the SandMatrix are fully functional with the three types of arrays.
For example, with "R22,Y22" in ALPHA and CLST , calling MMOVE will copy the array from standard registers to expanded registers.

Ahh the sweet taste of success... ;-)

[MikeOShea]

Thanks for the puzzle Ángel. It seems to make sense that the determinate of such a matrix would boil down to a lot of 1*1 - 1+1 (or zero) terms with a few 0 - 1*1 terms (-1). So, -29 is related to the dimension. I think it's something like norm(30), or 435 terms, most of which are zero.

Anyway, I thought about how this would look as a function of the diagonal value from 0 to 2. Interesting:



Outside that range, the numbers get big/small very quickly!
Nice how concise such an analysis turns into Python. I'm thinking of writing a python to RPN compiler. Maybe someone here has already done so?

[Didier Lachieze]

(05-31-2018 08:09 PM)MikeOShea Wrote:  I'm thinking of writing a python to RPN compiler. Maybe someone here has already done so?

Yes, see here.

[MikeOShea]

Yes, see here.
[/quote]

Ah! Thank you.

[Valentin Albillo]

.
Hi, Ángel:

(05-31-2018 07:40 AM)Ángel Martin Wrote:  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! ;-)

You can do better ... ;-) ... Try this one out instead with your new machine-code functionality:

Consider the determinant D(N) = determinant of this NxN matrix:
 
      3    1    1 ... 1
      1    4    1 ... 1
      1    1    5 ... 1
      ... ... ... ... ...
      1    1    1 ... N+2

For example, D(3) = 50. Now compute and time D(11), D(13),D(30).

If you succeed in getting the correct integer values, see if you can find an exact formula for D(N)    (Hint: after a while, D(n) ends with more and more 0's).

Congratulations on your costly (in terms of your time and effort) achievement. Regards and have a nice weekend.
V.
.

[Gerson W. Barbosa]

(06-02-2018 11:41 PM)Valentin Albillo Wrote:  Consider the determinant D(N) = determinant of this NxN matrix:
 
      3    1    1 ... 1
      1    4    1 ... 1
      1    1    5 ... 1
      ... ... ... ... ...
      1    1    1 ... N+2

For example, D(3) = 50. Now compute and time D(11), D(13),D(30).

If you succeed in getting the correct integer values, see if you can find an exact formula for D(N)    (Hint: after a while, D(n) ends with more and more 0's).

On my HP 50g:

D11 = 1486442880.
D13 = 283465647360.
D30 = 3.31153874629E34

« → n
« '(n+1)!*(Psi(n+2)-Psi(1))' EVAL
»
»

Probably not the formula you have in mind, but that’s what I could get so far.

D(n) = (n + 1)!*Hn+1

Regards,

Gerson.

[Gerson W. Barbosa]

(06-03-2018 03:56 AM)Gerson W. Barbosa Wrote:  D30 = 3.31153874629E34

D30 = 33115387462887740717065450291200000

« → n
« '(n+1)!' '(COMB((n+2)!+n+1,n+1)-1)/(n+2)!' DUP 'n+2' / EVAL FLOOR 'n+2' * - * EVAL
»
»

[Valentin Albillo]

.
Hi, Gerson:

(06-03-2018 03:56 AM)Gerson W. Barbosa Wrote:  On my HP 50g:

D11 = 1486442880.
D13 = 283465647360.
D30 = 3.31153874629E34

All Ok, of course.

Quote:Probably not the formula you have in mind, but that’s what I could get so far.

D(n) = (n + 1)!*Hn+1

Exactly the formula I had in mind:

      D(N) = (N+1)!*(1+1/2+1/3+...1/(N+1)) = (N+1)!*H(N+1)

where H(N) is the sum of the first N terms of the Harmonic series. For example:

      D(10) = FACT(11)*(1+1/2+1/3+1/4+1/5+1/6+1/7+1/8+1/9+1/10+1/11)
            = 120543840

The sum of the Harmonic series is thus:

      H(N) = D(N-1)/FACT(N)

which surely would be one of the most inefficient ways to compute the sum. :-D ... For example

      H(11) = D(10)/FACT(11) = 120543840/39916800 = 3.01987734488

           = 1+1/2+1/3+1/4+1/5+1/6+1/7+1/8+1/9+1/10+1/11 = 3.01987734488

This 2-line HP-71B program will create and display the NxN matrix, then compute and display its determinant for arbitrary N:

1 DESTROY ALL @ OPTION BASE 1 @ INPUT N @ DIM A(N,N) @ MAT A=CON
2 FOR I=1 TO N @ A(I,I)=I+2 @ NEXT I @ MAT DISP A; @ DISP @ DISP DET(A)

>RUN
? 11

3 1 1 1 1 1 1  1  1  1  1
1 4 1 1 1 1 1  1  1  1  1
1 1 5 1 1 1 1  1  1  1  1
1 1 1 6 1 1 1  1  1  1  1
1 1 1 1 7 1 1  1  1  1  1
1 1 1 1 1 8 1  1  1  1  1
1 1 1 1 1 1 9  1  1  1  1
1 1 1 1 1 1 1 10  1  1  1
1 1 1 1 1 1 1  1 11  1  1
1 1 1 1 1 1 1  1  1 12  1
1 1 1 1 1 1 1  1  1  1 13

1486442880

Still, the idea was that Ángel would use these determinants to check his new microcoded matrix handling functions to see what timings and accuracy they would achieve. Let's hope he sees these posts and obligues.

Are you listening, Ángel ? :-D

Thanks, Gerson, and regards.
V.
.

[Ángel Martin]

Yes, I'm listening :-D

in fact I've been busy fixing a few glitches in the MCODE as a result of the beta testing, which this test case has also contributed to...

Thanks for the example, it is indeed an interesting one (but of course it *had* to be given its source ;-). First off, there's a few remarks to be made:

• The execution times are pretty much independent from the actual element values, as there's no provision in the MCODE for special cases (such as 1 or 0) in the calculations. Therefore the D30 example takes about 11 seconds, same as the 30x30 anti-identity matrix took in the previous example.
  

• The accuracy is the same as using the Advantage module. Remember: I didn't write new code - only patched it to allow for CL Y-registers as well as the standard RAM registers or X-Mem registers allowed by the Advantage.
  

With those two remarks made, here's how I did it.

First, the program used to create the matrices and calculate the determinant:

01 LBL "DDN" - expects n in X-Reg
02 RCL X
03 E3
04 /
05 +
06 "Y" - matrix will use the Y-Registers, starting at Y_001
07 MATDIM
08 1
09 MCON - in the SandMatrix
10 CLX
11 MSIJA
12 2
13 LBL 00
14 1
15 +
16 MSC+
17 SF 25
18 J+
19 FS?C 25
20 GTO 00
21 MDET
22 END

And these are the results:

D(11) = 1,486,442,880.0 - in 1.8 seconds
D(13) = 2.834656472 E11 - in 2.01 seconds
D(30) = 3.311538747 E34 - in about 11 seconds

Not bad again, in fact surprisingly good. Using Gerson's values as the reference, these results are practically as accurate as those obtained using the equivalent formula D(n) = (n+1)!*H(n+1). - which can be easily programmed using the Harmonic function in the SandMath (step #4 in the trivial program below):

01 LBL "DN" - n is expected in X
02 1
03 +
04 HARM
05 LASTX
06 FACT
07 *
08 END

which returns respectively (in just a few milliseconds this time):

D(11) = 1,486,442,880.0
D(13) = 2.834656474 E11
D(30) = 3.311538746 E34

Thanks for the case study, it's helped with the beta testing and adds another interesting problem to the list.

Saludos,
Á.

[Ángel Martin]

ok, let's get crazy and pull all the stops as if we were really playing organ...

why stop at 1,024 registers if there are 3,072 available on the CL for the expanded RAM?

You guessed what's coming... ;-)

Hint: A 55x55 anti-Identity matrix determinant anyone?
or rather how about D(40) = ??

Let's try:

40, XEQ "DDN" = (.. and 79.6 seconds later ...)=> 1.439439902 E50

but wait, there's more - n=55 will require 3,025 registers !

55, XEQ "DDN"= (.. and 1 min 99 secs later ... ) = 3.278748200 E75

he, he... read my lips: it's *wor-king*!!

[Massimo Gnerucci]

Congrats Ángel!

Well done, as always...

(06-05-2018 01:40 PM)Ángel Martin Wrote:  55, XEQ "DDN"= (.. and 1 min 99 secs later ... ) = 3.278748200 E75

Does this translate to 2 min 39 secs, or is it just a typo?

[Gerson W. Barbosa]

(06-05-2018 01:40 PM)Ángel Martin Wrote:  40, XEQ "DDN" = (.. and 79.6 seconds later ...)=> 1.439439902 E50

On my HP 50g, in Exact mode:

40

« → n
  « '(n+1)!*((COMB((n+2)!+n+1,n+1)-1)/(n+2)!-(n+2)*FLOOR((COMB((n+2)!+n+1,n+1)-1)/((n+2)*(n+2)!)))' EVAL
  »
»

TEVAL

->

143943990158833766399457315931645601623572480000000

s: 54.2999

But for a moment I'd forgotten about the KISS principle:

40

« 1 + ! 0 1 LASTARG
  FOR i i INV +
  NEXT * EXPAND
»
TEVAL

->

143943990158833766399457315931645601623572480000000

s: 3.1465

[Valentin Albillo]

.
Hi, Ángel:

(06-05-2018 01:40 PM)Ángel Martin Wrote:  ok, let's get crazy and pull all the stops as if we were really playing organ... why stop at 1,024 registers if there are 3,072 available on the CL for the expanded RAM?

See ? I told you in post #47: "You can do better ... ;-)
And you did !

Quote:but wait, there's more - n=55 will require 3,025 registers !
55, XEQ "DDN"= (.. and 1 min 99 secs later ... ) = 3.278748200 E75
he, he... read my lips: it's *wor-king*!!

My most sincere congratulations, truly excellent speed and accuracy. Wish you'd apply your amazing microcoding talents to my beloved HP71B ... perhaps some day ... :-D

Best regards.
V.
.