[VA] SRC #015b - HP-15C & clones: COMPLEX Matrix Inverse up to 8x8

[Valentin Albillo]

      
Hi, all,

This thread's got a number of new posts so some comments ...

First Werner "2-byte" Proudly Wrote:I have to live up to my reputation, don't I?

Indeed you did, shaving off 2 bytes from my original code once again. That's the third or fourth time in a row so I think you've earned yourself the "2-byte" moniker.

Some remarks on your 18-step code vs. my 20-step routine:

• Your use of [÷] at step 007 plus [1/x] at steps 011 and 012 means that matrix A gets processed three times: (1) converting it to an LU-decomposed form, (2) then converting that LU form to A's inverse, and (3) finally inverting A's inverse to get back A (possibly with some negligible rounding errors) for further use.
        
  However, I'm not sure if [÷] (which returns the result of X^-1Y,) does indeed compute the inverse (of A) before the multiplication or, as I suspect, it simply computes A's LU form and then uses special multiplication rules as in the case of MATRIX 5, which computes Y^TX without actually transposing Y.
   
• On the other hand, my original code only subjects matrix A to two inversion procedures instead of your three procedures detailed above, so I feel that your code's accuracy might be slightly degraded (whether the matrix is degenerate or not,) and I also wonder if your routine's speed is negatively affected. Have you checked any of this ?

      
Well, I have !   These are the results:

• As for speed, using J-F's 5x5 complex matrix featured in my OP (Worked example) and calling our respective inversion routines 50 times in a loop I get these average timings:
  
        - my 20-step routine  : 0.89" per call
        - your 18-step routine: 0.91" per call
  
  so mine is negligibly faster (~2.2%).
        
• And as for accuracy, I've checked it (again using J-F's handy 5x5 example) by inverting the complex matrix twice in succession so that the original matrix ought to be returned if the inversion procedures were exact, then comparing each returned matrix A', B' vs. their respective originals A, B (conveniently saved to C, D because the inversion is done in-place and replaces the originals.)
  
  The comparison consists simply in subtracting the returned matrices from the saved originals (which if the inversions were error-free would result in the zero matrix,) and then computing the Row Norm (MATRIX 7) and the Frobenius Norm (MATRIX 8) of the respective A' - A (C) and B' - B (D) subtractions. The results are as follows:
  
    GSB C, GSBC, RESULT E,                        Mine       Yours  
    RCL MATRIX A, RCL MATRIX C, -, MATRIX 7 ->  3.8753E-8  2.9381E-8
                            LASTX, MATRIX 8 ->  3.0716E-8  2.2004E-8 
    RCL MATRIX B, RCL MATRIX D, -, MATRIX 7 ->  4.3000E-8  4.4401E-8
                            LASTX, MATRIX 8 ->  3.1621E-8  3.6484E-8 
  
  so you see, the differences are again very minor. Your errors are slightly better for matrix A while mine are slightly better for matrix B, no big deal. And remember, these are the errors after inverting the complex matrix twice, the errors for a single inversion should be significantly smaller.
  

Then J-F Garnier Wrote:I'm amazed by these little pieces of code posted by Valentin (congrats for your 1000 - now 1001 - posts !) and Werner.

Thanks, J-F. I was missing your comments.

Quote:I really appreciate short, powerful routines that expand the possibility of a machine. And this is exactly the case here: the first thread from Valentin was providing a workaround for the 8x8 limit for real matrix inversion.

With the HP-15C original and CE editions not having any sort of mass storage, it's vital that the programs are as short as possible lest no one will bother to key them in, not even to just try them out, as it would be such a lengthy, error-prone chore. And something had to be done about that annoying 8x8 limit, so partitioned matrices immediately sprang to mind.

Quote:Then this second thread goes further by getting the utmost from the limited memory of the HP-15c in the case of complex matrices. In the case of the 15c CE/192, this makes a difference and allows to manage 7x7 matrices instead of 6x6. An hypothetical 15c CE/224 could even handle 8x8 matrices.

On the other hand, a not-so-hypothetical-but-actually-extant DM15 clone with firmware M1B (229 regs) can handle 8x8 complex matrices (and perhaps even systems of 8 complex equations in as many complex unknowns.) And if the CE/192 user is not afraid of manually typing in just 16 or 18 instructions, he/she can invert an 8x8 complex matrix in the CE/192 as well, it'll take just a little care and likely less than 2 min.

BTW, it's a real pity HP didn't make the index registers R0, R1 and RI allocatable because that would've left just enough program space for the complex inversion routine itself so no need for manually typing in anything.

It's also a pity that Moravia didn't allow for the max. 229 registers in the CE, as SwissMicros did. The 37 extra regs might come in handy for everything, 8x8 complex matrix handling in particular.

Perhaps this capability was intentionally left for a future HP-15C CE 2.0 version including

• that many registers (229, or at least 224,)
• a cable and a Connectivity Kit of sorts for mass storage, backups, outputting results, easy and convenient program/data editing, sharing progs/data with an ...
• ... also-included software emulator for Windows/iOS/Android,
• a patched firmware which corrects the bug in the "nut" emulation layer that fatally broke the embedded HP-16C and might also affect the HP-15C ROM code in as-yet-unknown ways,
• a printed color copy of the "Advanced Functions Handbook",
• oh, and a "Thank you again !" card.

That way they'd easily milk yet another 130 € apiece from each of us ...  

Best regards.
V.