// Question 3
T=1; // an
S_0=100; 
r=0.05; // par an
sigma=0.3; // par racine d'annee sigma^2 * T  est sans dimension
K=100;

function [y]=put(x,K)
  y=max(K*ones(x)-x,0);
endfunction

function []=test_call_arbitrage(N)
// C-P = S_0 - K exp(-rT)
// On peut donc construire un nouvel estimateur
// S_0 - K exp(-rT) + exp(-rT) * (K-S_T)_+

  W_T=sqrt(T)*rand(1,N,"gauss");
  S_T=S_0*exp((r-sigma^2/2)*T + sigma*W_T);
  payoff= ... ; // utilisez le calcul Monte-Carlo du put

  estimation=mean(payoff);  // estimation de la moyenne
  ecart_type=stdev(payoff); // estimation de l'ecart type
  erreur=1.96*ecart_type/sqrt(N); // demi-largeur de l'intervalle de confiance

  printf("CallPut N=%d, %f +- %f\n",N, estimation, erreur);
endfunction

K=100;test_call(1000);test_call_arbitrage(1000);
K= 80;test_call(1000);test_call_arbitrage(1000);
K= 60;test_call(1000);test_call_arbitrage(1000);