#include <stdio.h> /* for input/output */
#include <math.h> /* for maths */

float f(float x) { 
  return 5.0*pow(x,4.0) + 4.0*pow(x,3.0) + 3.0*pow(x,2.0) + 2.0*x + 1.0;
}

main() {
  float x1, xN; /* first and last terms */
  float I; /* integral result */
  float h; /* trapezium bar width */
  float error; /* maximum error to accept */
  int i, lastN; /* term number and last term to do */
  printf("Enhanced Extended Trapezoidal Numeric Integrator\nBy Ian Hickson\nRecompile code to change equation.\n");
  printf("\nEnter the low bound: ");
  scanf("%f", &x1);
  printf("\nFYI, the equation at the low bound is: f(%f)=%f.", x1, f(x1));
  printf("\nEnter the high bound: ");
  scanf("%f", &xN);
  printf("\nFYI, the equation at the high bound is: f(%f)=%f.", xN, f(xN));
  printf("\nEnter the maximum error: ");
  scanf("%d", &error);
  printf("\nWorking...");
  I = 0.5 * (f(x1)+f(xN)); /* First and last terms */
  lastN = trunc(sqrt(1/error)+1);
  h = (xN - x1) / (lastN - 1); /* Find the trapezium width */
  for (i=1; i<lastN; i++) {
    I += f(x1 + i*h);
  }
  I *= h;
  printf ("\nThe integral is approximately %f, the error is of order %f.\n", I, 1/(pow(i,2)));
}



/* OUTPUT:
none yet
 */

/* COMPILATION ERRORS 
none yet
 */