alg032.txt

Numerical Anaysis 8th Edition Burden and Faires (Maple Source)

> restart;
> # NEWTON'S INTERPOLATORY DIVIDED-DIFFERENCE FORMULA ALGORITHM 3.2        
> #
> # To obtain the divided-difference coefficients of the                    
> # interpolatory polynomial P on the (n+1) distinct numbers x(0),          
> # x(1), ..., x(n) for the function f:                                    
> #
> # INPUT:   numbers x(0), x(1), ..., x(n); values f(x(0)), f(x(1)),          
> #          ..., f(x(n)) as the first column Q(0,0), Q(1,0), ...,          
> #    Q(N,0) of Q, or may be computed if function f is supplied.
> #                                                                           # OUTPUT:  the numbers Q(0,0), Q(1,1), ..., Q(N,N) where
> #          P(x) = Q(0,0) + Q(1,1)*(x-x(0)) + Q(2,2)*(x-x(0))*             
> #          (x-x(1)) +... + Q(N,N)*(x-x(0))*(x-x(1))*...*(x-x(N-1)).
> alg032 := proc() local OK, FLAG, N, I, X, Q, A, NAME, INP, F, OUP, J;
> printf(`Newtons form of the interpolation polynomial\n`);
> OK := FALSE;
> while OK = FALSE do
> printf(`Choice of input method:\n`);
> printf(`1. Input entry by entry from keyboard\n`);
> printf(`2. Input data from a text file\n`);
> printf(`3. Generate data using a function F\n`);
> printf(`Choose 1, 2, or 3 please\n`);
> FLAG := scanf(`%d`)[1];
> if FLAG = 1 or FLAG = 2 or FLAG = 3 then
> OK := TRUE;
> fi;
> od;
> if FLAG = 1 then
> OK := FALSE;
> while OK = FALSE do
> printf(`Input n\n`);
> N := scanf(`%d`)[1];
> if N > 0 then
> OK := TRUE;
> for I from 0 to N do
> printf(`Input X(%d) and F(X(%d)) `, I, I);
> printf(`separated by space\n`);
> X[I] := scanf(`%f`)[1];
> Q[I,0] := scanf(`%f`)[1];
> od;
> else 
> printf(`Number must be a positive integer\n`);
> fi;
> od;
> fi;
> if FLAG = 2 then
> printf(`Has a text file been created with the data in two columns ?\n`);
> printf(`Enter Y or N\n`);
> A := scanf(`\n%c`)[1];
> if A = "Y" or A = "y" then
> printf(`Input the file name in the form - `);
> printf(`drive:\\name.ext\n`);
> printf(`For example:   A:\\DATA.DTA\n`);
> NAME := scanf(`%s`)[1];
> INP := fopen(NAME,READ,TEXT);
> OK := FALSE;
> while OK = FALSE do
> printf(`Input n\n`); 
> N := scanf(`%d`)[1];
> if N > 0 then                  
> for I from 0 to N do
> X[I] := fscanf(INP, `%f`)[1];
> Q[I,0] := fscanf(INP, `%f`)[1];
> od;
> fclose(INP);
> OK := TRUE;
> else
> printf(`Number must be a positive integer\n`)
> fi;
> od;
> else
> printf(`Please create the input file in two column `);
> printf(`form with the X values and\n`);
> printf(`F(X) values in the corresponding columns.\n`);
> printf(`The program will end so the input file can `);
> printf(`be created.\n`);
> OK := FALSE;
> fi;
> fi;
> if FLAG = 3 then
> printf(`Input the function F(x) in terms of x\n`);
> printf(`For example: cos(x)\n`);
> F := scanf(`%a`)[1];
> F := unapply(F,x);
> OK := FALSE;
> while OK = FALSE do
> printf(`Input n\n`);
> N := scanf(`%d`)[1];
> if N > 0 then
> for I from 0 to N do
> printf(`Input X(%d)\n`, I);
> X[I] := scanf(`%f`)[1];
> Q[I,0] := F(X[I]);
> od;
> OK := TRUE;
> else
> printf(`Number must be a positive integer\n`);
> fi;
> od;
> fi;
> if OK = TRUE then
> printf(`Select output destination\n`);
> printf(`1. Screen\n`);
> printf(`2. Text file\n`);
> printf(`Enter 1 or 2\n`);
> FLAG := scanf(`%d`)[1];
> if FLAG = 2 then
> printf(`Input the file name in the form - drive:\\name.ext\n`);
> printf(`For example:   A:\\OUTPUT.DTA\n`);
> NAME := scanf(`%s`)[1];
> OUP := fopen(NAME,WRITE,TEXT);
> else
> OUP := default;
> fi;
> fprintf(OUP, `NEWTONS INTERPOLATION POLYNOMIAL\n\n`);
> # STEP 1
> for I from 1 to N do
> for J from 1 to I do
> Q[I,J] := (Q[I,J-1] - Q[I-1,J-1]) / (X[I] - X[I-J]);
> od;
> od;
> # STEP 2
> fprintf(OUP, `Input data follows:\n`);      
> for I from 0 to N do
> fprintf(OUP, `X(%d) = %12.8f F(X(%d)) = %12.8f\n`, I, X[I], I, Q[I,0]);      
> od;
> fprintf(OUP, `\nThe coefficients Q(0,0), ..., Q(N,N) are:\n`);      
> for I from 0 to N do
> fprintf(OUP, `%12.8f\n`, Q[I,I]);
> od;
> if OUP <> default then
> fclose(OUP);
> printf(`Output file %s created successfully`,NAME);
> fi;
> fi;
> RETURN(0);
> end;
> alg032();