%% class 20 October 2017

% numerical integrals
clear all
close all
clc

N=(10:1:100);
g=@(x)(exp(x));
for i=1:length(N)
    a=0; b=1; 
    x=linspace(a, b, N(i));
    dx(i)=x(2)-x(1);
    f=exp(x);
    % let's introduces the indices il and ir
    il=(1:1:N(i)-1); ir=(2:1:N(i));

    Il(i)=sum(f(il)*dx(i));
    Ir(i)=sum(f(ir)*dx(i));
    Itr(i)=sum(dx(i)*(f(il)+f(ir))/2);
%    Ic(i)=sum(exp((x(il)+x(ir))/2)*dx(i));
% or equivalently
    Ic(i)=sum( g((x(il)+x(ir))/2)*dx(i) );
end

hold on
plot(dx,ones(length(dx),1)*(exp(1)-1),'k--')
plot(dx,Il)
plot(dx,Ir)
plot(dx,Itr,'r','Linewidth',2)
plot(dx,Ic,'g','Linewidth',2)
xlabel('dx')
legend('exact integral','Il','Ir','Itr','Ic')
hold off

%% same as above by using trapz(x,y)
clear all
close all
clc
format long
a=0; b=1; N=100;
x=linspace(a, b, N);
dx=x(2)-x(1);
f=exp(x);
% let's introduces the indices il and ir
il=(1:1:N-1); ir=(2:1:N);

exact_int=exp(1)-1

Il=sum(f(il)*dx)
Ir=sum(f(ir)*dx)
Itr=sum(dx*(f(il)+f(ir))/2)
Ic=sum(exp((x(il)+x(ir))/2)*dx)

%using trapz
trapz(x,f)

%% we want to do now indefinite integrals
% by using cumsum: example
clear all
clc

x=(0:1:20);
y=cumsum(x)
length(y)

y(1)
y(2)
y(3)

%% excercise 8.0.17
clear all
clc
close all

% for a single value of x:
x=2.1*pi;

a=pi/2; b=x; N=500;
s=linspace(a, b, N);
ds=s(2)-s(1);
f=cos(s);
%Itr2=trapz(s,f)
% even better with our definition
il=(1:1:N-1);ir=(2:1:N);
Itr=sum( (f(il)+f(ir))/2 *ds)

%compare with the exact result
sin(x)-1

% for more values of x:
x=linspace(pi/2,5*pi,200);

for i=1:length(x)
    a=pi/2; b=x(i); N=500;
    s=linspace(a, b, N);
    ds=s(2)-s(1);
    f=cos(s);
%Itr2=trapz(s,f)
% even better with our definition
    il=(1:1:N-1);ir=(2:1:N);
    Itr(i)=sum( (f(il)+f(ir))/2 *ds);
end

hold on
plot(x,Itr)
plot(x,sin(x)-1,'rx')
legend('approx','exact','Location','Best')
xlabel('x')
hold off

% note that this way the accuracy for evaluating the integral at lower x is
% higher than for a larger x!!!!!

%% excercise 8.0.17 with cumsum
clear all
clc
close all

N=500;
x=linspace(pi/2,5*pi,N);
dx=x(2)-x(1);
il=(1:1:N-1); ir=(2:1:N);
f=cos(x);

Itr(1)=0
Itr(ir)=cumsum((f(il)+f(ir))/2 *dx);

hold on
plot(x,Itr)
plot(x,sin(x)-1,'rx')
legend('approx','exact','Location','Best')
xlabel('x')
hold off

% now Itr is evaluated with the same accuracy no matter of the value of
% x!!!

%% excercise 7 of week 2
clear all
close all
clc

f=@(x)( 2/5 - 1/5*exp(-x).*(2*cos(2*x)+sin(2*x)) );

g=@(s)(sin(2*s).*exp(-s));

a=0; b=4*pi; N=200;
x=linspace(a, b, N);

dx=x(2)-x(1);
il=(1:1:N-1); ir=(2:1:N);

Itr(1)=0
Itr(ir)=cumsum( (g(x(il))+g(x(ir)))/2 *dx);

figure()
hold on
plot(x,f(x),'r')
plot(x,Itr,'xg')
legend('exact integral','numerical integral')
xlabel('x')
hold off






