#include <math.h>
#include <stdio.h>

long fact(int n) {
    if (n == 0) {
        return 1;
    } else {
        return n * fact(n-1);
    }
}

double Differentiate(double f(double x), double x, double h, int grade) {
    double sum = 0;
    for (int i = 0; i <= grade; i++) {
        printf("%g %d %g %g\n", h, grade, pow(h, 4), pow(h, grade));
	    // printf("%d %f %f %f %g\n", i, sum, fact(grade) / (fact(i) * fact(grade - i)), f(x + (double) i * h), pow(h, grade));
        sum += pow(-1, grade - i) * (fact(grade) / (fact(i) * fact(grade - i))) * f(x + (double) i * h) / pow(h, grade);
    }
    return sum;
}

double f(double x) {
    return x*x*x*x*x - x/(2-x*x);
}

double SimpsonsRule(double f(double x), double a, double b, int errGrade) {
    int res = 0;
    double der4;
    double c, cmax;
    for (c = a; c <= b; c += 0.1) {
        der4 = Differentiate(f, c, sqrt(__DBL_EPSILON__), 4);
        if (res < der4) {
            res = der4;
            cmax = c;
        }
    }

    printf("%f\n", cmax);
    double limit = pow(10, -errGrade), n = 1;
    double h = fabs(b-a) / n;
    double err = (b-a) * pow(h, 4) * res / 180;
    while (err > limit) {
        printf("%g %g\n", err, limit);
        n++;
        h = fabs(b-a) / n;
        err = (b-a) * pow(h, 4) * res / (double) 180;
    }
    printf("%d", n);

// function s=simpson(f,a,b,n)
//   h=(b-a)/n;
//   y=[];for x=a:h:b;y=[y f(x)];end
//   E=y(1)+y(end);
//   P=sum(y(3:2:end-2)); Im=sum(y(2:2:end-1));
//   s=(h/3)*(E+2*P+4*Im);
// end

    double x, sum = 0;
    int i = 0; // de 2 a n-1
    while (i < n) {
        x = a + i * h;
        if (i % 2 == 0) {
            sum += 2 * f(x);
        } else {
            sum += 4 * f(x);
        }
        i++;
    }
    return (h/3) * (f(a) + f(b) + sum);
}

double g(double x) {
    return 1/(x*x*x) + x*x;
}

int main(void) {
    // printf("%f", )
    // printf("der %f\n", Differentiate(g, 1, pow(10, -1), 4));
    // printf("der %f\n", );
    // printf("culo %g", SimpsonsRule(g, 1, 8, 5));
}