# hp41programs

Clausen's Integral

# Clausen's Integral for the HP-41

Overview

-The following program computes   f(x) = - §0x Ln ( 2 sin (t/2) ) dt                     ( § = Integral )

-Integration by parts gives:   f(x) = -x Ln ( 2 sin (x/2) ) + 2 §0x/2 ( t / Tan t ) dt
-The series expansion   t / Tan t = 1 - x2/3 - x4/45 - 2 x6/945 - ..... - [ (-1)n-1 22n/(2n)! ] B2n x2n - .....................  is used  where  B2n are Bernoulli numbers.
-Bernoulli numbers B2n are calculated through  Zeta function for n > 3

Program Listing

Data Registers:    R00-R01  are used by "ZETA"

R02 = x    R03 thru R06: temp   (  R05 = §0x/2 ( t / Tan t ) dt )
Flags: /
Subroutine:   "ZETA"  ( cf "Riemann Zeta Function for the HP-41" )

-If x is small enough ( for example  x = 0.15 ), "ZETA" is not needed.

 01  LBL "CLAUS" 02  DEG 03  STO 02 04  2 05  / 06  ENTER^ 07  X^2 08  STO 04           09  LASTX 10  * 11  STO 03 12  9 13  / 14  - 15  RCL 03 16  RCL 04 17  * 18  STO 03  19  225 20  / 21  - 22  RCL 03 23  RCL 04  24  * 25  STO 03           26  3307.5 27  / 28  CHS 29  X<>Y 30  + 31  STO 05 32  LASTX 33  X=Y?  34  GTO 00 35  2 36  PI 37  6 38  STO 06  39  Y^X 40  / 41  CHS 42  ST* 03 43  FIX 9 44  LBL 01 45  RCL 06           46  2 47  + 48  STO 06 49  XEQ "ZETA"  50  RCL 03 51  RCL 04  52  * 53  PI 54  X^2 55  / 56  STO 03           57  * 58  RCL 06 59  1 60  + 61  / 62  RCL 05 63  + 64  STO 05 65  LASTX 66  X#Y? 67  GTO 01 68  LBL 00 69  ST+ X 70  RCL 02           71  2 72  / 73  R-D  74  SIN 75  ST+ X 76  LN 77  RCL 02 78  * 79  - 80  END

( 108 bytes / SIZE 007 )

 STACK INPUTS OUTPUTS Y / g(x) X x f(x)

where  f(x) = Clausen's Integral            ( 0 < x < 2 PI )
and    g(x) = §0x/2 ( t / Tan t ) dt

Examples:

PI  3  /   XEQ "CLAUS"  >>>>   f(PI/3) = 1.014941606                ---Execution time = 28s---
X<>Y   g(PI/3) = 0.507470803

0.15   R/S    >>>>   f(0.15) = 0.434614878       ( in 2 seconds )
X<>Y   g(0.15) = 0.074953114

6      R/S    >>>>   f(6) = -0.640782664           ( in 7mn49s )
X<>Y   g(6) = -4.115383667

Notes:

-Execution time tends to infinity as x tends to 2 PI
-This routine produces DATA ERROR if x = 0 but since f(0) = 0, simply add  X=0?  RTN  after line 01

-Clausen's integral is also equal to  SUMk=1,2,......  [ Sin (kx) ] / k2 but thousands of terms would be required to return accurate results if we used this formula.
-However, Clausen's integral Cl(x) is closely related to Lobachevsky function L(x), we have:

Cl(x) = 2 L(x/2)     at least if  x >= 0 & x <= 2.PI

-So, we can use the program "LOB" listed in "Lobachevsky Function for the HP-41"

Reference:

[1]   Abramowitz and Stegun - "Handbook of Mathematical Functions" -  Dover Publications -  ISBN  0-486-61272-4