% Plot f(x) and antiderivative
% fint(x) = \int_a^x f(t) dt
% On interval [a,b]

% Clear all previously defined variables
clear all
% Cloase all figures
close all

% Define interval [a,b], inital integration point x0
a=0;
b=1;
x0=0;

% Define function fder(x)
fder = @(x) x.^2; a=-4; b=4;
% fder = @(x) exp(-x.^2); a = -3; b=3;
% fder = @(x) sin(x.^2);  a=-5; b=5;
% fder = @(x) 1./(1+x.^2); b=10;
% fder = @(x) sin(x);

% define fint(x,y) = \int_x^y f(t) dt 
fint =  @(x,y) quadl(fder,x,y);

% Plot function
nplot = 1000;
dxplot = (b-a)/nplot;
xplot= a:dxplot:b;
yplot = fder(xplot);
hold off
plot(xplot,yplot);
xlabel('x');
ylabel('y');
pause

% Compute values of antiderivate fint(x)
% normalized so that fint(x0) = 0
if a==x0
    yintplot(1)=0;
else    
    yintplot(1) = fint(x0,a);
end
for i=1:nplot
  yintplot(i+1) = yintplot(i) + fint(xplot(i),xplot(i+1));
end

% Plot antiderivative 
hold on
plot(xplot,yintplot, 'r')
title('y=f(x) (blue) and y = F(x) (red),where F is an integral of f')
    
hold off