hp41programs

Kampé de Fériet

Kampé de Fériet Function for the HP-41


Overview
 

 1°)  Focal Program
 2°)  M-Code Routine
 

-Kampé de Fériet function is the generalized hypergeometric function of 2 variables.
-Given 6 non-negative integers:  A , B , C , P , Q , S ,  it is defined by the double series:
 

           A B C   /    a1 ,........, aA ;  b1 ,........, bB ;  c1 ,........, cC ;              \
        FP Q S   |                                                                             x , y    |      =    SUMm=0,1,2,...  SUMn=0,1,2,..    Km,n  xm yn / ( m! n! )
                      \    p1 ,........, pP ;  q1 ,........, qQ ;  s1 ,........, sS ;               /
 

      where      Km,n  =  [ (a1)m+n ......... (aA)m+n  (b1)m ......... (bB)m  (c1)n ......... (cC)n ] /  [ (p1)m+n ......... (pP)m+n  (q1)m ......... (qQ)m  (s1)n ......... (sS)n ]
 

       and    (a)n = Pochhammer's symbol:     (a)0 = 1
                                                                  (a)n = a(a+1)(a+2) ...... (a+n-1)    if  n > 0
 
 

1°)  Focal Program
 
 

Data Registers:              R00 = x                     ( Registers R01 thru R06+A+B+C+P+Q+S  are to be initialized before executing "KdF" )

                                      •  R01 = A     •  R07 =  a1   .......................................   •  R06+A = aA
                                      •  R02 = B     •  R07+A =  b1 ....................................   •  R06+A+B = bB
                                      •  R03 = C     •  R07+A+B =  c1 ................................  •  R06+A+B+C = cC
                                      •  R04 = P      •  R07+A+B+C =  p1 ..........................  •  R06+A+B+C+P = pP
                                      •  R05 = Q     •  R07+A+B+C+P =  q1 ......................  •  R06+A+B+C+P+Q = qQ
                                      •  R06 = S      •  R07+A+B+C+P+Q =  s1 .................  •  R06+A+B+C+P+Q+S = sS
Flag:  F10
Subroutines: /
 

-Line 169 is a three-byte  GTO 04  and  ""  is a TEXT0 ( NOP ) instruction
 
 

  01  LBL "KdF"
  02  CLA
  03  STO O
  04  X<>Y
  05  STO P
  06  CLX
  07  SIGN
  08  STO Q
  09  SIGN
  10  GTO 04
  11  LBL 01
  12  ISG 00
  13  FS? 30
  14  RTN
  15  RCL M
  16  RCL N
  17  +
  18  RCL IND 00
  19  +
  20  FC? 10
  21  1/X
  22  *
  23  GTO 01
  24  LBL 02
  25  ISG 00
  26  FS? 30
  27  RTN
  28  RCL M        
  29  RCL IND 00
  30  +
  31  FC? 10
  32  1/X
  33  *
  34  GTO 02
  35  LBL 03
  36  ISG 00
  37  FS? 30
  38  RTN
  39  RCL N
  40  RCL IND 00
  41  +
  42  FC? 10
  43  1/X
  44  *
  45  GTO 03
  46  LBL 04
  47  LASTX
  48  PI
  49  INT
  50  ST+ X
  51  STO 00
  52  RCL 01
  53  +
  54  PI
  55  INT
  56  10^X
  57  /
  58  ST+ 00
  59  RDN
  60  SF 10
  61  XEQ 01
  62  RCL 02        
  63  PI
  64  INT
  65  10^X
  66  /
  67  DSE 00
  68  ""
  69  ST+ 00
  70  RDN
  71  XEQ 02
  72  RCL 03
  73  ST+ 00
  74  RCL 04
  75  +
  76  PI
  77  INT
  78  10^X
  79  /
  80  DSE 00
  81  ""
  82  ST+ 00
  83  RDN
  84  CF 10
  85  XEQ 01
  86  RCL 05
  87  PI
  88  INT
  89  10^X
  90  /
  91  DSE 00
  92  ""
  93  ST+ 00
  94  RDN
  95  XEQ 02
  96  RCL O
  97  *
  98  RCL M        
  99  ISG X
100  ""
101  STO M
102  /
103  +
104  X#Y?
105  GTO 04
106  PI
107  INT
108  ST+ X
109  STO 00
110  RCL 01
111  +
112  PI
113  INT
114  10^X
115  /
116  ST+ 00
117  CLX
118  STO M        
119  X<> Q
120  SF 10
121  XEQ 03
122  RCL 02
123  ST+ 00
124  RCL 03
125  +
126  PI
127  INT
128  10^X
129  /
130  DSE 00
131  ""
132  ST+ 00
133  RDN
134  XEQ 03
135  RCL 04
136  PI
137  INT
138  10^X
139  /
140  DSE 00
141  ""
142  ST+ 00
143  RDN
144  CF 10
145  XEQ 03
146  RCL 05
147  ST+ 00
148  RCL 06
149  +
150  PI
151  INT
152  10^X
153  /
154  DSE 00
155  ""
156  ST+ 00
157  RDN
158  XEQ 03
159  RCL P
160  *
161  RCL N        
162  ISG X
163  ""
164  STO N
165  /
166  STO Q
167  +
168  X#Y?
169  GTO 04
170  CLA
171  END

 
   ( 255 bytes / SIZE 007+A+B+C+P+Q+S )
 
 

      STACK        INPUTS      OUTPUTS
           Y             y          F(x,y)
           X             x          F(x,y)

 
Example:

           A = 2   B = 3   C = 2                  a1 = 1.2    a2 = 2.3  ,   b1 = 1.1   b2 = 0.7   b3 = 1.4   ,   c1 = 1.6   c2 = 1.5
           P = 3   Q = 2   S = 4            p1 = 2    p2 =  3   p3 =  4  ,  q1 = 5   q2 = 6  ,  s1 = 7    s2 = 8   s3 = 9   s4 =  10.1

-Calculate  F(5,7)

-Store  these 22 numbers into

               R01  R02  R03                                R07  R08 ,  R09  R10  R11 ,  R12  R13
               R04  R05  R06                         R14  R15  R16 , R17  R18 , R19  R20  R21  R22

   7   ENTER^
   5   XEQ "KdF"  >>>>  F(5,7) = 1.022466709         ---Execution time = 129s---
 

Note:

-If  A = B = C = P = Q = S = 0 ,  "KdF"  returns  exp(x+y)
 

2°)  M-Code Routine
 

-This routine is structured like the program listed above.
 

Data Registers:              R00 = unused                 ( Registers R01 thru R06+A+B+C+P+Q+S  are to be initialized before executing "KdF" )

                                      •  R01 = A     •  R07 =  a1   .......................................   •  R06+A = aA
                                      •  R02 = B     •  R07+A =  b1 ....................................   •  R06+A+B = bB
                                      •  R03 = C     •  R07+A+B =  c1 ................................  •  R06+A+B+C = cC
                                      •  R04 = P      •  R07+A+B+C =  p1 ..........................  •  R06+A+B+C+P = pP
                                      •  R05 = Q     •  R07+A+B+C+P =  q1 ......................  •  R06+A+B+C+P+Q = qQ
                                      •  R06 = S      •  R07+A+B+C+P+Q =  s1 .................  •  R06+A+B+C+P+Q+S = sS
Flags: /
Subroutines: /
 

-There are 9 absolute jumps ( one  ?C GO  and  8  ?NC XQ ) written in red & green.
-You'll have to change these instructions according to the location of the corresponding lines ( written in the same colors ) in your own ROM.
 

086   "F"                                         at the address FF06 in my ROM
044   "d"
00B  "K"
378   READ 13(c)
03C  RCR 3
070   N=C ALL
106   A=C S&X
130   LDI S&X
006   006
146   A=A+C S&X
130   LDI S&X
200   200h                          this value is for an HP-41CV , CX or C with a Quad memory module or 4 memory modules
306   ?A<C S&X               if you have an HP-41C without memory modules, replace 200 by 100
381   ?NC GO
00A   02E0                        if register R06 doesn't exist, the routine stops after displaying "NONEXISTENT"
248   SETF 9
06E   A<>B ALL
130   LDI S&X
002   002
158   M=C ALL
0B0   C=N ALL
226   C=C+1 S&X
070   N=C ALL
270   RAMSLCT
038   READDATA
38D  ?NCXQ
008   02E3                         the different numbers in registers R01 thru R06 are converted to hexadecimal numbers
106   A=C S&X                and the corresponding addresses are stored in CPU-register N and in synthetic register Q
0CE  C=B ALL
1BC  RCR 11
0C6  C=B S&X
206   C=C+A S&X
0EE   C<>B ALL
198   C=M ALL
266   C=C-1 S&X
383   JNC -16d=-10h
046   C=0 S&X
270   RAMSLCT
0CE  C=B ALL
24C  ?FSET 9
023   JNC +04
244   CLRF 9
268   WRIT 9(Q)
333   JNC -26d=-1Ah
228   WRIT 8(P)
106   A=C S&X
130   LDI S&X
200   200h                          this value is for an HP-41CV , CX or C with a Quad memory module or 4 memory modules
306   ?A<C S&X               if you have an HP-41C without any additional memory module, replace 200 by 100
381   ?NC GO
00A   02E0                        if register R06+A+B+C+P+Q+S doesn't exist, the routine stops after displaying "NONEXISTENT"
2A0   SETDEC
0F8   READ 3(X)
128   WRIT 4(L)
04E   C=0 ALL
1A8   WRIT 6(N)
35C  PT=12                      C=
050   LD@PT- 1                1
1E8   WRIT 7(O)
0E8   WRIT 3(X)
068   WRIT 1(Z)
04E   C=0 ALL                Loop           this word is at the address  FF43  in my ROM
168   WRIT 5(M)
278   READ 9(Q)                                If you don't have a "newest" HP-41,  insert the 2 words  3CC   ?KEY   360   ?C RTN    in this area
17C  RCR 6                                        You'll simply have to press a key to stop an infinite loop
104   CLRF 8                                      if the series is divergent or converges too slowly.
288   CETF 7                                      In this case however, the address of the subroutine will have to be changed...
2BD  ?NXQ                   ?ncxq
3FC   FFAF                     sub              But with a "newest" HP-41, simply press  ENTER^   ON  to stop the routine.
278   READ 9(Q)
03C  RCR 3
244   CLRF 9
284   CLRF 7
2BD  ?NXQ                   ?ncxq
3FC   FFAF                     sub
238   READ 8(P)
17C  RCR 6
108   SETF 8
288   SETF 7
2BD  ?NXQ                   ?ncxq
3FC   FFAF                     sub
238   READ 8(P)
03C  RCR 3
284   CLRF 7
2BD  ?NXQ                   ?ncxq
3FC   FFAF                     sub
2A0   SETDEC
078    READ 1(Z)
10E   A=C ALL
138   READ 4(L)
135   ?NCXQ                 C=
060   184D                    A*C
070   N=C ALL
178   READ 5(M)
00E   A=0 ALL                A
35C  PT=12                     =
162   A=A+1@PT           1
01D   ?NCXQ                C=
060    1807                    A+C
168    WRIT 5(M)
0B0   C=N ALL
10E   A=C ALL
178   READ 5(M)
261   ?NCXQ                C=
060   1898                     A/C
068   WRIT 1(Z)
10E   A=C ALL
0F8   READ 3(X)
01D   ?NCXQ                C=
060    1807                    A+C
10E    A=C ALL
0F8   READ 3(X)
0AE  A<>C ALL
0E8   WRIT 3(X)
36E   ?A#C ALL
267   JC-34h=-52d
1F8   READ 7(O)
068   WRIT 1(Z)
278   READ 9(Q)
17C  RCR 6
104   CLRF 8
248   CLFR 9
2BD  ?NXQ                   ?ncxq
3FC   FFAF                     sub
278   READ 9(Q)
2BD  ?NXQ                   ?ncxq
3FC   FFAF                     sub
108   SETF 8
238   READ 8(P)
17C  RCR 6
2BD  ?NXQ                   ?ncxq
3FC   FFAF                     sub
238   READ 8(P)
2BD  ?NXQ                   ?ncxq
3FC   FFAF                     sub
2A0   SETDEC
078    READ 1(Z)
10E   A=C ALL
0B8   READ 2(Y)
135   ?NCXQ                 C=
060   184D                    A*C
070   N=C ALL
1B8   READ 6(N)
00E   A=0 ALL                A
35C  PT=12                     =
162   A=A+1@PT           1
01D   ?NCXQ                C=
060    1807                    A+C
1A8    WRIT 6(N)
0B0   C=N ALL
10E   A=C ALL
1B8   READ 6(N)
261   ?NCXQ                C=
060   1898                     A/C
068   WRIT 1(Z)
1E8   WRIT 7(O)
10E   A=C ALL
0F8   READ 3(X)
01D   ?NCXQ                C=
060    1807                    A+C
10E    A=C ALL
0F8   READ 3(X)
0AE  A<>C ALL
0E8   WRIT 3(X)
36E   ?A#C ALL
10D   ?C GO                 ?c go
3FF    FF43                    loop
345    ?NCGO              ?ncgo
042    10D1                   CLA                the routine stops after clearing the alpha register
260    SETHEX              sub                  The first word of this subroutine is at the address  FFAF in my ROM
23A   C=C+1 M
070   N=C ALL
106   A=C S&X
03C   RCR 3
306   ?A<C S&X
360   ?C RTN
270   RAMSLCT
038   READDATA
10E   A=C ALL
046   C=0 S&X
270   RAMSLCT
2A0  SETDEC
178   READ 5(M)
28C  ?FSET 7
033   JNC+06
01D   ?NCXQ                C=
060    1807                    A+C
10E    A=C ALL
1B8   READ 6(N)
023   JNC+04
24C  ?FSET 9
013   JNC+02
1B8  READ 6(N)
01D   ?NCXQ                C=
060    1807                    A+C
10E    A=C ALL
078   READ 1(Z)
0AE  A<>C ALL
10C  ?FSET 8
135   ?NCXQ                 C=
060   184D                    A*C
10C  ?FSET 8
261   ?C XQ                  C=
061   1898                     A/C
068   WRIT 1(Z)
0B0   C=N ALL
2DB  JNC-25h=-37d                          ( at the address FFD4 in my ROM )
 
 

      STACK        INPUTS      OUTPUTS
           Y             y             y
           X             x          F(x,y)
           L             /             x

 
Example1:     With the same parameters,

   7   ENTER^
   5   XEQ "KdF"   returns    F(5,7) = 1.022466709         ---Execution time = 15s---
 

Example2:     Appell hypergeometric function F1 is a special case of Kampé de Fériet functions with  A = B = C = P = 1 and  Q = S = 0

-For instance, to compute  F1( 1 , 2 , 3 , 7 ; 0.5 , 0.7 )

-Store  1  into registers  R01  R02  R03  R04
            0  into registers  R05  R06
   and   1  2  3  7   into   R07  R08  R09  R10

-Then,

  0.7  ENTER^
  0.5  XEQ "KdF"   gives     F1( 1 , 2 , 3 , 7 , 0.5 , 0.7 )  = 1.856874937      ---Execution time = 2mn20s
 

Warning:     This program does not check for alpha data.

-Therefore, be sure that registers R01 thru R06+A+B+C+P+Q+S only contain numbers!
 

Notes:

-If  A = B = C = P = Q = S = 0 ,  "KdF"  returns  exp(x+y)
-Register Z contains the first neglected term.
-Register T is unchanged.
-The alpha "register" is cleared.
 
 

Reference:

[1]   http://mathworld.wolfram.com/KampedeFerietFunction.html