hp41programs

Lauricella

Lauricella Functions for the HP-41


Overview
 

 1°)  Lauricella Function A

   a)  3 variables
   b)  n variables

 2°)  Lauricella Function B

   a)  3 variables
   b)  n variables

 3°)  Lauricella Function C

   a)  3 variables
   b)  n variables

 4°)  Lauricella Function D

   a)  3 variables: Focal Program
   b)  3 variables: M-Code Routine
   c)  Legendre Elliptic Integral of the 3rd kind
   d)  n variables
 

-These functions are defined by hypergeometric series of n variables.
-If n = 2, they correspond to the Appell hypergeometric functions F2 F3 F4 F1.
-If n = 1, all these 4 functions are equal to the ( Gauss ) hypergeometric series 2F1
 

1°)  Lauricella Function A
 

     a)  3 Variables
 

-This function is defined by

  FA( a ; b1 , b2 , b3 ; c1 , c2 , c3 ; x , y , z ) = SUMm,n,p=0 to infinity   (a)m+n+p (b1)m (b2)n (b3)p / ((c1)m (c2)n (c3)p)  xmynzp / (m! n! p!)

  where   (a)k = a(a+1)(a+2) ....... (a+k-1)   if  k > 0   and   (a)0 =  1    ( Pochhammer symbol )
 

Data Registers:               R00 , R08 to R10 , R14: temp               ( Registers R01 thru R07 are to be initialized before executing "LCFA3" )

                                      •  R01 = a       •  R02 = b1       •  R05 = c1      R11 = x
                                                            •  R03 = b2       •  R06 = c2      R12 = y
                                                            •  R04 = b3       •  R07 = c3      R13 = z
Flags: /
Subroutines: /
 
 

01  LBL "LCFA3"
02  STO 11         
03  RDN
04  STO 12
05  X<>Y
06  STO 13
07  CLX
08  STO 08
09  STO 09
10  STO 10
11  SIGN
12  STO 00
13  STO 14
14  SIGN
15  LBL 01
16  LASTX
17  RCL 01
18  RCL 08
19  +
20  RCL 09
21  +
22  RCL 10         
23  +
24  *
25  RCL 02
26  RCL 08
27  +
28  *
29  RCL 11
30  *
31  RCL 05
32  RCL 08
33  ST+ Y
34  ISG X
35  CLX
36  STO 08
37  *
38  /
39  +
40  X#Y?
41  GTO 01
42  CLX
43  STO 08 
44  X<> 00         
45  RCL 01
46  RCL 09
47  +
48  RCL 10
49  +
50  *
51  RCL 03
52  RCL 09
53  +
54  *
55  RCL 12
56  *
57  RCL 06
58  RCL 09
59  ST+ Y
60  ISG X
61  CLX
62  STO 09
63  *
64  /
65  STO 00         
66  +
67  X#Y?
68  GTO 01
69  CLX
70  STO 08
71  STO 09
72  X<> 14
73  RCL 01
74  RCL 10
75  +
76  *
77  RCL 04
78  RCL 10
79  +
80  *
81  RCL 13
82  *
83  RCL 07         
84  RCL 10 
85  ST+ Y
86  ISG X
87  CLX
88  STO 10
89  *
90  /
91  STO 00
92  STO 14
93  +
94  X#Y?
95  GTO 01
96  END

 
  ( 117 bytes / SIZE 015 )
 
 

      STACK        INPUTS      OUTPUTS
           Z             z       FA(x,y,z)
           Y             y       FA(x,y,z)
           X             x       FA(x,y,z)

 
Example:    With   a = 0.7 ,  b1 = 0.1  b2 = 0.2  b3 = 0.3 ,  c1 = 2.4  c2 = 2.5  c3 = 2.6       ( Store these 7 numbers into R01 thru R07 )

   0.13   ENTER^
   0.12   ENTER^
   0.11   XEQ "LCFA3"  >>>>>   FA( 0.11 , 0.12 , 0.13 ) = 1.021571587                ---Execution time = 104s---

   0.15   ENTER^
   0.25   ENTER^
   0.35      R/S        >>>>>   FA( 0.35 , 0.25 , 0.15 ) = 1.040641882                ---Execution time = 4mn28s---
 

Note:

-The series is convergent if  | x | + | y | + | z |  <  1
 

     b)  n Variables
 

-This function may be generalized to n variables as follows:
 

  FA( a ; b1 , ..... , bn ; c1 , ..... , cn ; x1 , ..... , xn ) = SUMi1,...,in=0 to infinity   (a)i1+...+in (b1)i1 ...... (bn)in / ((c1)i1 ... (cn)in)  x1i1 ... xnin / ( i1! .... in !)
 
 

Data Registers:           •  R00 = n               ( Registers R00 thru R3n+1 are to be initialized before executing "LCFA" )

                                      •  R01 = x1     •  Rn+1 = a       •  Rn+2 = b1       •  R2n+2 = c1               R3n+2 to R5n+1: temp
                                         .............                                ............................................

                                      •  Rnn = xn                              •  R2n+1 = bn     •  R3n+1 = cn

Flag:   F10  is cleared at the end
Subroutines: /
 

 2 M-code routines are used  X+1   to add 1 to register X
                                    and   X/E3  to divide X by 1000

  X+1  may be replaced by  ISG X  CLX
  X/E3 may be replaced by   PI   INT   10^X   /    ( we cannot use  E3  /   because  synthetic registers  P & Q are used )

-Of course, synthetic registers  M  N  O  P  Q  may be replaced by standard registers - for instance R95  R96  R97  R98  R99 - provided n is not too large
-The CX function  CLRGX  is also used.

-Line 113 is a three-byte GTO 12
 
 

  01  LBL "LCFA"
  02  RCL 00        
  03  STO N
  04  5
  05  *
  06  X+1 
  07  STO O
  08  .1
  09  %
  10  +
  11  RCL 00
  12  -
  13  CLRGX
  14  INT
  15  STO P
  16  RCL 00
  17  -
  18  STO Q
  19  RCL P
  20  X/E3
  21  ST+ O
  22  +
  23  X+1
  24  SIGN
  25  STO M
  26  LBL 00
  27  STO IND L
  28  ISG L
  29  GTO 00
  30  LBL 12
  31  SF 10
  32  LBL 01
  33  RCL 00        
  34  X+1
  35  RCL O
  36  RCL IND Y
  37  LBL 02
  38  RCL IND Y
  39  +
  40  DSE Y
  41  GTO 02
  42  RCL Q
  43  RCL 00
  44  -
  45  X<>Y
  46  RCL IND Y
  47  RCL IND O
  48  +
  49  *
  50  RCL IND Q
  51  RCL IND O
  52  ST+ Y
  53  X+1
  54  STO IND O
  55  *
  56  /
  57  RCL IND N
  58  *
  59  ST* IND P
  60  RCL IND P
  61  RCL M
  62  +
  63  STO M
  64  LASTX
  65  FC? 10
  66  GTO 04        
  67  X#Y?
  68  GTO 01
  69  LBL 03
  70  RCL 00
  71  ENTER^
  72  ST+ X
  73  ST+ X
  74  +
  75  X+1
  76  X/E3
  77  RCL O
  78  INT
  79  +
  80  CLRGX 
  81  DSE N
  82  X=0?
  83  GTO 06
  84  DSE O
  85  DSE P
  86  DSE Q
  87  CF 10
  88  GTO 01
  89  LBL 04
  90  RCL 00        
  91  ST+ X
  92  ST+ X
  93  X+1
  94  X/E3
  95  RCL P
  96  +
  97  RCL IND X
  98  LBL 05
  99  STO IND Y
100  ISG Y
101  GTO 05
102  R^
103  R^
104  X=Y?
105  GTO 03
106  RCL 00        
107  ENTER^
108  X<> N
109  -
110  ST+ O
111  ST+ P
112  ST+ Q
113  GTO 12
114  LBL 06
115  X<> M
116  CLA
117  END
 
 
    ( 194 bytes / SIZE 5n+2 )
 
 
      STACK        INPUTS     OUTPUTS
           X             /   FA(x1 , ... , xn)

 
Example:   the same one:    n = 3  STO 00

     x1 = 0.11   STO 01        a = 0.7   STO 04      b1 = 0.1   STO 05        c1 = 2.4    STO 08
     x2 = 0.12   STO 02                                        b2 = 0.2   STO 06        c2 = 2.5    STO 09
     x3 = 0.13   STO 03                                        b3 = 0.3   STO 07        c3 = 2.6    STO 10

  XEQ "LCFA"  >>>>   FA( 0.11 , 0.12 , 0.13 ) = 1.021571587                ---Execution time = 3m53s---

Note:

-The series is convergent if   | x1 | + ...... + | xn |  <  1
 

2°)  Lauricella Function B
 

     a)  3 Variables
 

-This function is defined by

  FB( a1 , a2 , a3 ; b1 , b2 , b3 ; c ; x , y , z ) = SUMm,n,p=0 to infinity    (a1)m (a2)n (a3)p (b1)m (b2)n (b3)p / (c)m+n+p   xmynzp / (m! n! p!)
 

Data Registers:               R00 , R08 to R10 , R14: temp               ( Registers R01 thru R07 are to be initialized before executing "LCFB3" )

                                      •  R01 = a1       •  R04 = b1       •  R07 = c         R11 = x
                                      •  R02 = a2       •  R05 = b2                                R12 = y
                                      •  R03 = a3       •  R06 = b3                                R13 = z
Flags: /
Subroutines: /
 
 

01  LBL "LCFB3"
02  STO 11         
03  RDN
04  STO 12
05  X<>Y
06  STO 13
07  CLX
08  STO 08
09  STO 09
10  STO 10
11  SIGN
12  STO 00
13  STO 14
14  SIGN
15  LBL 01
16  LASTX
17  RCL 01
18  RCL 08
19  +
20  *
21  RCL 04
22  RCL 08         
23  +
24  *
25  RCL 11
26  *
27  RCL 07
28  RCL 09
29  +
30  RCL 10 
31  +
32  RCL 08
33  ST+ Y
34  ISG X
35  CLX
36  STO 08
37  *
38  /
39  +
40  X#Y?
41  GTO 01
42  CLX
43  STO 08         
44  X<> 00
45  RCL 02
46  RCL 09
47  +
48  *
49  RCL 05
50  RCL 09 
51  +
52  *
53  RCL 12
54  *
55  RCL 07
56  RCL 10
57  +
58  RCL 09
59  ST+ Y
60  ISG X
61  CLX
62  STO 09         
63  *
64  /
65  STO 00
66  +
67  X#Y?
68  GTO 01
69  CLX
70  STO 08
71  STO 09 
72  X<> 14
73  RCL 03
74  RCL 10
75  +
76  *
77  RCL 06
78  RCL 10
79  +
80  *
81  RCL 13         
82  *
83  RCL 07
84  RCL 10
85  ST+ Y
86  ISG X
87  CLX
88  STO 10 
89  *
90  /
91  STO 00
92  STO 14
93  +
94  X#Y?
95  GTO 01
96  END

 
   ( 117 bytes / SIZE 015 )
 
 

      STACK        INPUTS      OUTPUTS
           Z             z       FB(x,y,z)
           Y             y       FB(x,y,z)
           X             x       FB(x,y,z)

 
Example:    With   a1 = 0.1  a2 = 0.2  a3 = 0.3 ,  b1 = 0.4  b2 = 0.5  b3 = 0.6 ,  c = 2.7      ( Store these 7 numbers into R01 thru R07 )

   0.43   ENTER^
   0.42   ENTER^
   0.41   XEQ "LCFB3"  >>>>>   FB( 0.41 , 0.42 , 0.43 ) = 1.057748024                ---Execution time = 2m38s---

   0.73   ENTER^
   0.72   ENTER^
   0.71      R/S        >>>>>   FB( 0.71 , 0.72 , 0.73 ) = 1.112251899                ---Execution time = 6mn35s---
 

Note:

-The series is convergent if   Max { | x | , | y | , | z |  }  <  1
 

     b)  n Variables
 

FB( a1 , ..... , an ; b1 , ..... , bn ; c ; x1 , ..... , xn ) = SUMi1,...,in=0 to infinity  (a1)i1 ...... (an)in (b1)i1 ... (bn)in / (c)i1+...+in     x1i1 ... xnin / ( i1! .... in!)
 
 

Data Registers:           •  R00 = n               ( Registers R00 thru R3n+1 are to be initialized before executing "LCFB" )

                                      •  R01 = x1     •  Rn+1 = a1       •  R2n+1 = b1         •  R3n+1 = c             R3n+2 to R5n+1: temp
                                         ..................................................................

                                      •  Rnn = xn      •  R2n = an        •  R3n = bn

Flag:   F10  is cleared at the end
Subroutines: /
 

 2 M-code routines are used  X+1   to add 1 to register X
                                    and   X/E3  to divide X by 1000

  X+1  may be replaced by  ISG X  CLX
  X/E3 may be replaced by   PI   INT   10^X   /

-Line 35 may be replaced by  ENTER^   ST+ X   +
-Line 115 is a three-byte GTO 12
 
 

  01  LBL "LCFB"
  02  RCL 00        
  03  STO N
  04  5
  05  *
  06  X+1 
  07  STO O
  08  .1
  09  %
  10  +
  11  RCL 00
  12  -
  13  CLRGX
  14  INT
  15  STO P
  16  RCL 00
  17  -
  18  STO Q
  19  RCL P
  20  X/E3
  21  ST+ O
  22  +
  23  X+1
  24  SIGN
  25  STO M
  26  ST- Q
  27  LBL 00
  28  STO IND L
  29  ISG L
  30  GTO 00        
  31  LBL 12
  32  SF 10
  33  LBL 01
  34  RCL 00
  35  3X 
  36  X+1
  37  RCL O
  38  RCL IND Y
  39  LBL 02
  40  RCL IND Y
  41  +
  42  DSE Y
  43  GTO 02
  44  RCL Q
  45  RCL 00
  46  -
  47  RCL IND X
  48  RCL IND O
  49  +
  50  R^
  51  /
  52  RCL IND Q
  53  RCL IND O
  54  ST+ Y
  55  X+1
  56  STO IND O
  57  /
  58  *
  59  RCL IND N
  60  *
  61  ST* IND P
  62  RCL IND P
  63  RCL M        
  64  +
  65  STO M
  66  LASTX
  67  FC? 10
  68  GTO 04
  69  X#Y?
  70  GTO 01
  71  LBL 03
  72  RCL 00
  73  ENTER^
  74  ST+ X
  75  ST+ X
  76  +
  77  X+1
  78  X/E3
  79  RCL O        
  80  INT
  81  +
  82  CLRGX 
  83  DSE N
  84  X=0?
  85  GTO 06
  86  DSE O
  87  DSE P
  88  DSE Q
  89  CF 10
  90  GTO 01
  91  LBL 04
  92  RCL 00
  93  ST+ X
  94  ST+ X
  95  X+1
  96  X/E3
  97  RCL P
  98  +
  99  RCL IND X
100  LBL 05
101  STO IND Y
102  ISG Y
103  GTO 05        
104  R^
105  R^
106  X=Y?
107  GTO 03
108  RCL 00
109  ENTER^
110  X<> N
111  -
112  ST+ O
113  ST+ P
114  ST+ Q
115  GTO 12
116  LBL 06
117  X<> M
118  CLA
119  END

 
   ( 198 bytes / SIZE 5n+2 )
 
 

      STACK        INPUTS      OUTPUTS
           X             /    FB(x1 , ... , xn)

 
Example:                 n = 3  STO 00

     x1 = 0.41   STO 01           a1 = 0.1   STO 04        b1 = 0.4    STO 07       c = 2.7   STO 10
     x2 = 0.42   STO 02           a2 = 0.2   STO 05        b2 = 0.5    STO 08
     x3 = 0.43   STO 03           a3 = 0.3   STO 06        b3 = 0.6    STO 09

  XEQ "LCFB"  >>>>   FB( 0.41 , 0.42 , 0.43 ) = 1.057748024                ---Execution time = 6m14s---

Note:

-The series is convergent if   Max { | x1 | ,......, | xn | }  <  1
 

3°)  Lauricella Function C
 

     a)  3 Variables
 

-This function is defined by

  FC( a ; b ; c1 , c2 , c3 ;  x , y , z ) = SUMm,n,p=0 to infinity    (a)m+n+p (b)m+n+p  / ( (c1)m(c2)n(c3)p )  xmynzp / (m! n! p!)
 

Data Registers:               R00 , R06 to R10: temp               ( Registers R01 thru R05 are to be initialized before executing "LCFC3" )

                                      •  R01 = a        •  R03 = b1                    R11 = x
                                      •  R02 = b        •  R04 = b2                    R12 = y
                                                              •  R05 = b3                    R13 = z
Flags: /
Subroutines: /
 
 

01  LBL "LCFC3"
02  STO 11         
03  RDN
04  STO 12
05  X<>Y
06  STO 13
07  CLX
08  STO 08
09  STO 09
10  STO 10
11  SIGN
12  STO 00
13  STO 06
14  SIGN
15  LBL 01
16  LASTX
17  RCL 08
18  RCL 09
19  +
20  RCL 10
21  +
22  STO 07         
23  RCL 01
24  +
25  *
26  RCL 02
27  RCL 07 
28  +
29  *
30  RCL 11
31  *
32  RCL 03
33  RCL 08
34  ST+ Y
35  ISG X
36  CLX
37  STO 08
38  *
39  /
40  +
41  X#Y?
42  GTO 01
43  CLX
44  STO 08         
45  X<> 00
46  RCL 09
47  RCL 10 
48  +
49  STO 07
50  RCL 01
51  +
52  *
53  RCL 02
54  RCL 07
55  +
56  *
57  RCL 12
58  *
59  RCL 04
60  RCL 09
61  ST+ Y
62  ISG X
63  CLX
64  STO 09         
65  *
66  /
67  STO 00 
68  +
69  X#Y?
70  GTO 01
71  CLX
72  STO 08
73  STO 09
74  X<> 06
75  RCL 01
76  RCL 10
77  +
78  *
79  RCL 02
80  RCL 10
81  +
82  *
83  RCL 13         
84  *
85  RCL 05
86  RCL 10 
87  ST+ Y
88  ISG X
89  CLX
90  STO 10
91  *
92  /
93  STO 00
94  STO 06
95  +
96  X#Y?
97  GTO 01
98  END

 
   ( 119 bytes / SIZE 014 )
 
 

      STACK        INPUTS      OUTPUTS
           Z             z       FC(x,y,z)
           Y             y       FC(x,y,z)
           X             x       FC(x,y,z)

 
Example:    With   a = 0.1 , b = 0.2 ,  c1 = 2.7   c2 = 2.8   c3 = 2.9         ( Store these 5 numbers into R01 thru R05 )

   0.12   ENTER^
   0.11   ENTER^
   0.10   XEQ "LCFC3"  >>>>>   FC( 0.10 , 0.11 , 0.12 ) = 1.002554253                ---Execution time = 2m34s---

Notes:

-The series is convergent if     | x |1/2 + | y |1/2 + | z |1/2   <  1
-This program may return a meaningless result if the above condition is not satisfied.
 

     b)  n Variables
 

FC( a ; b ; c1 , ..... , cn ; x1 , ..... , xn ) = SUMi1,...,in=0 to infinity  (a)i1+...+in (b)i1+...+in) / ( (c1)i1 ...... (cn)in)    x1i1 ... xnin / ( i1! .... in!)
 
 

Data Registers:           •  R00 = n               ( Registers R00 thru R3n+1 are to be initialized before executing "LCFC" )

                                      •  R01 = x1     •  Rn+1 = a       •  Rn+2 = b         •  Rn+3 = c1             R2n+3 to R4n+2: temp
                                         .............                                                            .................

                                      •  Rnn = xn                                                          •  R2n+2 = cn

Flag:   F10  is cleared at the end
Subroutines: /
 

 2 M-code routines are used  X+1   to add 1 to register X
                                    and   X/E3  to divide X by 1000

  X+1  may be replaced by  ISG X  CLX
  X/E3 may be replaced by  PI   INT   10^X   /

-Line 118 is a three-byte GTO 12
 
 

  01  LBL "LCFC"
  02  RCL 00        
  03  STO N
  04  4
  05  *
  06  2
  07  +
  08  STO O
  09  .1
  10  %
  11  +
  12  RCL 00
  13  -
  14  CLRGX
  15  INT
  16  STO P
  17  RCL 00
  18  -
  19  STO Q
  20  RCL P
  21  X/E3
  22  ST+ O
  23  +
  24  X+1
  25  SIGN
  26  STO M
  27  LBL 00
  28  STO IND L
  29  ISG L
  30  GTO 00
  31  LBL 12
  32  SF 10
  33  LBL 01        
  34  RCL O
  35  ENTER^
  36  CLX
  37  LBL 02
  38  RCL IND Y
  39  +
  40  DSE Y
  41  GTO 02
  42  RCL 00
  43  SIGN
  44  ST+ L
  45  RCL IND L
  46  X<>Y
  47  ST+ L
  48  RCL IND L
  49  R^
  50  ST+ T
  51  +
  52  R^
  53  *
  54  RCL IND Q
  55  RCL IND O
  56  ST+ Y
  57  X+1
  58  STO IND O
  59  *
  60  /
  61  RCL IND N
  62  *
  63  ST* IND P
  64  RCL IND P
  65  RCL M
  66  +
  67  STO M
  68  LASTX
  69  FC? 10
  70  GTO 04
  71  X#Y?
  72  GTO 01        
  73  LBL 03
  74  RCL 00
  75  ENTER^
  76  ST+ Y
  77  X+1
  78  ST+ X
  79  +
  80  X/E3
  81  RCL O
  82  INT
  83  +
  84  CLRGX 
  85  DSE N
  86  X=0?
  87  GTO 06
  88  DSE O
  89  DSE P
  90  DSE Q
  91  CF 10
  92  GTO 01        
  93  LBL 04
  94  RCL 00
  95  ENTER^
  96  X+1
  97  ST+ X
  98  +
  99  X/E3
100  RCL P
101  +
102  RCL IND X
103  LBL 05
104  STO IND Y
105  ISG Y
106  GTO 05        
107  R^
108  R^
109  X=Y?
110  GTO 03
111  RCL 00
112  ENTER^
113  X<> N
114  -
115  ST+ O
116  ST+ P
117  ST+ Q
118  GTO 12
119  LBL 06
120  X<> M
121  CLA
122  END

 
   ( 197 bytes / SIZE 4n+3 )
 
 

      STACK        INPUTS      OUTPUTS
           X             /    FC(x1 , ... , xn)

 
Example:               n = 3  STO 00

     x1 = 0.10   STO 01           a = 0.1   STO 04        b = 0.2    STO 05       c1 = 2.7   STO 06
     x2 = 0.11   STO 02                                                                                c2 = 2.8    STO 07
     x3 = 0.12   STO 03                                                                                c3 = 2.9    STO 08

  XEQ "LCFC"  >>>>   FC( 0.10 , 0.11 , 0.12 ) = 1.002554252                ---Execution time = 5m40s---

Notes:

-The series is convergent if     | x1 |1/2 + ........  + | xn |1/2   <  1
-This program may return a meaningless result if the above condition is not satisfied.
 

4°)  Lauricella Function D
 

     a)  3 Variables: Focal Program
 

-This function is defined by

  FD( a ; b1 , b2 , b3 ; c ;  x , y , z ) = SUMm,n,p=0 to infinity    (a)m+n+p (b1)m (b2)n (b3)p  /  (c)m+n+p   xmynzp / (m! n! p!)
 

Data Registers:               R00 , R06 to R10: temp                    ( Registers R01 thru R05 are to be initialized before executing "LCFD3" )

                                      •  R01 = a        •  R02 = b1       •  R05 = c            R11 = x
                                                              •  R03 = b2                                   R12 = y
                                                              •  R04 = b3                                   R13 = z
Flags: /
Subroutines: /
 
 

01  LBL "LCFD3"
02  STO 11         
03  RDN
04  STO 12
05  X<>Y
06  STO 13
07  CLX
08  STO 08
09  STO 09
10  STO 10
11  SIGN
12  STO 00
13  STO 06
14  SIGN
15  LBL 01
16  LASTX
17  RCL 08
18  RCL 09
19  +
20  RCL 10
21  +
22  STO 07         
23  RCL 01
24  +
25  *
26  RCL 05
27  RCL 07
28  +
29  /
30  RCL 11
31  *
32  RCL 02
33  RCL 08
34  ST+ Y
35  ISG X
36  CLX
37  STO 08 
38  /
39  *
40  +
41  X#Y?
42  GTO 01
43  CLX
44  STO 08         
45  X<> 00
46  RCL 09
47  RCL 10
48  +
49  STO 07
50  RCL 01
51  +
52  *
53  RCL 05
54  RCL 07
55  +
56  /
57  RCL 12 
58  *
59  RCL 03
60  RCL 09
61  ST+ Y
62  ISG X
63  CLX
64  STO 09         
65  /
66  *
67  STO 00
68  +
69  X#Y?
70  GTO 01
71  CLX
72  STO 08
73  STO 09
74  X<> 06
75  RCL 01
76  RCL 10 
77  +
78  *
79  RCL 05
80  RCL 10
81  +
82  /
83  RCL 13         
84  *
85  RCL 04
86  RCL 10
87  ST+ Y
88  ISG X
89  CLX
90  STO 10
91  /
92  *
93  STO 00
94  STO 06 
95  +
96  X#Y?
97  GTO 01
98  END

 
  ( 119 bytes / SIZE 014 )
 
 

      STACK        INPUTS      OUTPUTS
           Z             z       FD(x,y,z)
           Y             y       FD(x,y,z)
           X             x       FD(x,y,z)

 
Example:      With   a = 0.1 ,  b1 = 0.2   b2 = 0.3   b3 = 0.4 ,  c = 2.7             ( Store these 5 numbers into R01 thru R05 )

   0.41   ENTER^
   0.31   ENTER^
   0.21   XEQ "LCFD3"  >>>>>   FD( 0.21 , 0.31 , 0.41 ) = 1.012345372                ---Execution time = 2m53s---

    0.7    ENTER^
    0.6    ENTER^
    0.5       R/S        >>>>>   FD( 0.5 , 0.6 , 0.7 ) = 1.025992608                ---Execution time = 11m10s---
 

Notes:

-The series is convergent if    Max { | x | , | y | , | z |  }  <  1

-There are also 10 other Saran-Lauricella functions of 3 variables which may be programmed in a similar way.
-All these triple series are special cases of a generalized hypergeometric series of 3 variables.
-See for instance "Srivastava's Functions of Three Variables for the HP-41"
 

     b)  3 Variables: M-Code Routine
 

-We'll get faster - or at least less slow - results with M-Code.
-Registers R01 thru R05 are replaced by R00 thru R04.
-Their existence is tested but there is no check for alpha data.

-The synthetic registers M & N are denoted m & n in the simplified mnemonics below.
-The alpha "register" is cleared.
 

Data Registers:                                  ( Registers R00 thru R04 are to be initialized before executing "LCFD3" )

                                      •  R00 = a        •  R01 = b1       •  R04 = c
                                                              •  R02 = b2
                                                              •  R03 = b3
Flags: /
Subroutines: /
 

0B3  "3"
004  "D"
006  "F"
003  "C"
00C  "L"
3C8  CLRKEY        first executable word
378   C=c
03C  RCR 3
106  A=C S&X
130  LDI S&X
1FC 508d               correct value for an HP-41CV/CX or C with a quad memory module
306  ?A<C S&X
381  ?NCGO
00A  02E0              ( NONEXISTENT )
2A0  SETDEC
345   ?NCXQ
040   CLA
0F8   C=X
128   L=C
078   C=Z
028   T=C
04E   C
35C   =
050   1
0E8   X=C
068   Z=C
228   P=C
268   Q=C
178   C=M                LOOP  at the address  E48D  in my ROM - Change the 3 ?CGO E48D  written in red below according to your own ROM
10E  A=C ALL
1B8  C=n
01D  C=
060   A+C
1F8  C=O
025  C=
060  AB+C
070  N=C ALL
378  C=c
03C  RCR 3
270  RAMSLCT
038  READATA
025  C=
060  AB+C
046  C=0 S&X
270  RAMSLCT
078  C=Z
13D  C=
060   AB*C
068   Z=C
378   C=c
03C  RCR 3
260  SETHEX
226  C=C+1 S&X
226  C=C+1 S&X
226  C=C+1 S&X
226  C=C+1 S&X
270  RAMSLCT
038  READATA
10E  A=C ALL
0F0  C<>N ALL
2A0  SETDEC
01D  C=
060   A+C
10E   A=C ALL
046   C=0 S&X
270   RAMSLCT
078   C=Z
0AE  A<>C ALL
261   C=
060   A/C
138   C=L
13D  C=
060   AB*C
068   Z=C
178   C=m
10E   A=C ALL
378   C=c
03C  RCR 3
260   SETHEX
226   C=C+1 S&X
270  RAMSLCT
038  READATA
2A0  SETDEC
01D  C=
060   A+C
046   C=0 S&X
270   RAMSLCT
078   C=Z
13D  C=
060   AB*C
068   Z=C
178   C=m
00E   A
35C  =
162   1
01D  C=
060   A+C
168   m=C
10E   A=C ALL
078   C=Z
0AE  A<>C ALL
261   C=
060   A/C
068   Z=C
3CC  ?KEY                               These 2 words allow you to stop the routine
360   ?C RTN                             if it's too long: press any key to stop it.
0F8   C=X
025   C=
060   AB+C
10E   A=C ALL
0F8   C=X
0AE  A<>C ALL
0E8   X=C
36E   ?A#C ALL
235   ?CGO                  ?CGO
393    E48D                  LOOP
1B8   C=n
10E   A=C ALL
1F8   C=O
01D  C=
060   A+C
168   m=C
378   C=c
03C  RCR 3
270  RAMSLCT
038  READATA
025   C=
060   AB+C
046   C=0 S&X
270   RAMSLCT
238   C=P
13D  C=
060   AB*C
228   P=C
178   C=m
10E   A=C ALL
0B0   C=N ALL
01D  C=
060   A+C
10E   A=C ALL
238   C=P
0AE  A<>C ALL
261   C=
060   A/C
0B8  C=Y
13D  C=
060   AB*C
228   P=C
378   C=c
03C  RCR 3
260  SETHEX
226  C=C+1 S&X
226  C=C+1 S&X
270  RAMSLCT
038  READATA
10E  A=C ALL
04E  C=0 ALL
270   RAMSLCT
168   m=C
2A0  SETDEC
1B8  C=n
01D  C=
060   A+C
238   C=P
13D  C=
060   AB*C
228   P=C
1B8  C=n
00E   A
35C  =
162   1
01D  C=
060   A+C
1A8   n=C
10E   A=C ALL
238   C=P
0AE  A<>C ALL
261   C=
060   A/C
228   P=C
068   Z=C
0F8   C=X
025   C=
060   AB+C
10E   A=C ALL
0F8   C=X
0AE  A<>C ALL
0E8   X=C
36E   ?A#C ALL
235   ?CGO                  ?CGO
393    E48D                  LOOP
1F8   C=O
10E   A=C ALL
378   C=c
03C  RCR 3
270  RAMSLCT
038  READATA
01D  C=
060   A+C
046   C=0 S&X
270   RAMSLCT
278   C=Q
13D  C=
060   AB*C
268   Q=C
0B0  C=N ALL
10E   A=C ALL
1F8   C=O
01D  C=
060   A+C
10E   A=C ALL
278   C=Q
0AE  A<>C ALL
261   C=
060   A/C
04E   C=0 ALL
1A8  n=C
270   RAMSLCT
038   C=T
13D  C=
060   AB*C
268   Q=C
378   C=c
03C  RCR 3
260  SETHEX
226  C=C+1 S&X
226  C=C+1 S&X
226  C=C+1 S&X
270  RAMSLCT
038  READATA
10E  A=C ALL
046   C=0 S&X
270   RAMSLCT
2A0  SETDEC
1F8   C=O
01D  C=
060   A+C
278   C=Q
13D  C=
060   AB*C
268   Q=C
1F8   C=O
00E   A
35C  =
162   1
01D  C=
060   A+C
1E8   O=C
10E   A=C ALL
278   C=Q
0AE  A<>C ALL
261   C=
060   A/C
268   Q=C
228   C=P
068   Z=C
0F8   C=X
025   C=
060   AB+C
10E   A=C ALL
0F8   C=X
0AE  A<>C ALL
0E8   X=C
36E   ?A#C ALL
235   ?CGO                  ?CGO
393    E48D                  LOOP
046    C
270    =
038    T
068    Z=C
345    ?NCGO
042    CLA

( 271 words / SIZE 005 )
 
 

      STACK        INPUTS      OUTPUTS
           Z             z             z
           Y             y             y
           X             x       FD(x,y,z)
           L             /             x

 
Example:      With   a = 0.1 ,  b1 = 0.2   b2 = 0.3   b3 = 0.4 ,  c = 2.7                 ( Store these 5 numbers into R00 thru R04 )

   0.41   ENTER^
   0.31   ENTER^
   0.21   XEQ "LCFD3"  >>>>>   FD( 0.21 , 0.31 , 0.41 ) = 1.012345372                ---Execution time = 1m06s---

Note:

-The M-Code routine is almost 3 times as fast as the focal program
  but it remains relatively slow...
 

     c)  Legendre Elliptic Integral of the 3rd Kind
 

-Legendre elliptic integrals of the 3rd kind may be expressed in terms of Lauricella functions FD
-More explicitly:

      ¶ ( n ; phi | m ) = §0phi ( 1 + n sin2 t ) -1 ( 1 - m sin2 t ) -1/2 dt
                             = sin (phi) . FD ( 1/2 ; 1 , 1/2 , 1/2 ; 3/2 ;  - n sin2 (phi) , sin2 (phi) , m sin2 (phi) )
 
 

 01  LBL "LEI3"
 02  SIN
 03  STO 05
 04  X^2
 05  ST* Y
 06  ST* Z
 07  .5
 08  STO 00
 09  STO 02
 10  STO 03
 11  STO 04
 12  SIGN
 13  STO 01
 14  ST+ 04
 15  X<> T
 16  CHS
 17  LCFD3
 18  RCL 05
 19  *
 20  END

 
( 35 bytes / SIZE 006 )
 
 

      STACK       INPUTS    OUTPUTS
          Z            n           /
          Y            m           /
          X           phi   ¶( n ; phi | m )

 
Example:     in DEG mode:

       0.3  ENTER^
       0.7  ENTER^
       41   XEQ "LEI3"  >>>>  ¶ ( 0.3 ; 41° | 0.7 ) =  0.726768062                        ---Execution time = 3mn14s---

Notes:

-Obviously, this version remains very slow, even with the M-Code LCFD3.
-So, it cannot compete with Carlson elliptic integrals.
-It's just given here for completeness.

-Like in "Numerical Recipes", the sign convention for n is opposite that of Abramowitz & Stegun
-Otherwise, simply delete line 16 CHS
-The angle phi must be between  -90° and +90°
 

     d)  n Variables
 

 FD( a ; b1 , ..... , bn ; c ; x1 , ..... , xn ) = SUMi1,...,in=0 to infinity  (a)i1+....+in (b1)i1 ... (bn)in / (c)i1+...+in     x1i1 ... xnin / ( i1! .... in!)
 
 

Data Registers:           •  R00 = n               ( Registers R00 thru R2n+2 are to be initialized before executing "LCFD" )

                                      •  R01 = x1     •  Rn+1 = a        •  Rn+2 = b1         •  R2n+2 = c             R2n+3 to R4n+2: temp
                                         ...........                                      .................

                                      •  Rnn = xn                                •  R2n+1 = bn

Flag:   F10  is cleared at the end
Subroutines: /
 

 2 M-code routines are used  X+1   to add 1 to register X
                                    and   X/E3  to divide X by 1000

  X+1  may be replaced by  ISG X  CLX
  X/E3 -------------------   PI   INT   10^X   /

-Line 119 is a three-byte GTO 12
 
 

  01  LBL "LCFD"
  02  RCL 00        
  03  STO N
  04  4
  05  *
  06  2
  07  +
  08  STO O
  09  .1
  10  %
  11  +
  12  RCL 00
  13  -
  14  CLRGX
  15  INT
  16  STO P
  17  RCL 00
  18  -
  19  STO Q
  20  RCL P
  21  X/E3
  22  ST+ O
  23  +
  24  X+1
  25  SIGN
  26  STO M
  27  ST- Q
  28  LBL 00
  29  STO IND L
  30  ISG L
  31  GTO 00
  32  LBL 12
  33  SF 10
  34  LBL 01
  35  RCL O
  36  ENTER^
  37  CLX
  38  LBL 02
  39  RCL IND Y
  40  +
  41  DSE Y
  42  GTO 02
  43  RCL 00        
  44  SIGN
  45  ST+ L
  46  RCL IND L
  47  LASTX
  48  ST+ X
  49  RDN
  50  RCL IND T
  51  R^
  52  ST+ Z
  53  +
  54  /
  55  RCL IND Q
  56  RCL IND O
  57  ST+ Y
  58  X+1
  59  STO IND O
  60  /
  61  *
  62  RCL IND N
  63  *
  64  ST* IND P
  65  RCL IND P
  66  RCL M
  67  +
  68  STO M
  69  LASTX
  70  FC? 10
  71  GTO 04        
  72  X#Y?
  73  GTO 01
  74  LBL 03
  75  RCL 00
  76  ENTER^
  77  ST+ Y
  78  X+1
  79  ST+ X
  80  +
  81  X/E3
  82  RCL O
  83  INT
  84  +
  85  CLRGX
  86  DSE N
  87  X=0?
  88  GTO 06 
  89  DSE O
  90  DSE P
  91  DSE Q
  92  CF 10
  93  GTO 01
  94  LBL 04
  95  RCL 00        
  96  ENTER^
  97  X+1
  98  ST+ X
  99  +
100  X/E3
101  RCL P
102  +
103  RCL IND X
104  LBL 05
105  STO IND Y
106  ISG Y
107  GTO 05
108  R^
109  R^
110  X=Y?
111  GTO 03
112  RCL 00        
113  ENTER^
114  X<> N
115  -
116  ST+ O
117  ST+ P
118  ST+ Q
119  GTO 12
120  LBL 06
121  X<> M
122  CLA
123  END

 
   ( 199 bytes / SIZE 4n+3 )
 
 

      STACK        INPUTS      OUTPUTS
           X             /    FD(x1 , ... , xn)

 
Example:               n = 3  STO 00

     x1 = 0.21   STO 01           a = 0.1   STO 04        b1 = 0.2   STO 05       c = 2.7   STO 08
     x2 = 0.31   STO 02                                              b2 = 0.3    STO 06
     x3 = 0.41   STO 03                                              b3 = 0.4    STO 07

  XEQ "LCFD"  >>>>   FD( 0.21 , 0.31 , 0.41 ) = 1.012345372                ---Execution time = 5m47s---

Note:

-The series is convergent if   Max { | x1 | ,......, | xn | }  <  1
 

Reference:

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