
x = 0;
L = 2*pi;
h = L/2;
for j=1:8
  y1(j) = d_mid(@f,x,h);
  H1(j) = h;
  h = h/2;
end

for j=1:8
  N = 2^j;
  L = 2*pi;
  h = L/N;
  x = 0:h:L-h;
  F = f(x);
  Fx = fdiff(F,L);
  y2(j) = Fx(1);
  H2(j) = h;
%   [k,amp] = find_spec(F);  
%   A(j) = amp(length(k));
%   figure(j)
%   plot(k,amp)
end
% f(x) = sin(4x) --> f'(x) = 4*cos(4*x) --> f'(0) = 4
err1 = abs(y1-4)
err2 = abs(y2-4)
% f(x) = exp(sin(x)) --> f'(x) = cos(x)*exp(sin(x)) --> f'(0) = 1
err1 = abs(y1-1);
err2 = abs(y2-1);

format compact
format short e
display('meshwidth')
H1
display('finite difference errors')
err1
display('ratio of subsequent errors')
err1(1:7)./err1(2:8)
display(' ')
display('spectral errors')
err2 


% to test u(x) = tanh((x-x1)/eps) - tanh((x-x2)/eps)
% where x1 = .3, x2 = .7, eps = .05 for u'(1/4)
x = 1/4;
L = 1;
N = 4*16;
h = L/N;
for j=1:8
  y1(j) = d_mid(@f,x,h);
  H1(j) = h;
  h = h/2;
end

for j=1:8
  N = 2^(j+5);
  h = L/N;
  x = 0:h:L-h;
  F = f(x);
  Fx = fdiff(F,L);
  y2(j) = Fx(1+N/4);
  H2(j) = h;
%   [k,amp] = find_spec(F);  
%   A(j) = amp(length(k));
%   figure(j)
%   plot(k,amp)
end
% f'(1/4) = (-20*tanh(-1)^2+20*tanh(-9)^2))
fd_true =  (-20*tanh(-1)^2+20*tanh(-9)^2);
err1 = abs(y1-fd_true);
err2 = abs(y2-fd_true);


format compact
format short e
display('meshwidth')
H1
display('finite difference errors')
err1
display('ratio of subsequent errors')
err1(1:7)./err1(2:8)
display(' ')
display('spectral errors')
err2 




