hp41programs

Anger-Weber

Anger & Weber Functions for the HP-41


Overview
 

 1°)  Ascending Series
 2°)  Asymptotic Expansions
 

-Anger's function is defined by     Jn(x) = (1/PI) §0pi  cos ( n.t - x.sin t ) dt              § = integral
-Weber's function is defined by    En(x) = (1/PI) §0pi  sin ( n.t - x.sin t ) dt
 

1°)  Ascending Series
 

Formulae:

   Jn(x) = + (x/2) sin ( 90°n )  1F2( 1 ; (3-n)/2 , (3+n)/2 ; -x2/4 ) / Gam((3-n)/2) / Gam ((3+n)/2)
               + cos ( 90°n ) 1F2( 1 ; (2-n)/2 , (2+n)/2 ; -x2/4 ) / Gam((2-n)/2) / Gam ((2+n)/2)

   En(x) = - (x/2) cos ( 90°n )  1F2( 1 ; (3-n)/2 , (3+n)/2 ; -x2/4 ) / Gam((3-n)/2) / Gam ((3+n)/2)
               + sin ( 90°n ) 1F2( 1 ; (2-n)/2 , (2+n)/2 ; -x2/4 ) / Gam((2-n)/2) / Gam ((2+n)/2)
 

Data Registers:  R00 to R07: temp
Flags: /
Subroutines:   "1F2"  ( cf "Hypergeometric Functions" ) ,  "GAM"  ( cf "Gamma Function" )
 
 

01  LBL "ANWEB"
02  2
03  ST/ Z
04  /
05  STO 04
06  X<>Y
07  STO 03
08  STO 05
09  CHS
10  1
11  STO 01
12  ST+ 03
13  +
14  STO 02
15  RCL 04
16  X^2
17  CHS
18  XEQ "1F2"
19  STO 06
20  RCL 02
21  XEQ "GAM"    
22  ST/ 06
23  RCL 03
24  XEQ "GAM"
25  ST/ 06
26  2
27  1/X
28  ST+ 02
29  ST+ 03
30  RCL 02
31  XEQ "GAM"    
32  STO 07
33  RCL 03
34  XEQ "GAM"
35  ST* 07
36  RCL 00
37  XEQ "1F2"      
38  RCL 07
39  /
40  RCL 04
41  *
42  1
43  CHS
44  ACOS
45  RCL 05
46  *
47  STO 05
48  X<>Y
49  P-R
50  RCL 05           
51  RCL 06
52  P-R
53  X<> Z
54  -
55  X<> Z
56  +
57  END

 
    ( 100 bytes / SIZE 008 )
 
 

      STACK        INPUTS      OUTPUTS
           Y             n         En(x)
           X             x         Jn(x)

 
Example:

   2   SQRT
  PI  XEQ "ANWEB"  >>>>    Jsqrt(2)(PI) =  0.366086558                                ---Execution time = 26s---
                                  X<>Y    Esqrt(2)(PI) = -0.315594385
 

Notes:

-If n is an integer,  Jn(x) = Jn(x) = Bessel function of the first kind.
-However, this routine doesn't work if n is an integer, except if n = -1 , 0 , +1
 but it will work if we use the regularized hypergeometric function  1F~2
-The variant hereunder calls "HGF+"
 

Data Registers:  R00 to R06: temp
Flags: /
Subroutine:   "HGF+"  ( cf "Hypergeometric Functions" )

-If you want to use the M-Code routine HGF+ , replace line 30 by  HGF+  and line 20 by  STO 00   HGF+
 
 

01  LBL "ANWEB"
02  2
03  ST/ Z
04  /
05  STO 04
06  X<>Y
07  STO 03
08  STO 05
09  CHS
10  STO 02
11  1
12  STO 01
13  ST+ 02
14  ST+ 03
15  LASTX
16  CHS
17  RCL 04
18  X^2
19  CHS
20  XEQ "HGF+"   
21  STO 06
22  2
23  1/X
24  ST+ 02
25  ST+ 03
26  1
27  LASTX
28  CHS
29  RCL 00
30  XEQ "HGF+"   
31  RCL 04            
32  *
33  1
34  CHS
35  ACOS
36  RCL 05
37  *
38  STO 05
39  X<>Y
40  P-R
41  RCL 05            
42  RCL 06
43  P-R
44  X<> Z
45  -
46  X<> Z
47  +
48  END

 
     ( 75 bytes / SIZE 007 )
 
 

      STACK        INPUTS      OUTPUTS
           Y             n         En(x)
           X             x         Jn(x)

 
Examples:

   2   SQRT
  PI  XEQ "ANWEB"  >>>>    Jsqrt(2)(PI) =  0.366086558                                ---Execution time = 19s---
                                  X<>Y    Esqrt(2)(PI) = -0.315594386

   5   ENTER^
  PI      R/S     >>>>    J5(PI) =  0.052141184                                ---Execution time = 19s---
                     X<>Y    E5(PI) =  0.207000292
 

2°)  Asymptotic Expansions ( large x-values )
 

Formulae:

    Jn(x)  ~   Jn(x) + [ ( sin n.PI ) / ( x.PI ) ] ( F - n G / x )

    En(x)  ~  -Yn(x) - [ ( 1 + cos (n.PI) ) F - n ( 1 - cos n.PI ) G / x ] / ( x.PI )

  where   F = 1 + SUM k=1 to infinity  (n2-12) ....... (n2-(2k-1)2 / x2k
    and     G = 1 + SUM k=1 to infinity  (n2-22) ....... (n2-(2k)2 / x2k
 

Data Registers:   R00 = x ,  R01 = n

                             R02 = Jn(x) , R03 = Yn(x)      R04 to R06: temp
Flags: /
Subroutine:   "AEJY"  Asymptotic series for Bessel Functions ( cf "Bessel Functions for the HP-41" §6°) )
 
 
 

 01  LBL "AEANW"
 02  XEQ "AEJY"
 03  STO 02
 04  X<>Y
 05  STO 03
 06  RCL 01
 07  X^2
 08  1
 09  STO 05
 10  -
 11  RCL 00
 12  X^2
 13  /
 14  STO 06
 15  RCL 05
 16  +
 17  XEQ 01           
 18  STO 04
 19  CLX
 20  STO 05
 21  SIGN
 22  STO 06
 23  XEQ 01
 24  RCL 01
 25  *
 26  RCL 00
 27  /
 28  RCL 04
 29  X<>Y
 30  ST+ 04
 31  -
 32  GTO 02           
 33  LBL 01
 34  RCL 01
 35  X^2
 36  RCL 05
 37  2
 38  +
 39  STO 05
 40  X^2
 41  -
 42  RCL 06
 43  *
 44  RCL 00
 45  X^2
 46  /
 47  STO 06           
 48  +
 49  X#Y?
 50  GTO 01
 51  RTN
 52  LBL 02
 53  RCL 01
 54  PI
 55  R-D
 56  *
 57  X<>Y
 58  RCL 00
 59  PI
 60  *
 61  ST/ 04
 62  /
 63  P-R
 64  RCL 04           
 65  +
 66  CHS
 67  RCL 03
 68  -
 69  X<>Y
 70  RCL 02
 71  +
 72  END

 
   ( 95 bytes / SIZE 007 )
 
 

      STACK        INPUTS      OUTPUTS
           Y             n         En(x)
           X             x         Jn(x)

 
Examples:

    PI    ENTER^
  24.4  XEQ "AEANW"  >>>>    JPI(24.4) =   0.157155848               ---Execution time = 26s---
                                       X<>Y    EPI(24.4) = -0.008969347

Notes:

-With n = PI , the asymptotic series diverges too early if x = 24.3 and we fall into an infinite loop.
-We can also express F & G in terms of hypergeometric functions and we have:

    F = 3F0( (1+n)/2 , (1-n)/2 , 1 ; -4/x2 )   and   G = 3F0( (2+n)/2 , (2-n)/2 , 1 ; -4/x2 )

-So, the M-Code function HGF+ may be used to get faster results:
 

Data Registers:   R00 = x ,  R10 = n

                             R04 = Jn(x) , R05 = Yn(x)
Flags: /
Subroutine:   "AEJY"  ( cf "Bessel Functions for the HP-41" §6°) )  - the 2nd variant that uses HGF+
 
 

 01  LBL "AEANW"
 02  X<>Y
 03  STO 10
 04  X<>Y
 05  XEQ "AEJY"
 06  STO 04
 07  X<>Y
 08  STO 05
 09  1
 10  STO 03
 11  RCL 10
 12  +
 13  2
 14  /
 15  STO 01            
 16  RCL 03
 17  RCL 10
 18  -
 19  2
 20  /
 21  STO 02
 22  3
 23  0
 24  LASTX
 25  RCL 00
 26  /
 27  X^2
 28  CHS
 29  HGF+
 30  STO 06            
 31  .5
 32  ST+ 01
 33  ST+ 02
 34  3
 35  0
 36  LASTX
 37  HGF+
 38  RCL 10
 39  *
 40  RCL 00
 41  /
 42  RCL 06
 43  X<>Y
 44  ST+ 06
 45  -
 46  RCL 10            
 47  PI
 48  R-D
 49  *
 50  X<>Y
 51  RCL 00
 52  PI
 53  *
 54  ST/ 06
 55  /
 56  P-R
 57  RCL 06
 58  +
 59  CHS
 60  RCL 05            
 61  -
 62  X<>Y
 63  RCL 04
 64  +
 65  END

 
     ( 89 bytes / SIZE 011 )
 
 

      STACK        INPUTS      OUTPUTS
           Y             n         En(x)
           X             x         Jn(x)

 
Examples:

    PI    ENTER^
  24.4  XEQ "AEANW"  >>>>    JPI(24.4) =   0.157155848               ---Execution time = 16s---
                                       X<>Y    EPI(24.4) = -0.008969347
 

References:

[1]   Abramowitz and Stegun - "Handbook of Mathematical Functions" -  Dover Publications -  ISBN  0-486-61272-4
[2]   http://mathworld.wolfram.com/AngerFunctions.html
[3]   http://mathworld.wolfram.com/WeberFunctions.html
[4]   http://dlmf.nist.gov/11.11