//TP8: Options Américaines
//Exercice 1
clear;
n=50;T=1;delta_t=T/n;S0=140;K=S0;sigma=0.4;r=0.25;
up=exp(sigma*sqrt(delta_t));down=exp(-sigma*sqrt(delta_t));
R=exp(r*delta_t);p=(R-down)/(up-down);
//Calcul des valeurs de l'actif sous jacent S
function s=S(i,j);s=S0*up.^j.*down.^(i-j);endfunction;
// Calcul par récurrence de la matrice CC du Call
CC=zeros(n+1,n+1);
for j=0:n CC(1+n,1+j)=max(S(n,j)-K,0);
end;
for i=n-1:-1:0
    for j=0:i CC(1+i,1+j)=(p*CC(1+i+1,1+j+1)+(1-p)*        CC(1+i+1,1+j))/R;
    end;
end;
// Calcul par récurrence de la matrice PP du Put
PP=zeros(n+1,n+1);
for j=0:n PP(1+n,1+j)=max(K-S(n,j),0);
end;
for i=n-1:-1:0
    for j=0:i PP(1+i,1+j)=(p*PP(1+i+1,1+j+1)+(1-p)*        PP(1+i+1,1+j))/R;
    end;
end;

//Exercice 2
//Calcul par récurrence retrograde de la matrice CCamer du Call américain
CCamer=zeros(n+1,n+1); 
for j=0:n CCamer(n+1,j+1)=max(S(n,j)-K,0); end;
for i=n-1:-1:0
    for j=0:i
      CCamer(i+1,j+1)=max(max(S(i,j)-K,0),(p*CC(i+2,j+2)+(1-p)*CC(i+2,j+1))/R);
    end;
end;

//Exercice 3 
CC(1,1)
CCamer(1,1)

//Exercice 4
//Calcul par récurrence retrograde de la matrice PPamer du Put américain
PPamer=zeros(n+1,n+1); 
for j=0:n, PPamer(n+1,j+1)=max(K-S(n,j),0); end;
for i=n-1:-1:0
    for j=0:i
      PPamer(i+1,j+1)=max(max(K-S(i,j),0),(p*PPamer(i+2,j+2)+(1-p)*PPamer(i+2,j+1))/R);
    end;
end;
PPamer(1,1)-PP(1,1)

PP=zeros(n+1,n+1);
for j=0:n PP(1+n,1+j)=max(K-S(n,j),0);
end;
for i=n-1:-1:0
    for j=0:i PP(1+i,1+j)=(p*PP(1+i+1,1+j+1)+(1-p)*        PP(1+i+1,1+j))/R;
    end;
end;
PPamer=zeros(n+1,n+1); 
for j=0:n, PPamer(n+1,j+1)=max(K-S(n,j),0); end;
for i=n-1:-1:0
    for j=0:i
      PPamer(i+1,j+1)=max(max(K-S(i,j),0),(p*PPamer(i+2,j+2)+(1-p)*PPamer(i+2,j+1))/R);
    end;
end;
PPamer(1,1)-PP(1,1)

//Exercices 5 et 6 
K=140;
for r=0:0.05:0.25
    R=exp(r*delta_t);p=(R-down)/(up-down);
PP=zeros(n+1,n+1);
for j=0:n PP(1+n,1+j)=max(K-S(n,j),0);
end;
for i=n-1:-1:0
    for j=0:i PP(1+i,1+j)=(p*PP(1+i+1,1+j+1)+(1-p)*        PP(1+i+1,1+j))/R;
    end;
end;
PPamer=zeros(n+1,n+1); 
for j=0:n PPamer(n+1,j+1)=max(K-S(n,j),0); end;
for i=n-1:-1:0
    for j=0:i
      PPamer(i+1,j+1)=max(max(K-S(i,j),0),(p*PPamer(i+2,j+2)+(1-p)*PPamer(i+2,j+1))/R);
    end;
end;
//r
PP(1,1)
PPamer(1,1)
end;

