alg027.txt

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

> restart;
> #  HORNER'S ALGORITHM 2.7
> #
> #  To evaluate the polynomial
> #  p(x) = a(n) * x^n + a(n-1) * x^(n-1) + ... + a(1) * x + a(0)
> #  and its derivative p'(x) at x = x0;
> #
> #  INPUT:   degree n; coefficients aa(0),aa(1),...,aa(n);
> #           value of x0.
> #
> #  OUTPUT:  y = p(x0), z = p'(x0).
> alg027 := proc() local OK, N, I, AA, X0, Y, Z, MM, J;
> printf(`This is Horners Method\n`);
> OK := FALSE;
> while OK = FALSE do
> printf(`Input degree n of polynomial - no decimal point\n`);
> N := scanf(`%d`)[1];
> if N < 0 then
> printf(`Incorrect input - degree must be nonnegative.\n`);
> else
> OK := TRUE;
> fi;
> od;
> printf(`Input coefficients of P(X) in ascending order\n`);
> for I from 0 to N do
> printf(`Input coefficient of X^%d\n`, I);
> AA[I] := scanf(`%f`)[1];
> od;
> printf(`Input argument X0 at which to evaluate P(X)\n`);
> X0 := scanf(`%f`)[1];
> # Step 1
> # Compute b(n) for p(x)
> Y := AA[N];
> # Compute b(n-1) for p'(x)
> if N = 0 then
> Z := 0;
> else
> Z := AA[N];
> fi;
> MM := N-1;
> # Step 2
> for I from 1 to MM do
> J := N-I;
> # Compute b(j) for p(x)
> Y := Y*X0+AA[J];
> # Compute b(j-1) for p'(x)
> Z := Z*X0+Y;
> od;
> # Step 3
> # Compute b(0) for p(x)
> if N <> 0 then
> Y := Y*X0+AA[0];
> fi;
> printf(`Coefficients of polynominal P :\n`);
> # Step 4
> for I from 0 to N do
> printf(`Exponent = %3d Coefficient = %12.8f\n`, I, AA[I]);
> od;
> printf(`\n P ( %.10e ) = %12.8f\n`, X0, Y);
> printf(` P' ( %.10e ) = %12.8f\n`, X0, Z);
> RETURN(0);
> end;
> alg027();