// OPTION COLLAR
n=50
T=1
delta_t=T/n
r=0.25
sigma=0.4
S0=120
K=120
up=exp(sigma*sqrt(delta_t))
down=1/up
R=exp(r*delta_t)
p=(R-down)/(up-down)

function s=S(i,j);
        s=S0*up^j.*down^(i-j);
endfunction
function phi=phi(S,K);
    phi=max(S-K,0);
endfunction
function moncall=Call(K);
    moncall=phi(S(n,0:n),K)*binomial(p,n)'/R^n
endfunction
function psi=psi(S,K);
    psi=max(K-S,0);
endfunction
function monput=Put(K)
    monput=Call(K)-S0+K*exp(-r*T);
endfunction
function dif=difference(K)
    dif=Call(S0)-Put(K);
endfunction

Kcollar=fsolve(S0,difference);
disp(Kcollar,"Kcollar=");

// FRONTIERE D'EXERCICE DU PUT AMERICAIN
// EXERCICE 1

n=20
T=1
delta_t=T/n
r=0.25
sigma=0.4
S0=120
K=S0
up=exp(sigma*sqrt(delta_t))
down=1/up
R=exp(r*delta_t)
p=(R-down)/(up-down)


function psu=psi(x)
    psu=max(0,K-x)
endfunction

SS=zeros(n+1,n+1)
SS(1,1)=S0
for i=1:n
    SS(i+1,1)=SS(i,1)*down
    for j=1:i
        SS(i+1,j+1)=SS(i,j)*up
    end
end

PP=zeros(n+1,n+1)
for j=1:(n+1)
	PP(1+n,j)=psi(SS(1+n,j))
end
for i=n-1:-1:0
    for j=0:i
        PP(i+1,j+1)=(p*PP(i+2,j+2)+(1-p)*PP(i+2,j+1))/R;
	end
end

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

PPsi=zeros(n+1,n+1)
for i=1:n
    for j=0:i
        PPsi(i+1,j+1)=psi(SS(i+1,j+1))
    end
end

// EXERCICE 2

EPA=-ones(n+1,n+1)
for j=0:n
    if SS(n+1,j+1)<=K then EPA(n+1,j+1)=1; else EPA(n+1,j+1)=0; end
end;
for i=0:n-1
    for j=0:i
        if psi(SS(i+1,j+1))>(p*PPamer(i+2,j+2)+(1-p)*PPamer(i+2,j+1))/R then EPA(i+1,j+1)=1;
        else EPA(i+1,j+1)=0; end,
    end
end
            

// EXERCICE 3

//for i=0:n
//    for j=0:i
//        if SS(i+1,j+1)<2*S0 then plot2d(i+1,SS(i+1,j+1),-2*EPA(i+1,j+1)-2); end,
//    end;
//end;

// EXERCICE 5
k=0;
Frontieret=zeros(n+1)
FrontiereS=zeros(n+1)
for ii=1:n+1
    if max(EPA(ii,1:ii))==1 then 
        k=k+1; 
        Frontieret(k)=ii;
        FrontiereS(k)=SS(ii,sum(EPA(ii,1:ii)));
    end,
end;
plot(Frontieret,FrontiereS,'g-')
