beta=5;
N=1000;

function [z]=g(x);
  y=0;
  if (x>0) then
    y=exp(beta*x);
 end,
 z=y;
endfunction
    
function [z]=ftilde(x);
  z=exp(-(x-beta)^2/2)/sqrt(2*%pi);
endfunction

function [z]=f(x);
  z=exp(-x*x/2)/sqrt(2*%pi);
endfunction

//méthode directe  
s=0;
for i=1:N
  x=grand(1,1,'nor',0,1);
  s=s+g(x);
end,
s=s/N;
disp('estimation par methode directe');
s
s=0;
for i=1:N
  x=grand(1,1,'nor',0,1);
  y=grand(1,1,'nor',0,1);
  s=s+(g(x)-g(y))^2;
end,
s=s/(2*N);
disp('variance par methode directe');
s

//méthode d'échantillonage préférentiel
s=0;
for i=1:N
  x=grand(1,1,'nor',beta,1);
  s=s+g(x)*f(x)/ftilde(x);
end,
s=s/N;
disp('estimation par echantillonage preferentiel');
s
s=0;
for i=1:N
  x=grand(1,1,'nor',beta,1);
  y=grand(1,1,'nor',beta,1);
  x1=g(x)*f(x)/ftilde(x);
  y1=g(y)*f(y)/ftilde(y);
  s=s+(x1-y1)^2;
end,
s=s/(2*N);
disp('variance de la methode');
s

//méthode de variable de contrôle
s=0;
for i=1:N
  x=grand(1,1,'nor',0,1);
  s=s+g(x)-exp(beta*x);
end;
s=s/N+exp(beta^2/2); //à vérifier
disp('estimation par variable de controle');
s
s=0;
for i=1:N
  x=grand(1,1,'nor',0,1);
  y=grand(1,1,'nor',0,1);
  s=s+(g(x)-exp(beta*x)-g(y)+exp(beta*y))^2;
end,
s=s/(2*N);
disp('variance de la methode');
s

//méthode de variable de contrôle + variables antithétiques
s=0;
for i=1:N
  x=grand(1,1,'nor',0,1);
  y=-x;
  s=s+g(x)-exp(beta*x)+g(y)-exp(beta*y);
end,
s=s/(2*N)+exp(beta^2/2);
disp('estimation par variable de controle + variables antithetiques');
s
s=0;
for i=1:N
  x=grand(1,1,'nor',0,1);
  y=-x;
  x1=grand(1,1,'nor',0,1);
  y1=-x1;
  s=s+(g(x)-exp(beta*x)+g(y)-exp(beta*y)-(g(x1)-exp(beta*x1)+g(y1)-exp(beta*y1)))^2;
end,
s=s/(8*N);
disp('variance de la methode');
s



  
   

    
