% EX 1.3 c)

% Beräkna derivata y'(7) mha differensapprox
% a) y'(x) = 1/h (y(x + h) - y(x))
% b) y'(x) = 1/(2h) (-y(x + 2h) + 4y(x + h) - 3y(x)) 
% för h = 0.04, 0.02 0.01 och 0.005

% Plotta funktionen
x = 7:0.005:8;
y = ((x - 7).^(5/2) + 2*sin(pi*sqrt(x)))./...
  (sqrt(x + 4*log(x - 2*pi)) - 1);
plot(x, y)

% Sant y'(7) och y(7)
y_prim_7 = -1.68043;
y_7 = 2*sin(pi*sqrt(7))/(sqrt(7 + 4*log(7 - 2*pi)) - 1); 


% Tabell
h = [0.04, 0.02, 0.01, 0.005];
disp('       |    approx a)         approx b)')
disp('   h   |  värde    fel      värde     fel')
disp('-------|------------------------------------')
for i = 1:length(h)
  y_7_plus_h = (h(i)^(5/2) + 2*sin(pi*sqrt(7 + h(i))))/...
    (sqrt(7 + h(i) + 4*log(7 + h(i) - 2*pi)) - 1);
  y_7_plus_2h = ((2*h(i))^(5/2) + 2*sin(pi*sqrt(7 + 2*h(i))))/...
    (sqrt(7 + 2*h(i) + 4*log(7 + 2*h(i) - 2*pi)) - 1);

  y_prim_7_a = 1/h(i)*(y_7_plus_h - y_7);
  y_prim_7_b = 1/(2*h(i))*(-y_7_plus_2h + 4*y_7_plus_h - 3*y_7);
  
  disp([' ', num2str(h(i), '%0.3f'), ' | ', num2str(y_prim_7_a, '%0.4f'), ' ', ...
    num2str(y_prim_7_a - y_prim_7, '%0.5f'), '   ', ...
    num2str(y_prim_7_b, '%0.4f'), ' ', num2str(y_prim_7_b - y_prim_7, '%0.7f')])
end



% Körning:
% >> EX_13
%        |    approx a)         approx b)
%    h   |  värde    fel      värde     fel
% -------|------------------------------------
%  0.040 | -1.5784 0.10207   -1.6707 0.0096856
%  0.020 | -1.6282 0.05219   -1.6781 0.0023037
%  0.010 | -1.6541 0.02632   -1.6800 0.0004567
%  0.005 | -1.6672 0.01319   -1.6804 0.0000644
