% INTSIN0 Integration relative to the weight function w(t)=sin(1/t)
%
f1='%5.0f %21.13e %21.13e %21.13e\n';
f2='%5.0f %21.13e\n';
fprintf('\n')
disp('    n         1+sin0               Legendre               sin0')
load -ascii absin0;
ab=absin0; abl=r_jacobi01(40);
%delta=.5;
%delta=.3;
delta=.1;
%for n=2:2:16
%for n=3:3:21
for n=4:4:36
  xwl=gauss(n,abl); xw=gauss(n,ab);
  intl=sum(xwl(:,2).*tan((pi/2-delta).*xwl(:,1)));
  ints=sum(xw(:,2).*tan((pi/2-delta).*xw(:,1)));
  int=ints-intl;
  fprintf(f1,n,ints,intl,int)
end
fprintf('\n')
abL=r_jacobi01(1000);
for N=200:200:1000
  xwL=gauss(N,abL);
  intL=sum(xwL(:,2).*tan((pi/2-.1)*xwL(:,1)).*sin(1./xwL(:,1)));
  fprintf(f2,N,intL)
end