Series for the HP-41

Overview
 

 1°)  Simple Series

   a)  U_n = f(n)
   b)  U_n+1 = f(U_n)
   c)  Euler Transformation

 2°)  Double Series
 3°)  Triple Series
 4°)  Multiple Series
 

Warning:

-In these programs, the iterations stop when 2 consecutive partial sums are equal.
-Therefore, wrong results may be obtained if one term is negligible and the next one is not !
 

1°)   Simple Series
 

          a)  U_n = f(n)
 

-The following program calculates  S = u_k + u_k+1 +  u_k+2 + u_k+3 + .........    =  Sigma  u_n    for  n >= k
  defined by    u_n = f(n)  where f is a known function.
 

Data Registers:           •  R00 =  f name              ( Register R00 is to be initialized before executing "SIGMA1" )

                                          R01 = S
                                          R02 = n
Flags: /
Subroutine:  a program which computes  u_n = f(n)  assuming n is in X-register upon entry
 
 

01  LBL "SIGMA1"  02  STO 02  03  CLX  04  STO 01  05  LBL 01  06  RCL 02  07  XEQ IND 00  08  ISG 02  09  LBL 02  10  RCL 01  11  +  12  STO 01  13  LASTX  14  X#Y?  15  GTO 01  16  END

 
( 30 bytes / SIZE 003 )
 
 

STACK	INPUTS	OUTPUTS
Y	/	S
X	k	S

 
 
Example:   Calculate  S = 1 + 1/4 + 1/27 + ....... + 1/nⁿ + .........

-Here we have  f(n) = 1/nⁿ  and  k = 1 and this function may be obtained by

  LBL "T"
  ENTER^
  CHS
  Y^X
  RTN

-Then,

    "T"  ASTO 00
     1    XEQ "SIGMA1"  yields   X = 1.291285997 = R01   ( in 9 seconds )
 

          b) U_n+1 = f(U_n)
 

-The following program calculates  S = u_k + u_k+1 +  u_k+2 + u_k+3 + .........    =  Sigma  u_n    for  n >= k
  when  u_k is given  and   u_n+1 =  f(u_n;n)   where f is a known function.
 

Data Registers:           •  R00 =  f name              ( Register R00 is to be initialized before executing "SIGMAR" )

                                          R01 = S
                                          R02 = n
                                          R03 = u_n
Flags: /
Subroutine:  a program which computes  u_n+1 = f(u_n;n)  assuming u_n is in X-register and n is in Y-register upon entry
 
 

01  LBL "SIGMAR"  02  STO 01  03  STO 03  04  X<>Y  05  STO 02  06  LBL 01  07  RCL 02  08  RCL 03  09  XEQ IND 00  10  STO 03  11  ISG 02  12  LBL 02  13  RCL 01  14  +  15  STO 01  16  LASTX  17  X#Y?  18  GTO 01  19  END

 
( 33 bytes / SIZE 004 )
 
 

STACK	INPUTS	OUTPUTS
Y	k	S
X	uk	S

 
Example:    Calculate  S = 1 + 1/2! + 1/3! + ....... + 1/n! + .........

-It may be obtained with the first program but we can also define  u₁ = 1 and u_n+1 = u_n/(n+1)

  LBL "T"
  X<>Y
  1
  +
  /
  RTN

  "T"  ASTO 00
   1    ENTER^   ( k = 1 = u₁ )
   XEQ "SIGMAR"    produces    X = 1.718281830 = R01   in 7 seconds  ( the exact result is e-1 )
 

          c) The Euler Transformation
 

-The Euler transformation is an acceleration method developed for alternating series:    S = u₀ - u₁ +  u₂ - u₃ + ...... + (-1)ⁿ  u_n + .......

-The sum is rewritten:   S = a₀/2 + ( C₁¹ a₀ - C₁⁰ a₁ ) / 2² + ( C₂² a₀ - C₂¹ a₁ + C₂⁰ a₂  ) / 2³ +  ( C₃³ a₀ - C₃² a₁ + C₃¹ a₂ - C₃⁰ a₃  ) / 2⁴ + ................

     where   C_n^p =  n! / ( p! ( n-p )! )   are the binomial coefficients.

-This may produce superb acceleration but it can also fail...
 

Data Registers:           •  R00 =  f name              ( Register R00 is to be initialized before executing "EULER" )

                                          R01 = S
                                          R02 thru R05: temp
Flags: /
Subroutine:  a program which computes  u_n = f(n)  assuming n is in X-register upon entry
 
 

01  LBL "EULER" 02  CLX 03  XEQ IND 00 04  2 05  / 06  STO 01 07  1 08  STO 02 09  LBL 01 10  RCL 02 11  STO 03 12  SIGN

13  STO 04 14  CLX 15  STO 05 16  LBL 02 17  RCL 02 18  RCL 03 19  - 20  XEQ IND 00  21  RCL 04 22  * 23  ST+ 05 24  LASTX

25  CHS 26  RCL 03 27  * 28  RCL 02 29  LASTX           30  - 31  1 32  + 33  / 34  STO 04 35  DSE 03 36  GTO 02

37  RCL 02 38  XEQ IND 00  39  RCL 04 40  * 41  RCL 05 42  + 43  2 44  RCL 02 45  1 46  + 47  STO 02 48  Y^X

49  / 50  RCL 01 51  + 52  STO 01           53  LASTX 54  VIEW X 55  X#Y? 56  GTO 01 57  END

 
( 75 bytes / SIZE 006 )
 
 

STACK	INPUTS	OUTPUTS
X	/	S

 
  Example:   Compute  S = 1 - 2^-1/2 +  3^-1/2 - 4^-1/2 + .......  + (-1)ⁿ (n+1)^-1/2 + ......

-Here, we have   f(n) = (n+1)^-1/2

  LBL "T"
  SIGN
  LASTX
  +
  SQRT
  1/X
  RTN

 "T"  ASTO 00
  XEQ "EULER"  the successive approximations are displayed and  6mn27s later,   X = 0.6048986431 = R01

-The error is  -3 10^-10 (!)
-Actually, only 28 terms are calculated
-Without an acceleration method, more than 10¹⁸ terms would be necessary to achieve the same accuracy...
  ( execution time would be much greater than the age of the Universe! )

Notes:

-This program computes the different f(n) several times and this may be time-consuming if f is a complicated function.
-An alternative is to store the successive terms into data registers:

   Add   ST+ 06           after line 45
   Add   STO IND 06  after line 38
   Replace line 20  by   7  +  RCL IND X
   Add   8  STO 06      after line 08
   Add   STO 07          after line 03

-The same result is now obtained in 5mn19s ( but the SIZE must be at least 035 in this example )

-It is often better to sum the first terms and apply the Euler transformation to the rest:
-In the example above, if you sum  u₀ - u₁ +  u₂ - u₃ + ...... -  u₉ =  1 - 2^-1/2 +  3^-1/2 - 4^-1/2 + .......  - 10^-1/2
  and execute "EULER" with

  LBL "T"
  11
  +
  SQRT
  1/X
  RTN

-Execution time becomes:  1mn22s
 

2°)   Double Series
 

-The following program calculates  S = Sigma  u_n;m    for  n >= n₀ ; m >= m₀
 

Data Registers:           •  R00 =  f name              ( Register R00 is to be initialized before executing "SIGMA2" )

                                          R01 = n₀    R03 = S    R05 = m
                                          R02 = m₀   R04 = n

Flags: /
Subroutine:  a program which computes  u_n;m = f(n;m)  assuming n is in X-register and m is in Y-register upon entry
 
 

01  LBL "SIGMA2" 02  STO 01 03  STO 04 04  X<>Y 05  STO 02 06  STO 05 07  CLX

08  STO 03 09  DSE 05             10  LBL 01 11  ISG 05 12  CLX 13  RCL 05 14  RCL 04

15  XEQ IND 00    16  RCL 03 17  + 18  STO 03 19  LASTX 20  X#Y? 21  GTO 01

22  ISG 04 23  CLX 24  RCL 02 25  STO 05 26  RCL 04 27  XEQ IND 00    28  RCL 03

29  + 30  STO 03 31  LASTX             32  X#Y? 33  GTO 01 34  END

 
   ( 52 bytes / SIZE 006 )
 
 

STACK	INPUTS	OUTPUTS
Y	m0	S
X	n0	S

 
Example:    Calculate  S = Sigma  1 / ( nⁿ m! )      with  n >= 1 and m >= 1

-Here, we have   f(n;m) = 1 / ( nⁿ m! )  and this function may be obtained by

  LBL "T"
  ENTER^
  Y^X
  X<>Y
  FACT
  *
  1/X
  RTN

-Then

   "T"  ASTO 00
    1    ENTER^
    XEQ "SIGMA2"    and   1mn18s later     X = 2.218793264 = R03
 

Note:

-If the series converges too slowly, execution time becomes prohibitive.
-Lines 16 to 19 and 28 to 31 may then be replaced by    X<> 03   ST+ 03  RND  RCL 03  RND
  and the accuracy will be determined by the display format.
-Similar modifications apply to the other programs listed in this page.
 

3°) Triple Series
 

-The following program calculates  S = Sigma  u_n;m;p    for  n >= n₀ ; m >= m₀ ; p >= p₀
 

Data Registers:           •  R00 =  f name              ( Register R00 is to be initialized before executing "SIGMA3" )

                                          R01 = n₀        R04 = S       R07 = p
                                          R02 = m₀       R05 = n
                                          R03 = p₀        R06 = m

Flags: /
Subroutine:  a program which computes  u_n;m;p = f(n;m;p)  assuming n is in X-register ; m is in Y-register and p is in Z-register upon entry.
 
 

01  LBL "SIGMA3" 02  STO 01 03  STO 05 04  RDN 05  STO 02 06  STO 06 07  X<>Y 08  STO 03 09  STO 07 10  CLX 11  STO 04

12  DSE 07 13  LBL 01 14  ISG 07 15  CLX 16  RCL 07 17  RCL 06 18  RCL 05 19  XEQ IND 00   20  RCL 04 21  + 22  STO 04

23  LASTX 24  X#Y? 25  GTO 01 26  ISG 06 27  CLX 28  RCL 03 29  STO 07 30  RCL 06 31  RCL 05 32  XEQ IND 00   33  RCL 04

34  + 35  STO 04            36  LASTX 37  X#Y? 38  GTO 01 39  ISG 05 40  CLX 41  RCL 03 42  STO 07 43  RCL 02 44  STO 06

45  RCL 05 46  XEQ IND 00   47  RCL 04 48  + 49  STO 04 50  LASTX 51  X#Y? 52  GTO 01 53  END

 
    ( 74 bytes / SIZE 008 )
 
 

STACK	INPUTS	OUTPUTS
Z	p0	/
Y	m0	S
X	n0	S

 
Example:    Calculate  S = Sigma  1 / ( nⁿ m! (p!)² )      with  n >= 1 ; m >= 1 ; p >= 1

-Here, we have   f(n;m;p) = 1 / ( nⁿ m! (p!)² )  and this function may be obtained by

  LBL "T"
  ENTER^
  Y^X
  X<>Y
  FACT
  *
  X<>Y
  FACT
  X^2
  *
  1/X
  RTN

-Then,

  "T"  ASTO 00
    1   ENTER^  ENTER^
   XEQ "SIGMA3"   and  6mn53s  later,   X = 2.839135243 = R04   ( error = -7 10^-9 )
 

4°) Multiple Series
 

-The following program calculates  S = Sigma  u(n₁;n₂; ..... ;n_k)   with  n₁>= 0 ;  n₂>= 0 ; ........ ;  n_k>= 0
 

Data Registers:           •  R00 =  f name              ( Register R00 is to be initialized before executing "SIGMA" )

                                          R01 = n₁   R02 = n₂  .........  Rkk = n_k
Flags: /
Subroutine:  a program which computes  u(n₁;n₂; ..... ;n_k)  = f( n₁;n₂; ..... ;n_k)  assuming n₁is in R01 ; n₂ is in R02 ; ..... ; n_k is in Rkk  upon entry.

-Synthetic registers  M  N  O  may be replaced by any unused data registers.
-If you don't have an HP-41CX , replace line 08 with   0   LBL 00   STO IND Y   ISG Y   GTO 00
 
 

01  LBL "SIGMA" 02  CLA 03  STO N 04  STO O 05   E3 06  / 07  ISG X 08  CLRGX  09  DSE IND N

10  LBL 01 11  ISG IND N 12  CLX 13  XEQ IND 00  14  RCL M 15  + 16  STO M 17  LASTX 18  X#Y?

19  GTO 01 20  LBL 02 21  CLX 22  STO IND N   23  DSE N 24  X<0? 25  GTO 03 26  ISG IND N 27  CLX

28  XEQ IND 00  29  RCL M 30  + 31  STO M 32  LASTX 33  X=Y? 34  GTO 02 35  RCL O 36  STO N

37  GTO 01         38  LBL 03 39  X<>Y 40  CLA 41  END

 
   ( 73 bytes / SIZE kkk+1 )
 
 

STACK	INPUTS	OUTPUTS
X	k	S

   With  k = 1 for simple series ; 2 for double series ; 3 for triple series ; 4 for quadruple series ... etc ...

Example:   Let's take once again the triple series   S = Sigma  1 / ( nⁿ m! (p!)² )      for  n >= 1 ; m >= 1 ; p >= 1

-We make a change of arguments to reduce to the standard:   n >= 0 ; m >= 0 ; p >= 0   by replacing  n with (n+1) ; m with (m+1) ; p with (p+1)

  LBL "T"
  RCL 01
  1
  +
  ENTER^
  Y^X
  RCL 02
  1
  +
  FACT
  *
  RCL 03
  1
  +
  FACT
  X^2
  *
  1/X
  RTN

-Then

  3  XEQ "SIGMA"  yields   X = 2.839135243    ( in 8mn39s )