home *** CD-ROM | disk | FTP | other *** search
- PROGRAM d4r8(input,output);
- (* driver for routine GAULEG *)
- CONST
- npoint=10;
- x1=0.0;
- x2=1.0;
- x3=10.0;
- nval=10;
- TYPE
- double = real;
- darray=ARRAY [1..npoint] OF double;
- VAR
- i,j : integer;
- xx : real;
- x,w : darray;
-
- FUNCTION func(x: real): real;
- BEGIN
- func := x*exp(-x)
- END;
-
- FUNCTION sngl(x : real): real;
- BEGIN
- sngl := x
- END;
-
- (*$I GAULEG.PAS *)
-
- BEGIN
- gauleg(x1,x2,x,w,npoint);
- writeln;
- writeln('#':2,'x[i]':11,'w[i]':12);
- FOR i := 1 to npoint DO BEGIN
- writeln(i:2,x[i]:12:6,w[i]:12:6)
- END;
- (* demonstrate the use of gauleg for an integral *)
- gauleg(x1,x3,x,w,npoint);
- xx := 0.0;
- FOR i := 1 to npoint DO BEGIN
- xx := xx+sngl(w[i])*func(sngl(x[i]))
- END;
- writeln;
- writeln('Integral from GAULEG: ',xx:12:6);
- writeln('Actual value: ',((1.0+x1)*exp(-x1)-(1.0+x3)*exp(-x3)):12:6)
- END.
-