/*
*   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).
*/

#include<stdio.h>
#include<math.h>
#define true 1
#define false 0

static void  horner(double,double*, double*);

int main()
{
    double  Y, Z;
    double  X0;

    // set X0 value;
    X0 = -2.0;

    horner(X0, &Y,&Z);
  
    // print out the value of the P(x) and P'(x)
    printf("\n P ( %.10e ) = %12.8f\n", X0, Y);
    printf(" P' ( %.10e ) = %12.8f\n", X0, Z);
}

static void horner(
        double   X0,
	double   *outy,
        double   *outz)
{
   double AA[51];
   double Y,Z;
   int N,MM,I,J;

   // Degree of polynomial
   N = 4;
   // coefficients of P(X) in ascending order
   AA[0] = -4.0;
   AA[1] = 3.0;
   AA[2] = -3.0;
   AA[3] = 0.0;
   AA[4] = 2.0;

   /* STEP 1 */ 
   /* compute b(n) and for p(x) */ 
   Y = AA[N];
   /* compute b(n-1) for q(x) = p'(x) */
   if (N == 0)  Z = 0;
   else Z = AA[N];
   MM = N-1;
   /* STEP 2 */ 
   for (I=1; I<=MM; I++) 
   {
      J = N-I;
      /* compute b(j) for p(x) */
      Y = Y * X0 + AA[J];
      /* compute b(j-1) for q(x) */
      Z = Z * X0 + Y;
   }
   /* STEP 3 */ 
   /* compute b(0) for p(x) */ 
   if (N != 0) 
      Y = Y * X0 + AA[0];

   /* STEP 4 */ 
   printf("Coefficients of polynominal P :\n");
   for (I=0;I<=N;I++)   
      printf("Exponent = %3d    Coefficient = %12.8f\n", I, AA[I]);

   /* STEP 5*/
   // use pointer to return the computed Y and Z so
   // that we are able to use these values outside the 
   // function horner().

   *outy = Y;
   *outz = Z;
}