+- HP Forums (https://www.hpmuseum.org/forum)
+-- Forum: HP Software Libraries (/forum-10.html)
+--- Forum: HP-41C Software Library (/forum-11.html)
+--- Thread: SandMatrix Routine of the week: Pseudo-Inverse (/thread-2718.html)
SandMatrix Routine of the week: Pseudo-Inverse - Ángel Martin - 12-28-2014 11:37 AM
Even if the convergence is finicky, using the BenIsrael-Cohen algorithm it takes a relatively short program to calculate the Moore-Penrose matrix pseudo-inverse. Once the initial estimation is calculated it's a matter of a few iterations until the final result is obtained - although each iteration also takes its time.
For additional references see:
http://benisrael.net/COHEN-BI-ITER-GI.pdf
and JM Baillard's page:
http://hp41programs.yolasite.com/pseudoinverse.php
The code is listed below; lots of Alpha stuff as they are the input control strings that drive the matrix operations. Also note three auxiliary matrices are needed - mostly due to the M*M instruction that requires a result matrix different from the operands.
Code:
01 LBL "MPINV"
02 ASTO O1
03 "|-,#"
04 MAT=
05 DIM?
06 INT
07 I<>J
08 LASTX
09 +
10 "P"
11 MATDIM
12 "#"
13 TRNPS
14 XEQ 00
15 "P"
16 MTRACE
17 "#,X,#"
18 MAT/
19 FNRM
20 STO 00
21 DIM?
22 "Q"
23 MATDIM
24 LBL 01
25 XEQ 00
26 "#,P,Q"
27 M*M
28 "#"
29 MAT+
30 "#,Q,#"
31 MAT-
32 FNRM
33 X<> 00
34 RCL- (00)
35 ABS
36 VIEW X
37 E-7
38 X<=Y?
39 GTO 01
40 "P"
41 PURFL
42 "Q"
43 PURFL
44 MNAME?
45 RTN
46 LBL 00
47 CLA
48 ARCL 01
49 "|-,#,P"
50 M*M
51 ENDYou'll need the SandMatrix for functions I<>J, MTRACE and MAT=
Auxiliary matrices "P" and "Q" are purged upon completion.
The result matrix "#" will not overwrite the initial one.
Cheers,
'AM
RE: SandMatrix Routine of the week: Pseudo-Inverse - Gene - 12-28-2014 08:48 PM
Thanks for this and all the work you have done on the modules for the 41. They make the 41CL a joy to use.
I'm sure the surface hasn't even been scratched with the abilities of the roms you have created. Many thanks!