% Sample file for left-endpoint singular integral using Simpson's

%Input data
% clear all;

a = 0;
b = 1;
q=1;
s = ones(1,3);
nmin = 128; % needs to be even
for k=1:3
  n = nmin*2^k; 
  h = (b-a)/n;
  x = a:h:b;
  xr = a+h:h:b;
  fr = xr.^(3*q-1).*sin(1./(xr.^q));
  f = [0 fr];

%weight=(1,4,2,4,2,...,4,2,4,1)
  w = ones(1,n+1);
  w(2:2:n)=4;
  w(3:2:n-1)=2;
  s(k) = q*(h/3)*w*f'; 
end
%Sum
log2( (s(1)-s(2))/(s(2)-s(3)) )

%s =  3.10933712650887