# hp41programs

Anionic Functions4

# Anionic Special Functions(IV) for the HP-41

"Fractional-Integro-Differentiation"

Overview

1°)  A few Elementary Functions
2°)  Kummer's Functions
3°)  Error Function
4°)  Sine & Hyperbolic Sine Integral
5°)  Cosine & Hyperbolic Cosine Integral
6°)  Exponential Integral Ei
7°)  Generalized Laguerre's Functions
8°)  Bessel & Modified Bessel Functions of the 1st Kind
9°)  Bessel & Modified Bessel Functions of the 2nd Kind - Non-Integer Order
10°) Spherical Bessel Functions of the 1st Kind
11°) Hermite Functions
12°) Fresnel Integrals
13°) Airy Functions

-The fractional-integro-differentiation of a complex function is quite well defined.
-For quaternions, octonions and so on, this is more doubtful:
-Even the 1st derivative of  a2 is not  2.a  because the multiplication is not commutative in general.
-So, these functions may be regarded as other special functions not always related to differentiation of special functions... except for complex arguments !

-I've used M-Code routines which may be replaced by focal programs with small modifications in the listings below

1°)  A few Elementary Functions

Formulae:

•  Real Power                          Dµ ax  = [ x ! / Gam(x-µ+1) ] ax-µ           assuming  x # -1 , -2 , -3 , .....

•  Hyperbolic Sine                   Dµ Sinh a = 2µ-1 sqrt(PI) a1-µ1F~2 ( 1 ; (2-µ)/2 , (3-µ)/2 ; a2/4 )

•  Hyperbolic Cosine               Dµ Cosh a = (2/a)µ sqrt(PI)  1F~2 ( 1 ; (1-µ)/2 , (2-µ)/2 ; a2/4 )

•  Sine                     Dµ Sin a = 2µ-1 sqrt(PI) a1-µ1F~2 ( 1 ; (2-µ)/2 , (3-µ)/2 ; -a2/4 )

•  Cosine                 Dµ Cos x = (2/a)µ sqrt(PI)  1F~2 ( 1 ; (1-µ)/2 , (2-µ)/2 ; -a2/4 )

•  Logarithm            Dµ Ln a = a FC(µ)log (a)

where     FC(µ)log (a) = (-1)µ-1 (µ-1) !                                              if  µ is a positive integer
and       FC(µ)log (a) = [ Ln a - Psi(1-µ) - gamma ] / Gam(1-µ)        otherwise

Psi = Digamma Function , gamma = Euler's constant = 0.5772156649...  and  Gam = Gamma Function.

•  Exponential          Dµ Exp a = a 1F~1 ( 1 ; 1-µ ; a )

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing these routines )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R4n:  temp ( sometimes )

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flags:  F09 & F10 ( for "DSHA"  "DCHA"  "DSINA"  "DCOSA" )

Subroutines:  A*A1  "1/A"  "A^X"   "LNA"                            ( cf "Anions for the HP41" )
ST*A   ST/A   AMOVE   ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
"HGFB"                                              ( cf "Anionic Special Functions(II) for the HP41" )
GAM   1/G+   PSI                              ( cf "Gamma, Digamma Functions for the HP41" )

 01  LBL "DA^X"     02  STO M   03  X<>Y   04  -   05  STO N   06  XEQ "A^X"   07  RCL M   08  1   09  +   10  GAM   11  RCL N   12  1   13  +   14  1/G+   15  *   16  ST*A   17  X<>Y   18  CLA   19  RTN   20  LBL "DSHA"   21  CF 09   22  CF 10   23  GTO 00   24  LBL "DCHA"   25  SF 09   26  CF 10   27  GTO 00   28  LBL "DSINA"   29  CF 09   30  SF 10 31  GTO 00   32  LBL "DCOSA"   33  SF 09   34  SF 10   35  LBL 00   36  STO M   37  12   38  AMOVE   39  14   40  AMOVE   41  A*A1   42  4   43  FS?C 10   44  CHS   45  ST/A   46  CLX   47  SIGN   48  ENTER^   49  FC?C 09   50  ST+ X   51  RCL M   52  -   53  RCL X   54  1   55  +   56  2   57  ST/ Z   58  /   59  1.002   60  CHS 61  XEQ "HGFB"    62  X<> Z   63  ST+ X   64  1   65  -   66  STO M   67  41   68  AMOVE   69  X<>Y   70  XEQ "A^X"   71  2   72  RCL M   73  CHS   74  Y^X   75  PI   76  SQRT   77  *   78  ST*A   79  32   80  AMOVE   81  A*A1   82  CLA   83  RTN   84  LBL "DLNA"   85  STO M   86  12   87  AMOVE   88  X<>Y   89  CHS   90  XEQ "A^X" 91  STO N   92  RCL M   93  FRC   94  X#0?   95  GTO 01             96  LASTX   97  X<=0?   98  GTO 01   99  SIGN 100  CHS 101  RCL M 102  1 103  - 104  Y^X 105  LASTX 106  FACT 107  * 108  ST*A 109  RCL N 110  CLA 111  RTN 112  LBL 01 113  12 114  ASWAP 115  XEQ "LNA" 116  .5772156649 117  1 118  RCL M 119  - 120  PSI 121  ST+ Y 122  X<>Y 123  ST- 01 124  LASTX 125  1/G+ 126  ST*A 127  A*A1 128  CLA 129  RTN 130  LBL "DE^A" 131  14 132  AMOVE 133  SIGN 134  1 135  RCL Z 136  - 137  1.001 138  CHS 139  XEQ "HGFB"  140  R^ 141  41 142  AMOVE 143  X<>Y 144  CHS 145  XEQ "A^X" 146  32 147  AMOVE 148  A*A1 149  END

( 315 bytes )

 STACK INPUTS OUTPUTS INPUT OUTPUT Y µ / / / X x 1.nnn µ 1.nnn

Where   1.nnn   is the control number of the result f(a)     ;    x is required for "DA^X" only

Examples:    With  a = 1.1 + 1.2 i + 1.3 j + 1.4 k  and  µ = -PI                                        ( quaternion ->  4  STO 00 )

•  "DA^X"       if  x = e     f(a) =  1.630330082 + 0.233604201 i + 0.253071217 j + 0.272538234 k          ---Execution time = 6s---

•  "DSHA"       f(a) =  0.093634336 - 0.712248302 i - 0.771602327 j - 0.830956352 k          ---Execution time = 25s---

•  "DCHA"       f(a) =  -1.721108064 - 0.629274549 i - 0.681714094 j - 0.734153640 k          ---Execution time = 27s---

•  "DSINA"      f(a) =  -0.415438235 - 0.886075931 i - 0.959915592 j - 1.033755253 k          ---Execution time = 25s---

•  "DCOSA"    f(a) =  -2.980759335 - 0.217479859 i - 0.235603180 j - 0.253726502 k          ---Execution time = 26s---

•  "DLNA"       f(a) =  3.234305049 - 0.932622695 i - 1.010341253 j - 1.088059811 k          ---Execution time = 10s---

with µ = 3    f(a) = -0.123704743 + 0.014013365 i +0.015181146 j + 0.016348926 k          ---Execution time = 4s---

•  "DE^A"       f(a) =  -1.627473729 - 1.341522850 i - 1.453316422 j - 1.565109992 k          ---Execution time = 29s---

2°)  Kummer's Functions

Formula:               Dµ M(p;q;a) = a Gam(q) 2F~2 ( 1 , p ; 1-µ , q ; a )

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DKUMA" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R2n:  temp       R2n+1 ..... R3n = a

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flags:  /
Subroutines:   A*A1                                                              ( cf "Anions for the HP41" )
ST*A   DSA  AMOVE  ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
1/G+                                                ( cf "Gamma, Digamma Functions for the HP41" )
X+1   ST<>A                                    ( cf "A few M-Code Routines for the HP-41" )

-Lines 45-46 may be replaced by  RCL O  RCL N  RCL M
-Line 18 may be replaced by  PI  SIGN  +
-Line 02 may be replaced by  STO M  RDN  STO N  X<>Y  STO O

 01  LBL "DKUMA"  02  ST<>A  03  12  04  AMOVE  05  CLX  06  ST*A  07  SIGN  08  STO 01  09  13  10  AMOVE 11  CLX  12  STO P  13  LBL 01              14  A*A1  15  RCL N  16  RCL O  17  CHS  18  X+1  19  RCL M  20  RCL P 21  ST+ T  22  ST+ Z  23  +  24  *  25  /  26  ISG P                 27  CLX  28  ST*A  29  DSA  30  X#0? 31  GTO 01  32  21  33  MOVE  34  RCL O  35  CHS  36  XEQ "A^X"       37  RCL O  38  CHS  39  X+1  40  1/G+ 41  ST*A  42  23  43  ASWAP             44  A*A1  45  RDN  46  ST<>A  47  R^  48  CLA  49  END

( 93 bytes / SIZE 3n+1 )

 STACK INPUTS OUTPUTS T / µ Z µ p Y p q X q 1.nnn

Where   1.nnn   is the control number of the result

Examples:      µ = PI     p = sqrt(2)   q = sqrt (3)     a = 0.4 + 0.5 i + 0.6 j + 0.7 k                        ( quaternion ->  4  STO 00 )

0.4  STO 01
0.5  STO 02
0.6  STO 03
0.7  STO 04

PI    2   SQRT
3   SQRT    XEQ "DKUMA"  >>>>   1.004                                ---Execution time = 38s---

R01 = 0.630946964
R02 = 0.307230135
R03 = 0.368676162
R04 = 0.430122188

-Whence      F(a) =  0.630946964 + 0.307230135 i + 0.368676162 j + 0.430122188 k

Note:

-This program does not work if  µ = 1 , 2 , 3 , .....

3°)  Error Function

Formula:            Dµ Erf (a) = 2µ a1-µ2F~2 [ 1/2 , 1 ; (2-µ)/2 , (3-µ)/2 ; -a2 ]

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DERFA" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R3n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flags:  /
Subroutines:   A*A1   "A^X"                                                   ( cf "Anions for the HP41" )
ST*A   DSA  AMOVE  ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
1/G+                                                ( cf "Gamma, Digamma Functions for the HP41" )

-Line 62 may be replaced with  E3   /

 01  LBL "DERFA"  02  STO M  03  .5  04  STO N  05  *  06  STO O  07  12  08  AMOVE  09  SIGN  10  CHS  11  ST*A  12  A*A1  13  12 14  ASWAP  15  SIGN  16  RCL M  17  -  18  XEQ "A^X"     19  SIGN  20  RCL O  21  -  22  1/G+  23  RCL N  24  LASTX  25  +  26  1/G+ 27  *  28  2  29  RCL M  30  Y^X  31  *  32  ST*A  33  13  34  AMOVE         35  CLX  36  STO P  37  LBL 01  38  A*A1  39  RCL N 40  ENTER^  41  SIGN  42  RCL O  43  -  44  RCL X  45  RCL N  46  +  47  RCL P              48  ST+ T  49  ST+ Z  50  +  51  *  52  / 53  ST*A  54  ISG P  55  CLX  56  DSA  57  X#0?  58  GTO 01  59  31  60  AMOVE         61  RCL 00  62  X/E3  63  ISG X  64  CLA  65  END

( 115 bytes / SIZE 3n+1 )

 STACK INPUT OUTPUT X µ 1.nnn

Where   1.nnn   is the control number of the result

Examples:            µ = PI   ,    a =  0.4 + 0.5 i + 0.6 j + 0.7 k                        ( quaternion ->  4  STO 00 )

0.4  STO 01
0.5  STO 02
0.6  STO 03
0.7  STO 04

PI   XEQ "DERFA"  >>>>   1.004                                ---Execution time = 36s---

R01 =  0.924128302
R02 = 10.81672321
R03 = 12.98006783
R04 = 15.14341247

-Whence,  Dµ Erf (a) = 0.924128302 + 10.81672321 i + 12.98006783 j + 15.14341247 k

Note:

-This program does not work if  µ = 2 , 3 , 4 , .....

4°)  Sine & Hyperbolic Sine Integral

Formulae:

•  Sine Integral             Dµ Si a = 2µ-2  PI  a1-µ2F~3 ( 1/2 , 1 ; 3/2 , (2-µ)/2 , (3-µ)/2 ; -a2/4 )

•  Hyperbolic Sine Integral         Dµ Shi x = 2µ-2  PI  a1-µ2F~3 ( 1/2 , 1 ; 3/2 , (2-µ)/2 , (3-µ)/2 ; a2/4 )

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DSIA" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R3n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flag:    F01        CF 01  <->  Sine Integral
SF 01   <->  Hyperbolic Sine Integral

Subroutines:   A*A1   "A^X"                                                               ( cf "Anions for the HP41" )
ST*A  ST/A   DSA  AMOVE  ASWAP       ( cf "Anionic Special Functions(I) for the HP41" )
1/G+                                                             ( cf "Gamma, Digamma Functions for the HP41" )
X+1                                                                ( cf "A few M-Code Routines for the HP-41" )

Line 66 may be replaced with  E3  /
X+1  ( line 49 ) may be replaced by  ISG X   CLX

 01  LBL "DSIA"  02  CLA  03  STO M  04  12  05  AMOVE  06  4  07  FC? 01  08  CHS  09  ST/A  10  A*A1  11  12  12  ASWAP  13  SIGN  14  RCL M 15  -  16  XEQ "A^X"  17  2  18  RCL M  19  .5  20  CHS  21  ST* M  22  STO N  23  SIGN  24  +  25  Y^X  26  PI  27  SQRT  28  * 29  1  30  RCL M         31  +  32  1/G+  33  *  34  1.5  35  RCL M  36  +  37  1/G+  38  *  39  ST*A  40  13  41  AMOVE  42  LBL 01 43  A*A1  44  RCL N  45  ENTER^       46  CHS  47  RCL M  48  RCL O  49  X+1  50  STO O  51  ST+ T  52  ST+ Z  53  +  54  ST* Y  55  RCL N  56  - 57  *  58  /  59  ST*A  60  DSA  61  X#0?  62  GTO 01  63  31  64  AMOVE      65  RCL 00  66  X/E3  67  ISG X  68  CLA  69  END

( 121 bytes / SIZE 3n+1 )

 STACK INPUT OUTPUT X µ 1.nnn

Where   1.nnn   is the control number of the result

Examples:            µ = PI   ,    a =  1 + 2 i + 3 j + 4 k                        ( quaternion ->  4  STO 00 )

CF 01     Dµ Si(a) = -6.216067679 + 4.558034704 i + 6.837052053 j + 9.116069403 k                     ---Execution time = 32s---

SF 01     Dµ Shi(a) = -0.078381628 - 0.067284573 i - 0.100926858 j - 0.134569145 k                      ---Execution time = 33s---

Note:

-This routine does not work for  µ = 2 , 3 , 4 , .....

5°)  Cosine & Hyperbolic Cosine Integral

Formulae:

•  Cosine Integral

Dµ Ci a =  [ (-1)µ-1 (µ-1) ! ] a - 2µ-3  sqrt(PI)  a2-µ 2F~3 ( 1 , 1 ; 2 , (3-µ)/2 , (4-µ)/2 ; -a2/4 )       if  µ is a positive integer

Dµ Ci a =  [ Ln a - Psi(1-µ) ] a - 2µ-3  sqrt(PI)  a2-µ  2F~3 ( 1 , 1 ; 2 , (3-µ)/2 , (4-µ)/2 ; -a2/4 )       otherwise

•  Hyperbolic Cosine Integral

Dµ Chi a = [ (-1)µ-1 (µ-1) ! ] a + 2µ-3  sqrt(PI)  a2-µ 2F~3 ( 1 , 1 ; 2 , (3-µ)/2 , (4-µ)/2 ; a2/4 )       if  µ is a positive integer

Dµ Chi a = [ Ln a - Psi(1-µ) ] a + 2µ-3 sqrt(PI)  a2-µ  2F~3 ( 1 , 1 ; 2 , (3-µ)/2 , (4-µ)/2 ; a2/4 )      otherwise

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DCIA" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R4n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flag:   F01        CF 01  <->  Cosine Integral
SF 01   <->  Hyperbolic Cosine Integral

Subroutines:   A*A1   A+A   "A^X"  "LNA"                                       ( cf "Anions for the HP41" )
ST*A   ST/A   DSA  AMOVE  ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
1/G+     PSI                                                ( cf "Gamma, Digamma Functions for the HP41" )

-Line 35 may be replaced by  ISG X   CLX
-Line 103 ------------------   1  -

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

( 219 bytes / SIZE 4n+1 )

 STACK INPUT OUTPUT X µ 1.nnn

Where   1.nnn   is the control number of the result

Examples:            µ = PI   ,    a =  1.1 + 1.2 i + 1.3 j + 1.4 k                        ( quaternion ->  4  STO 00 )

CF 01     Dµ Ci(a) = 0.770099307 + 0.231600213 i + 0.250900230 j + 0.270200248 k                     ---Execution time = 43s---

SF 01     Dµ Chi(a) = -0.154706962 + 0.244319821 i + 0.264679806 j + 0.285039791 k                      ---Execution time = 48s---

Note:

-This program does not work for  µ = 3 , 4 , 5 , .....

6°)  Exponential Integral Ei

Formulae:

Dµ Ei a = [ (-1)µ-1 (µ-1) ! ] a + a1-µ  2F~2 ( 1 , 1 ; 2 , 2-µ ; a )       if µ is a positive integer

Dµ Ei a = [ Ln a - Psi(1-µ) ] a + a1-µ  2F~2 ( 1 , 1 ; 2 , 2-µ ; a )     otherwise

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DEIA" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R3n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flags:  /
Subroutines:   A*A1   A+A   "A^X"                                        ( cf "Anions for the HP41" )
ST*A   DSA  AMOVE  ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
1/G+     PSI                                     ( cf "Gamma, Digamma Functions for the HP41" )

-Line 63 may be replaced by  1  -
-Lines 79 to 81 may be replaced by  LASTX

 01  LBL "DEIA"  02  CLA  03  STO M  04  12  05  AMOVE  06  CLX  07  ST*A  08  SIGN  09  STO 01  10  13  11  AMOVE  12  LBL 01  13  A*A1  14  SIGN  15  ENTER^  16  ST+ X  17  RCL X  18  RCL M 19  -  20  RCL N  21  ST+ T  22  ST+ Z  23  +  24  *  25  /  26  ST*A  27  ISG N  28  CLX  29  DSA  30  X#0?  31  GTO 01  32  21  33  AMOVE      34  SIGN  35  RCL M  36  - 37  XEQ "A^X"  38  23  39  ASWAP  40  A*A1  41  2  42  RCL M  43  -  44  1/G+  45  ST*A  46  13  47  ASWAP       48  12  49  AMOVE  50  RCL M  51  CHS  52  XEQ "A^X"  53  RCL M  54  FRC 55  X#0?  56  GTO 02  57  LASTX  58  X<=0?  59  GTO 02  60  1  61  CHS  62  LASTX  63  X-1  64  Y^X  65  LASTX  66  FACT  67  *  68  ST*A  69  GTO 03  70  LBL 02  71  12  72  ASWAP 73  XEQ "LNA"  74  SIGN  75  RCL M  76  -  77  PSI  78  ST- 01  79  1  80  RCL M  81  -  82  1/G+  83  ST*A  84  A*A1  85  LBL 03  86  32  87  AMOVE  88  A+A  89  CLA  90  END

( 160 bytes / SIZE 3n+1 )

 STACK INPUT OUTPUT X µ 1.nnn

Where   1.nnn   is the control number of the result

Example:            µ = PI   ,    a =  1.1 + 1.2 i + 1.3 j + 1.4 k                        ( quaternion ->  4  STO 00 )

Dµ Ei(a) = -0.272860689 + 0.388792411 i + 0.421191778 j + 0.453591146 k                     ---Execution time = 51s---

Note:

-"DEIA" does not work for  µ = 2 , 3 , 4 , .....

7°)  Generalized Laguerre's Functions

Formula:

Dµ Lpq (a) = [ Gam(p+q+1)/Gam(q+1) ] a2F~2 ( 1 , -q ; p+1 , 1-µ ; a )

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DLANA" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R3n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flags:  /
Subroutines:   A*A1   "A^X"                                     ( cf "Anions for the HP41" )
ST*A   DSA  AMOVE      ( cf "Anionic Special Functions(I) for the HP41" )
1/G+                                  ( cf "Gamma, Digamma Functions for the HP41" )

-Line 52 may be replaced with  E3  /
-Line 38 may be replaced by  ISG P  CLX
-Line 03 may be replaced with  STO M  RDN  STO N  X<>Y  STO O
-Lines 02 & 22  may be replaced by  1  +

 01  LBL "DLANA"  02  X+1  03  ST<>A  04  12  05  AMOVE  06  RCL O  07  CHS  08  XEQ "A^X"  09  RCL M  10  1/G+  11  LASTX 12  RCL N  13  +  14  1/G+  15  /  16  1  17  RCL O             18  -  19  1/G+  20  *  21  RCL N  22  X+1 23  1/G+  24  *  25  ST*A  26  13  27  AMOVE           28  CLX  29  STO P  30  LBL 01  31  A*A1  32  RCL M  33  CHS 34  RCL O             35  CHS  36  RCL N  37  RCL P  38  X+1  39  STO P  40  ST+ T  41  ST+ Z  42  +  43  *  44  / 45  ST*A  46  DSA  47  X#0?  48  GTO 01  49  31  50  AMOVE           51  RCL 00  52  X/E3  53  ISG X  54  CLA  55  END

( 103 bytes / SIZE 3n+1 )

 STACK INPUTS OUTPUTS Z µ / Y p / X q 1.nnn

Where   1.nnn   is the control number of the result

Examples:      µ = -PI     p = sqrt(2)   q = sqrt (3)     a = 1 + 2 i + 3 j + 4 k                        ( quaternion ->  4  STO 00 )

1  STO 01
2  STO 02
3  STO 03
4  STO 04

PI   CHS   2  SQRT  3  SQRT  XEQ "DLANA"  >>>>  1.004                             ---Execution time = 34s---

R01 = -108.7414887
R02 =  -11.05061335
R03 =  -16.57592004
R04 =  -22.10122675

-Whence,   Dµ Lpq (a) = -108.7414887 - 11.05061335 i - 16.57592004 j - 22.10122675 k

Note:

µ must be different from  1 , 2 , 3 , .....

8°)  Bessel & Modified Bessel Functions of the 1st Kind

Formulae:

•    Dµ Jm (a) = 2µ-2m sqrt(PI)  a-m-µ  Gam(m+1)  2F~3 [ (m+1)/2 , (m+2)/2 ; (m+1-µ)/2 , (m+2-µ)/2 , m+1 ; -a2/4 ]

•    Dµ Im (a) = 2µ-2m sqrt(PI) am-µ  Gam(m+1)  2F~3 [ (m+1)/2 , (m+2)/2 ; (m+1-µ)/2 , (m+2-µ)/2 , m+1 ; a2/4 ]

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DJNA" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R3n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flag:  F01    CF 01  <->   Bessel functions of the 1st kind
SF 01  <->   Modified Bessel Functions of the 1st kind

Subroutines:   A*A1   "A^X"                                                             ( cf "Anions for the HP41" )
ST*A   ST/A   DSA  AMOVE  ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
1/G+                                                           ( cf "Gamma, Digamma Functions for the HP41" )

-Line 82 may be replaced with  E3  /
-Line 70 may be replaced by  ISG X  CLX
-Line 29 may be replaced by  1  +

 01  LBL "DJNA"  02  CLA  03  STO N  04  X<>Y  05  STO M  06  12  07  AMOVE  08  4  09  FC? 01  10  CHS  11  ST/A  12  A*A1  13  12  14  ASWAP  15  RCL N  16  RCL M  17  - 18  XEQ "A^X"   19  PI  20  SQRT  21  2  22  RCL N  23  ST+ X  24  RCL M  25  -  26  Y^X  27  /  28  RCL N  29  X+1  30  STO N  31  RCL M  32  -  33  2  34  / 35  STO M  36  1/G+  37  *  38  RCL M  39  .5  40  +  41  1/G+  42  *  43  ST*A  44  13  45  AMOVE        46  .5  47  STO P  48  LBL 01  49  A*A1  50  RCL N  51  RCL P 52  *  53  STO Y  54  LASTX  55  +  56  RCL M  57  RCL O  58  ST+ T  59  ST+ Z  60  +  61  /  62  *  63  RCL M          64  RCL P  65  +  66  RCL N  67  RCL O  68  ST+ Z 69  ST+ Y  70  X+1  71  STO O  72  *  73  *  74  /  75  ST*A  76  DSA  77  X#0?  78  GTO 01  79  31  80  AMOVE       81  RCL 00  82  X/E3  83  ISG X  84  CLA  85  END

( 149 bytes / SIZE 3n+1 )

 STACK INPUTS OUTPUTS Y µ / X p 1.nnn

Where   1.nnn   is the control number of the result

Examples:      µ = PI     p = sqrt(2)        a = 1.1 + 1.2 i + 1.3 j + 1.4 k                        ( quaternion ->  4  STO 00 )

•  CF 01   Dµ Jm (a) = -1.115289298 + 0.405457582 i + 0.439245714 j + 0.473033847 k                 ---Execution time = 30s---

•  SF 01    Dµ Im (a) = -0.213651981 + 0.178633476 i + 0.193519599 j + 0.208405722 k                 ---Execution time = 30s---

9°)  Bessel & Modified Bessel Functions of the 2nd Kind - Non-Integer Order

Formulae:

•   Dµ Ym (a) = 2µ-2m (PI)1/2  a-µ-m csc(m.PI) { -16m Gam(1-m)  2F~3 [ (1-m)/2 , (2-m)/2 ; (1-µ-m)/2 , (2-µ-m)/2 , 1-m ; -a2/4 ]

+ x2m Cos(m.PI) Gam(m+1)  2F~3 [ (m+1)/2 , (m+2)/2 ; (m+1-µ)/2 , (m+2-µ)/2 , m+1 ; -a2/4 ] }    where m is not an integer.

•   Dµ Km (a) = 2µ-2m-1 (PI)3/2  a-µ-m csc(m.PI) { 16m Gam(1-m)  2F~3 [ (1-m)/2 , (2-m)/2 ; (1-µ-m)/2 , (2-µ-m)/2 , 1-m ; x2/4 ]

- x2m Gam(m+1)  2F~3 [ (m+1)/2 , (m+2)/2 ; (m+1-µ)/2 , (m+2-µ)/2 , m+1 ; a2/4 ] }    where m is not an integer.

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DYNA" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R5n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flag:  F01    CF 01  <->   Bessel functions of the 2nd kind
SF 01  <->   Modified Bessel Functions of the2nd kind

Subroutines:   A*A1   A-A   "A^X"                                                     ( cf "Anions for the HP41" )
ST*A   ST/A   DSA  AMOVE  ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
1/G+                                                           ( cf "Gamma, Digamma Functions for the HP41" )

-Line 124 may be replaced with  2  /
-Line 54 may be replaced by  ISG X  CLX
-Line 34 may be replaced by  PI  SIGN  +

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

( 263 bytes / SIZE 5n+1 )

 STACK INPUTS OUTPUTS Y µ / X m 1.nnn

Where   1.nnn   is the control number of the result  and  µ is not an integer

Examples:      µ = PI     m = sqrt(2)        a = 1.1 + 1.2 i + 1.3 j + 1.4 k                        ( quaternion ->  4  STO 00 )

•  CF 01   Dµ Ym (a) = -0.710530741 - 0.501187143 i - 0.542952739 j - 0.584718334 k                 ---Execution time = 72s---

•  SF 01    Dµ Km (a) =  0.019073486 - 0.197769632 i - 0.214250435 j - 0.230731237 k                 ---Execution time = 70s---

10°) Spherical Bessel Functions of the 1st Kind

Formula:

Dµ jn (a) = 2µ-2m-1 PI  am-µ  Gam(m+1)  2F~3 [ (m+1)/2 , (m+2)/2 ; (m+1-µ)/2 , (m+2-µ)/2 , m+3/2 ; -a2/4 ]

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DBS1A" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R3n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flags:  /
Subroutines:   A*A1   "A^X"                                                           ( cf "Anions for the HP41" )
ST*A   ST/A   DSA  AMOVE  ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
1/G+                                                           ( cf "Gamma, Digamma Functions for the HP41" )

-Line 96 may be replaced by  E3  /
-Line 91 may be replaced with  1  +
-Lines 62 to 67  may be replaced with  SIGN  RCL N  +  RCL M  -  RCL X  1  +
-Line 49 may be replaced by   ISG X   CLX
-Lines 25 to 27 may be replaced by  PI  SIGN  +  RCL X   LASTX   +

 01  LBL "DBS1A"  02  CLA  03  STO N  04  X<>Y  05  STO M  06  12  07  AMOVE  08  A*A1  09  4  10  CHS  11  ST/A  12  12  13  ASWAP  14  RCL N  15  RCL M  16  -  17  XEQ "A^X"  18  13  19  AMOVE  20  .5 21  STO P  22  LBL 01            23  A*A1  24  RCL N  25  X+1  26  ENTER^  27  X+1  28  RCL Y  29  RCL M  30  -  31  RCL P  32  ST* T  33  ST* Z  34  *  35  RCL O  36  ST+ T  37  ST+ Z  38  +  39  /  40  * 41  RCL N  42  RCL M            43  -  44  RCL N  45  RCL P  46  ST* Z  47  +  48  RCL O  49  X+1  50  STO O  51  ST+ Z  52  ST+ Y  53  *  54  *  55  /  56  ST*A  57  DSA  58  X#0?  59  GTO 01  60  31 61  AMOVE          62  RCL N  63  X+1  64  RCL M  65  -  66  ENTER^  67  X+1  68  2   69  ST/ Z  70  /  71  1/G+  72  X<>Y  73  1/G+  74  *  75  RCL N  76  1.5  77  +    78  1/G+  79  *  80  PI 81  *  82  2  83  ST/ Y  84  RCL M  85  RCL N  86  ST+ X  87  -  88  Y^X  89  *  90  RCL N            91  X+1  92  1/G+  93  /  94  ST*A  95  RCL 00  96  X/E3  97  ISG X  98  CLA  99  END

( 175 bytes / SIZE 3n+1 )

 STACK INPUTS OUTPUTS Y µ / X m 1.nnn

Where   1.nnn   is the control number of the result

Example:      µ = PI     m = sqrt(2)        a = 1.1 + 1.2 i + 1.3 j + 1.4 k                        ( quaternion ->  4  STO 00 )

Dµ jm (a) = -0.513818255 + 0.162870796 i + 0.176443362 j + 0.190015928 k                   ---Execution time = 37s---

11°) Hermite Functions

Formula:

Dµ Hm (a) = [ 2m+µ (PI) a / Gam((1-m)/2) ] 2F~2 [ 1 , -m/2 ; (1-µ)/2 , (2-µ)/2 ; a2 ]
- [ 2m+µ (PI) a1-µ / Gam((-m)/2) ] 2F~2 [ 1 , (1-m)/2 ; 1-µ/2 , (3-µ)/2 ; a2 ]

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DHMTA" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R4n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flag:   F10
Subroutines:   A*A1   A-A   "A^X"                                                 ( cf "Anions for the HP41" )
ST*A   ST/A   DSA  AMOVE  ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
1/G+                                                           ( cf "Gamma, Digamma Functions for the HP41" )

-Line 66 may be replaced by  1  -
-Lines 82-72-67-58-55 may be replaced with  2  /
-Lines 28 to 30 may be replaced by  PI  SIGN  +  RCL X   LASTX   +

 01  LBL "DHMTA"  02  CHS  03  STO N  04  X<>Y  05  CHS  06  STO M  07  12  08  AMOVE  09  14  10  AMOVE  11  A*A1  12  12  13  ASWAP  14  SF 10  15  LBL 01  16  RCL M  17  XEQ "A^X"  18  13  19  AMOVE  20  CLX 21  STO O  22  .5  23  STO P  24  LBL 02  25  A*A1  26  RCL N  27  RCL M              28  X+1  29  ENTER^  30  X+1  31  RCL P  32  ST* T  33  ST* Z  34  *  35  RCL O  36  ST+ T  37  ST+ Z  38  +  39  *  40  / 41  ST*A  42  ISG O  43  CLX  44  DSA  45  X#0?  46  GTO 02  47  31  48  AMOVE           49  FC?C 10  50  RTN  51  SIGN  52  ST+ M  53  ST+ N  54  RCL M  55  X/2  56  1/G+  57  RCL N  58  X/2  59  1/G+  60  * 61  ST*A  62  14  63  ASWAP            64  XEQ 01  65  RCL N  66  X-1  67  X/2  68  1/G+  69  2  70  RCL M  71  +  72  X/2  73  1/G+  74  *  75  ST*A  76  42  77  AMOVE  78  A-A  79  STO O  80  SIGN 81  RCL M             82  +  83  X/2  84  1/G+  85  2  86  RCL M  87  RCL N  88  +  89  2  90  -  91  Y^X  92  /  93  PI  94  *  95  CHS  96  ST*A  97  RCL O  98  CLA  99  END

( 179 bytes / SIZE 4n+1 )

 STACK INPUTS OUTPUTS Y µ / X m 1.nnn

Where   1.nnn   is the control number of the result

Example:      µ = -PI     m = sqrt(2)        a = 1.1 + 1.2 i + 1.3 j + 1.4 k                        ( quaternion ->  4  STO 00 )

Dµ Hm (a) = 2.770233592 - 1.881749624 i - 2.038562092 j - 2.195374561 k              ---Execution time = 100s---

12°) Fresnel Integrals

Formulae:

•  Fresnel Cosine Integral     Dµ C(a) = 22µ-3/2 PI3/2 a1-µ3F~4 [ 1/4 , 3/4 , 1 ; (2-µ)/4 , (3-µ)/4 , (4-µ)/4 , (5-µ)/4 ; -(PI)2 a4/16 ]

•  Fresnel Sine Integral    Dµ S(a) = 22µ-11/2 PI5/2 a3-µ3F~4 [ 3/4 , 1 , 5/4 ; (4-µ)/4 , (5-µ)/4 , (6-µ)/4 , (6-µ)/4 ; -(PI)2 a4/16 ]

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DCSA" )

•  R01   ......  •  Rnn = the n components of the anion

Rn+1 ..........  R3n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of the result

Flag:  F01    CF 01  <->   Fresnel Cosine Integral     Dµ C(a)
SF 01  <->   Fresnel Sine Integral          Dµ S(a)

Subroutines:   A*A1   A-A   "A^X"                                                     ( cf "Anions for the HP41" )
ST*A   ST/A   DSA  AMOVE  ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
1/G+                                                           ( cf "Gamma, Digamma Functions for the HP41" )

Lines 61 to 65 may be replaced by  RCL X  PI  SIGN  +  ENTER^  X<> L  ST+ L

 01  LBL "DCSA"   02  CLA   03  STO M   04  12   05  AMOVE   06  13   07  AMOVE   08  A*A1   09  12   10  AMOVE   11  A*A1   12  PI   13  4   14  STO N   15  /   16  X^2   17  CHS   18  ST*A   19  12   20  AMOVE   21  31   22  AMOVE   23  3   24  FC? 01   25  SIGN   26  RCL M   27  -   28  XEQ "A^X" 29  13   30  AMOVE         31  LBL 01   32  A*A1   33  SIGN   34  PI   35  INT   36  2   37  FS? 01   38  X^2   39  RCL M   40  -   41  RCL N   42  FS? 01   43  ST+ T   44  ST/ T   45  ST/ Z   46  /   47  RCL O   48  ST+ T   49  ST+ Z   50  +   51  /   52  *   53  PI   54  INT   55  2   56  FC? 01 57  CLX   58  +   59  RCL M           60  -   61  ENTER^   62  X+1   63  ENTER^   64  X+1   65  SIGN   66  CLX   67  RCL N   68  ST/ Z   69  ST/ Y   70  ST/ L   71  CLX   72  RCL O   73  ST+ Z   74  ST+ Y   75  ST+ L   76  X<> L   77  *   78  *   79  /   80  ST*A   81  ISG O   82  CLX   83  DSA    84  X#0? 85  GTO 01   86  31   87  AMOVE         88  PI   89  X^2   90  LASTX   91  3   92  Y^X   93  FS? 01   94  *   95  SQRT   96  2   97  RCL M   98  2   99  FC? 01 100  CLX 101  - 102  ST+ X 103  1.5 104  - 105  Y^X 106  * 107  2 108  FS? 01 109  X^2 110  RCL M 111  - 112  RCL N 113  / 114  1/G+ 115  ST* Y 116  X<> L 117  RCL N         118  1/X 119  + 120  1/G+ 121  ST* Y 122  X<> L 123  RCL N 124  1/X 125  + 126  1/G+ 127  ST* Y 128  X<> L 129  RCL N 130  1/X 131  + 132  1/G+ 133  * 134  ST*A 135  RCL 00 136   E3 137  / 138  ISG X 139  CLA 140  END

( 227 bytes / SIZE 3n+1 )

 STACK INPUT OUTPUT X µ 1.nnn

Where   1.nnn   is the control number of the result

Example:      µ = -PI    ,     a = 1.1 + 1.2 i + 1.3 j + 1.4 k                        ( quaternion ->  4  STO 00 )

•  CF 01   Dµ C (a) =   0.418951593 - 0.420163597 i - 0.455177231 j - 0.490190864 k                 ---Execution time = 62s---

•  SF 01    Dµ S (a) =  -0.233950973 + 0.363367171 i + 0.393647768 j + 0.423928366 k                 ---Execution time = 63s---

Notes:

-With  Dµ C (a) ,  µ must be different from  2 , 3 , 4 , ....
-With  Dµ S (a) ,  µ ----------------------   4 , 5 , 6 , ....

13°) Airy Functions

Formulae:

Dµ Ai(a) = 3µ-4/3 a  { 32/3 Gam(1/3) 2F~3 [ 1/3 , 1 ; (1-µ)/3 , (2-µ)/3 , (3-µ)/3 ; a3/9 ] - a Gam(2/3) 2F~3 [ 2/3 , 1 ; (4-µ)/3 , (2-µ)/3 , (3-µ)/3 ; a3/9 ] }

Dµ Bi(a) = 3µ-5/6 a  { 32/3 Gam(1/3) 2F~3 [ 1/3 , 1 ; (1-µ)/3 , (2-µ)/3 , (3-µ)/3 ; a3/9 ] + a Gam(2/3) 2F~3 [ 2/3 , 1 ; (4-µ)/3 , (2-µ)/3 , (3-µ)/3 ; a3/9 ] }

"DAIRYA"  calculates  Dµ Ai(a)  &  Dµ Bi(a)  simultaneously and stores the results in R01 thru R2n

Data Registers:           •  R00 = n > 1         ( Registers R00 thru Rnn are to be initialized before executing "DAIRYA" )

•  R01   ......  •  Rnn = the n components of the anion

R2n+1 ..........  R4n:  temp

>>>  When the programs stop:     R01   ......  Rnn = the n components of  Dµ Ai(a)
Rn+1 .....  R2n =  -------------------    Dµ Bi(a)

Flag:  F10
Subroutines:   A*A1   A-A   A+A   "A^X"                                         ( cf "Anions for the HP41" )
ST*A   ST/A   DSA  AMOVE  ASWAP    ( cf "Anionic Special Functions(I) for the HP41" )
1/G+                                                           ( cf "Gamma, Digamma Functions for the HP41" )

-Lines 125 to 127 may be replaced by  4  3  /   RCL M   X<>Y
-Lines 73-74 may be replaced by   SIGN  ST* X  ST- L
-Lines 47 to 50 may be replaced by  CLX  FS? 10  SIGN  RCL M  -
-Lines 43-33-27-21 may be replaced by  3  /

 01  LBL "DAIRYA"   02  STO M   03  3   04  STO N   05  12   06  AMOVE   07  A*A1   08  A*A1   09  9   10  ST/A   11  12   12  ASWAP   13  14     14  AMOVE   15  SIGN   16  ENTER^   17  SF 10   18  LBL 01   19  FS? 10   20  ST+ X   21  X/3   22  1/G+   23  /   24  2   25  RCL M   26  -   27  X/3 28  1/G+   29  *   30  3   31  RCL M   32  -   33  X/3   34  1/G+   35  *   36  1   37  RCL N   38  FC? 10   39  CLX   40  +   41  RCL M   42  -   43  X/3   44  1/G+     45  *   46  STO O   47  RCL M   48  CHS   49  FS? 10   50  X+1   51  XEQ "A^X"         52  CLX   53  X<> O   54  ST*A 55  13   56  AMOVE             57  LBL 02   58  A*A1   59  SIGN   60  FS? 10   61  ST+ X   62  1   63  RCL N   64  FC? 10   65  CLX   66  +   67  RCL M   68  -   69  RCL N   70  LASTX   71  -   72  ENTER^   73  X-1   74  SIGN   75  CLX   76  RCL N   77  ST/ T   78  ST/ Z   79  ST/ Y   80  ST/ L   81  CLX 82  RCL O   83  ST+ T   84  ST+ Z   85  ST+ Y   86  ST+ L   87  X<> L   88  *   89  *   90  /   91  ST*A   92  ISG O   93  CLX   94  DSA   95  X#0?    96  GTO 02   97  31   98  AMOVE             99  FC?C 10 100  RTN 101  14 102  ASWAP 103  RCL N 104  2 105  RCL Y 106  / 107  Y^X 108  1 109  XEQ 01 110  42 111  AMOVE 112  A+A 113  3 114  RCL M 115  5 116  6 117  / 118  - 119  Y^X 120  ST*A 121  13 122  ASWAP           123  A-A 124  3 125  RCL M 126  4 127  X/3 128  - 129  Y^X 130  ST*A 131  32 132  AMOVE 133  R^ 134  CLA 135  END

( 231 bytes / SIZE 4n+1 )

 STACK INPUT OUTPUT X µ 1.nnn

Where   1.nnn   is the control number of Dµ Ai(a)                       [  Dµ Bi(a)  is in registers  Rn+1 thru R2n ]

Example:      µ = PI    ,     a = 1.1 + 1.2 i + 1.3 j + 1.4 k                        ( quaternion ->  4  STO 00 )

•  R01 thru R04  ->     Dµ Ai (a) =  -0.264935707 + 0.693505589 i + 0.751297721 j + 0.809089853 k                 ---Execution time = 86s---

•  R05 thru R08  ->     Dµ Bi (a) =  -2.400100467 - 0.665074361 i - 0.709663892 j - 0.764253423 k

Note:

"DAIRYA" does not work if  µ = 1 , 2 , 3 , ....

Reference: