// TP1_08.sce 
clear;
//Exercice 1 : La marche de Wiener  ///////
n=16;T=1;delta_t=T/n;WW0=0;
WW=zeros(n+1,n+1);
WW(1,1)=WW0;
for i=1:n
  WW(i+1,1)=WW(i,1)-sqrt(delta_t); // ici j=0
  for j=1:i
    WW(i+1,j+1)=WW(i,j)+sqrt(delta_t);
  end;
end;

//Exercice 2 : idem pour la marche de Cox
n=100;T=1;delta_t=T/n;
SS=zeros(n+1,n+1);
S0=140;sigma=0.2;
up=exp(sigma*sqrt(delta_t));down=1/up;
SS(1,1)=S0;
for i=1:n
  SS(i+1,1)=SS(i,1)*down; // ici j=0
  for j=1:i
    SS(i+1,j+1)=SS(i,j)*up;
  end;
end;

//Exercice 3 : Complement: tracé de l'arbre de Wiener /////////
Abs=zeros(n+1,n+1);Ord=zeros(n+1,n+1);//initialisation
for k=0:n
  for l=0:n-k
    Abs(l+1,k+1)=(k+l)*delta_t;
    Ord(l+1,k+1)=WW(k+l+1,k+1);
  end;
  for l=1:k
    Abs(n-k+l+1,k+1)=(n-l)*delta_t;
    Ord(n-k+l+1,k+1)=WW(n-l+1,k-l+1);
  end;
end;
xset("window",0);
plot2d(Abs,Ord)
//  tracé de l'arbre de Cox /////////
Abs=zeros(n+1,n+1);Ord=zeros(n+1,n+1);//initialisation
for k=0:n
  for l=0:n-k
    Abs(l+1,k+1)=(k+l)*delta_t;
    Ord(l+1,k+1)=SS(k+l+1,k+1);
  end;
  for l=1:k
    Abs(n-k+l+1,k+1)=(n-l)*delta_t;
    Ord(n-k+l+1,k+1)=SS(n-l+1,k-l+1);
  end;
end;
xset("window",1);
plot2d(Abs,Ord)//xs2eps(1,'ArbreCox.eps');
//xs2gif(1,'ArbreCox.gif');

// Exercice 6 : Simulations de M trajectroires de Wiener /////////////////
p=.5;M=100;
deltasJ=int(p+rand(n,M,'uniform'));
J=cumsum(deltasJ,"r");
Traj=zeros(1+n,M);
for m=1:M
  Traj(1,m)=0;
  for i=1:n
    Traj(1+i,m)=WW(i+1,J(i,m)+1);
  end;
end;
xset("window",6);
plot2d(0:n,Traj)

//Exercice 5 : idem pour M trajectoires CRR
deltasJ=int(p+rand(n,M,'uniform'));
J=cumsum(deltasJ,"r");
Traj=zeros(1+n,M);
for m=1:M
  Traj(1,m)=S0;
  for i=1:n
    Traj(1+i,m)=SS(i+1,J(i,m)+1);
  end;
end;
xset("window",5);
plot2d(0:n,Traj)
xs2eps(5,'TrajAlea3.eps');

//Exercice 7 : :  histogramme des valeurs finales /////
xset("window",2);
histplot(int(sqrt(M)),Traj(n+1,1:M))
//xs2eps(2,'Histog.eps');