function y = BlackScholes_AV_MC_call(S0,K,T,r,sigma,M)
stacksize('max')
    tic();
    Y=rand(1,M/2,'normal');
    X=[Y -Y];
    S=S0*exp(sigma*X*sqrt(T)+(r-1/2*sigma^2)*T);
    C=exp(-r*T)*max(S-K,0);
    Price=1/M*sum(C);
    VarEst=1/(M-1)*sum((C-Price)^2);
    CI95=[Price-1.96*sqrt(VarEst)/sqrt(M),Price+1.96*sqrt(VarEst)/sqrt(M)];
    CI99=[Price-2.58*sqrt(VarEst)/sqrt(M),Price+2.58*sqrt(VarEst)/sqrt(M)];
    time=toc();
    y=[Price 0; CI95; CI99; time 0]
endfunction