function [z] = trap(f,a,b,n)
    s = f(a);
    c = a;
    dx = (b-a)/n
    for i=1:1:n-1
        d = a+i*dx
        s = s + 2*f(d)
    end
    s = s + f(b)
    z = s/2;
endfunction


function y = fun(x)
    y = x^2
endfunction

a = trap(fun, 0, -5 , 5)
disp(a)