import <Utilities/Conversion.sl>;
import <Utilities/Sequence.sl>;

integrateLeft(f, a, b, n) :=
    let
        h := (b - a) / n;
        vals[x] := f(x) foreach x within (0 ... (n-1)) * h + a;
    in
        h * sum(vals);

integrateRight(f, a, b, n) :=
    let
        h := (b - a) / n;
        vals[x] := f(x+h) foreach x within (0 ... (n-1)) * h + a;
    in
        h * sum(vals);

integrateMidpoint(f, a, b, n) :=
    let
        h := (b - a) / n;
        vals[x] := f(x+h/2.0) foreach x within (0 ... (n-1)) * h + a;
    in
        h * sum(vals);

integrateTrapezium(f, a, b, n) :=
    let
        h := (b - a) / n;
        vals[i] := 2.0 * f(a + i * h) foreach i within 1 ... n-1;
    in
        h * (sum(vals) + f(a) + f(b)) / 2.0;

integrateSimpsons(f, a, b, n) :=
    let
        h := (b - a) / n;
        vals1[i] := f(a + h * i + h / 2.0) foreach i within 0 ... n-1;
        vals2[i] :=  f(a + h * i) foreach i within 1 ... n-1;
    in
        h / 6.0 * (f(a) + f(b) + 4.0 * sum(vals1) + 2.0 * sum(vals2));

xCubed(x) := x^3;
xInverse(x) := 1/x;
identity(x) := x;

tests[method] :=
    [method(xCubed, 0.0, 1.0, 100),
     method(xInverse, 1.0, 100.0, 1000),
     method(identity, 0.0, 5000.0, 5000000),
     method(identity, 0.0, 6000.0, 6000000)]
     foreach method within [integrateLeft, integrateRight, integrateMidpoint, integrateTrapezium, integrateSimpsons];

//String manipulation for ouput display.
main :=
    let
        heading := [["Func", "Range\t", "L-Rect\t", "R-Rect\t", "M-Rect\t", "Trapezium", "Simpson"]];
        ranges := [["0 - 1\t", "1 - 100\t", "0 - 5000", "0 - 6000"]];
        funcs := [["x^3", "1/x", "x", "x"]];
    in
        delimit(delimit(heading ++ transpose(funcs ++ ranges ++ trimEndZeroes(floatToString(tests, 8))), '\t'), '\n');

trimEndZeroes(x(1)) := x when size(x) = 0 else x when x[size(x)] /= '0' else trimEndZeroes(x[1...size(x)-1]);