function [ e ] = midpoint_error( a,b,f,h )

syms x
fdp = diff(f,2);
c = (b-a)/2;
e = (((b-a) * abs(vpa(subs(fdp,x,c)))) / 24)*h^2;
end