HP Forums
(49G) Function of a Matrix - Printable Version

+- HP Forums (https://www.hpmuseum.org/forum)
+-- Forum: HP Software Libraries (/forum-10.html)
+--- Forum: General Software Library (/forum-13.html)
+--- Thread: (49G) Function of a Matrix (/thread-3786.html)



(49G) Function of a Matrix - Gerald H - 05-05-2015 03:47 PM

MFUN applies a function to a square matrix, calculator should be in approx & complex modes. For example, to find the natural log of matrix

M1=[[ -9. 0. -1. -4. ]
[ 4. 9. 3. 3. ]
[ 0. 4. -9. 9. ]
[ 5. -1. -6. -4. ]]

place M1 on stack level Y & 'LN(X)' on stack level X. MFUN returns the complex matrix

M2=[[ (2.17927367373,1.99350072846) (-.185449038469,-.25062812471) (.12978172209,1.13653167793) (-1.21072226166,-.24561453301) ]
[ (-.349990278973,-.725022350755) (2.17006519946,9.11517059164E-2) (.316994726521,-.413348667032) (-8.29926914019E-3,8.93282974816E-2) ]
[ (.499972969813,1.87320392796) (.360844076572,-.235504096367) (1.81399293871,1.06794824449) (2.26427803714,-.230793047339) ]
[ (1.5738586231,8.93453089839E-2) (.311847366461,-1.12327258886E-2) (-1.58567642989,5.09374150115E-2) (2.97368489567,-.011008025242) ]]

to stack level X.

As a check, you could then put 'e^X' on the stack & actuate MFUN again, to return

[[ (-8.99999999989,-1.62595437076E-10) (-7.789692811E-12,2.87950512152E-11) (-.999999999977,-1.80997814311E-10) (-3.99999999999,-1.76112719874E-11) ]
[ (3.99999999998,4.21474277251E-11) (9.00000000017,7.44661562123E-11) (2.99999999996,3.96498859739E-11) (3.00000000004,1.34469661894E-12) ]
[ (2.98413096713E-11,-1.53814112532E-10) (3.99999999998,4.04329178333E-11) (-8.99999999998,-1.19988781132E-10) (8.99999999994,1.82967287016E-10) ]
[ (4.99999999998,6.94233268647E-11) (-1.00000000004,6.827569037E-12) (-5.99999999997,-5.4121217567E-11) (-4.00000000003,1.19154531002E-10) ]]

where the imaginary parts are v small, so RE 0 RND looks more like the original M1.

I am sure the programme has enormous inefficiencies & would welcome improvements.

TI fans - the TI-92 does this out of the box, just use LN(M1).

Code:
::
  CK2&Dispatch
  # 49
  ::
    {
      Z1_
    }
    Z0_
    NULL{}
    DUPDUP
    BINT1
    BINT1
    BINT0
    FPTR2 ^RCLVX
    '
    NULLLAM
    BINT11
    NDUPN
    DOBIND
    BINT105
    BINT110
    ClrUserFlag
    ClrUserFlag
    BINT110
    BINT110
    SysITE
    SetUserFlag
    SetSysFlag
    BINT105
    BINT105
    SysITE
    ::
      SetUserFlag
      11GETLAM
      FPTR2 ^CKNUMARRY
      11PUTLAM
    ;
    DROP
    11GETLAM
    DUP
    MDIMSDROP
    DROP
    4PUTLAM
    FPTR2 ^MATEGVL
    FPTR2 ^MATRIX2LIST
    DUPLENCOMP
    4GETLAM
    #<
    IT
    ::
      BINT105
      SetSysFlag
      DROP
      11GETLAM
      FPTR2 ^CKNUMARRY
      DUP
      11PUTLAM
      FPTR2 ^MATEGVL
      FPTR2 ^MATRIX2LIST
    ;
    7PUTLAM
    BEGIN
    4GETLAM
    #0<>
    WHILE
    ::
      6GETLAM
      7GETLAM
      CARCOMP
      >TCOMP
      DUP
      6PUTLAM
      5GETLAM
      BINT1
      >TCOMP
      5PUTLAM
      3GETLAM
      NTHCOMPDROP
      7GETLAM
      CDRCOMP
      DUP
      7PUTLAM
      SWAP
      EQUALPOSCOMP
      2PUTLAM
      BEGIN
      2GETLAM
      #0<>
      WHILE
      ::
        7GETLAM
        INNERDUP
        #2+
        2GETLAM
        #-
        ROLLDROP
        #1-{}N
        DUP
        7PUTLAM
        5GETLAM
        DUP
        3GETLAM
        NTHCOMPDROP
        #1+
        3GETLAM
        ROT
        PUTLIST
        5PUTLAM
        6GETLAM
        3GETLAM
        NTHCOMPDROP
        EQUALPOSCOMP
        2PUTLAM
      ;
      REPEAT
      3GETLAM
      #1+
      3PUTLAM
      7GETLAM
      LENCOMP
      4PUTLAM
    ;
    REPEAT
    3GETLAM
    #1-
    3PUTLAM
    11GETLAM
    MDIMSDROP
    DROP
    4PUTLAM
    9GETLAM
    4GETLAM
    ONE_DO
    1GETLAM
    INDEX@
    FPTR2 ^#>Z
    PTR 2EF1F
    >TCOMP
    LOOP
    10GETLAM
    >TCOMP
    DUP
    9PUTLAM
    1GETLAM
    6GETLAM
    CARCOMP
    PTR 2EF27
    FPTR2 ^FLAGSYMBEXEC
    FPTR2 ^LIST2MATRIX
    5GETLAM
    CARCOMP
    BINT1
    #>
    IT
    ::
      5GETLAM
      CARCOMP
      ONE_DO
      9GETLAM
      INDEX@
      ZERO_DO
      1GETLAM
      FPTR2 ^DERIVext
      LOOP
      1GETLAM
      6GETLAM
      CARCOMP
      PTR 2EF27
      xSUBST
      xAXL
      LOOP
    ;
    3GETLAM
    BINT1
    #>
    IT
    ::
      3GETLAM
      #1+
      BINT2
      DO
      9GETLAM
      1GETLAM
      6GETLAM
      INDEX@
      NTHCOMPDROP
      PTR 2EF27
      FPTR2 ^FLAGSYMBEXEC
      FPTR2 ^LIST2MATRIX
      5GETLAM
      INDEX@
      NTHCOMPDROP
      BINT1
      #>
      IT
      ::
        5GETLAM
        INDEX@
        NTHCOMPDROP
        ONE_DO
        9GETLAM
        INDEX@
        ZERO_DO
        1GETLAM
        FPTR2 ^DERIVext
        LOOP
        1GETLAM
        6GETLAM
        JINDEX@
        NTHCOMPDROP
        PTR 2EF27
        xSUBST
        xAXL
        LOOP
      ;
      LOOP
    ;
    4GETLAM
    FPTR 3 68
    FPTR2 ^MATRREF
    4GETLAM
    #1+
    FPTR2 ^MATRIX-COL
    SWAPDROP
    FPTR2 ^MATRIX2LIST
    DUP
    8PUTLAM
    CARCOMP
    4GETLAM
    UNCOERCE
    PTR 2F28F
    FPTR2 ^SCL*MAT
    4GETLAM
    #1+
    BINT2
    DO
    8GETLAM
    INDEX@
    NTHCOMPDROP
    11GETLAM
    INDEX@
    #1-
    FPTR2 ^#>Z
    FPTR2 ^MAT^
    FPTR2 ^SCL*MAT
    FPTR2 ^MAT+
    LOOP
    ABND
    FPTR2 ^FLAGEXPAND
    BINT110
    BINT110
    UserITE
    ClrUserFlag
    ClrSysFlag
    BINT105
    BINT105
    UserITE
    ClrUserFlag
    ClrSysFlag
  ;
;



RE: HP 49G: Function of a Matrix - Gerald H - 05-15-2015 03:09 PM

Here a second attempt, shorter & faster.

Input is a square matrix & a programme,

eg

[Matrix]
<< LN >>

Code:
::
  CK2&Dispatch
  # 48
  ::
    2DUP
    BINT0
    FPTR2 ^2LAMBIND
    FPTR2 ^PUSHFLAGS_
    FPTR2 ^CLREXACT
    BINT103
    SetSysFlag
    FPTR2 ^MATEGV
    FPTR2 ^XEQARRY>
    INCOMPDROP
    COERCEDUP
    1PUTLAM
    ZERO_DO
    2GETLAM
    EVAL
    1GETLAM
    ROLL
    LOOP
    1GETLAM
    TYPEMATRIX_
    COMPN_
    1GETABND
    DUP
    TWO{}N
    FPTR 3 6B
    OVER
    FPTR2 ^MATINV
    FPTR2 ^MMMULT
    FPTR2 ^MMMULT
    FPTR2 ^POPFLAGS_
    3UNROLL
    xDROP2
  ;
;



RE: (49G) Function of a Matrix - Gilles59 - 10-31-2017 05:16 PM

Hi Gerald

I think you can do this with the integrated DIAGMAP function. With your example :

Code:
[[ -9 0 -1 -4 ]
[ 4 9 3 3 ]
[ 0 4 -9 9 ]
[ 5 -1 -6 -4 ]]

'LN(X)' DIAGMAP XNUM

takes ~11,6s in approx mode on a 50g (i dont test with the 49g)

Then :

Code:
 'EXP(X)' DIAGMAP XNUM
and your retrieve exactly the real parts and x.xxxxE-12 for the imaginary parts.


RE: (49G) Function of a Matrix - Gerald H - 10-31-2017 10:04 PM

Well done, Gilles.

I've never used nor heard of DIAGMAP, how ignorant I am & how good a device the 50g is.

My programme above requires 6.7s to transform the example matrix, so remains of some use.


RE: (49G) Function of a Matrix - Gerald H - 10-31-2017 10:15 PM

If you try DIAGMAP with a your matrix * 1. &

<< LN >>

the time comes down to 6.7s.


RE: (49G) Function of a Matrix - Gilles59 - 10-31-2017 11:03 PM

(10-31-2017 10:15 PM)Gerald H Wrote:  If you try DIAGMAP with << LN >>
the time comes down to 6.7s.

Interesting... And no need of XNUM in this case!
But the precision is less. Curious. It seems that for example 'LN(X)' syntax use partially a symbolic calculation (so need an XNUM). << LN >> is full numerical


RE: (49G) Function of a Matrix - Gerald H - 11-01-2017 01:19 PM

I shall use the programme below in future as I only require numerical answers. The programme also has the advantage of returning LASTARG if required.

FPTR 7 134 is internal DIAGMAP
FPTR 7 1DC is internal PUSH
FPTR 7 1DD is internal POP

All three addresses are the same in OS 1.19-6 & 2.10-8.

Size: 54.

CkSum: # D173h

Code:
::
  CK2&Dispatch
  # 48
  ::
    2DUP
    FPTR 7 1DC
    FPTR2 ^CLREXACT
    BINT103
    SetSysFlag
    FPTR 7 134
    FPTR 7 1DD
    3UNROLL
    xDROP2
  ;
;



RE: (49G) Function of a Matrix - pier4r - 11-01-2017 06:06 PM

(10-31-2017 10:04 PM)Gerald H Wrote:  I've never used nor heard of DIAGMAP, how ignorant I am & how good a device the 50g is.

I agree with this a lot (and thanks Gilles for sharing!).

I estimate that I know the 50g way less than Gerald H, so I wonder how many years of problems I should solve before having a grasp of it.

PS: and thanks Gerald H to share your solutions, they provide nice inputs. I wish people would share their work often (I say this because in my experience many people don't share their work unless they know it is really good for many and that is a pity. See Gauss and his works, or Riemann. The same applies to small but potentially interesting works like little programs.)


RE: (49G) Function of a Matrix - Gerald H - 11-01-2017 06:27 PM

I agree with you wholeheartedly, pier4r.

Who knows what gems of programmes are never published.