//paramètres
T=1;//horizon
sig=0.2;//volatilité
r=0.05;//intérêt
S0=100;//valeur initiale du sous-jacent
N=2;//nombre initial de pas de discrétisation
M=10000;//nombre de simulations indépendantes

erreul=[];//vecteur des erreurs fortes Euler
errmil=[];//vecteur des erreurs fortes Milshtein
liceul=[];//largeur des intervalles de confiance Euler
licmil=[];//largeur des intervalles de confiance Milshtein
Npas=[];// vecteur des nombres de pas

for j=1:5, //boucle sur le nombre de pas
  //paramètres utiles pour la discrétisation avec N pas
  
  //////////////////////////////////////////////////
  //A compléter avec le calcul des paramètres utiles
  //////////////////////////////////////////////////
  
  //variables de stockage
  someul=0;
  careul=0;
  sommil=0;
  carmil=0;
  
  
  for i=1:M, //simulations indépendantes
	S=S0;
	Se=S0;
	Sm=S0;
	for k=1:N, //boucle sur les pas de temps
	  g=rand(1,'g');//génération d'une gaussienne centrée réduite
	  
	  ////////////////////////////////////////////////////////////
	  //À compléter avec l'evolution du sous-jacent et des schémas
	  ////////////////////////////////////////////////////////////
	 
	end;
	someul=someul+(S-Se)^2;
	careul=careul+(S-Se)^4;
	sommil=sommil+(S-Sm)^2;
	carmil=carmil+(S-Sm)^4;
  end;
  erreul=[erreul,someul/M];
  liceul=[liceul,1.96*sqrt((careul/M-(someul/M)^2)/M)];
  errmil=[errmil,sommil/M];
  licmil=[licmil,1.96*sqrt((carmil/M-(sommil/M)^2)/M)];
  Npas=[Npas,N];
  
  N=N*2;//multiplication du nombre N de pas par 2
end;

//Affichage des vecteurs d'erreurs et des demi-largeurs des IC à 95%
erreul
liceul
errmil
licmil


//representation graphique de (1/N,erreul) et (1/(N*N),errmil)
xbasc();
subplot(2,1,1);
plot2d(1../Npas',[erreul;erreul-liceul;erreul+liceul]');
subplot(2,1,2);
plot2d(1../Npas^2',[errmil;errmil-licmil;errmil+licmil]');