(71B) Gauss-Chebyshev Quadrature

[Namir]

This is an HP-71B program that performs the Gauss-Chebyshev Quadrature for a function f(x). This algorithm is of teh same class as the Gauss-Legendre Qudrature, except calculating the roots of the Chebyshev polynomials is MUCH easier than calculating the roots of the Legendre polynomial. Thus, you can easily specify the number of points (which is also the order of the Chebyshev polynomial) used to calculate the integral.

The program prompts you to enter A, B which define the integral limit and also prompts you for N the Chebyshev polynomial order (which is also the number of points used in the quadrature. The program displays the value of the integral, stored in variable R. The listing is:

Code:

10 REM GAUSS-CHEBYSHEV QUADRATURE
20 INPUT "A? ";A
30 INPUT "B? ";B
40 INPUT "N? ";N
50 T1=(B-A)/2
60 T2=(B+A)/2
70 T3=PI/2/N
80 S=0
90 RADIANS
100 FOR I=1 TO N
110 Y=COS((2*I-1)*T3)
120 X=T1*Y+T2
130 GOSUB 1000
140 S=S+F*SQR(1-Y*Y)
150 NEXT I
160 R=2*T1*T3*S
170 DISP "AREA=";R @ PAUSE
180 GOSUB 2000
190 END
1000 REM F(X)
1010 F=1/X
1020 RETURN
2000 REM COMPARE WITH ACTUAL INTEGRAL
2010 Y=LOG(B/A)
2020 DISP "EXACT AREA=";Y @ PAUSE
2030 DISP "% ERR=";100*(Y-R)/Y
2040 RETURN

Line 1010 starts the calculation for the function f(x) which is stored in the variable F. Line 2010 is optional and allows you to calculate or assign the exact value of the integral and store it in the variable Y. The subroutine displays that value and also displays the % error. If calculating the exact integral is too difficult or not possible, then make line 2010 to simply contain a RETURN statement.

Enjoy!

Namir