function mcm09unifnorm % Example Output: FINM 331 Demo for Uniform Monte Carlo; % for uniform rand-deviates on (a,b); clc clear % global a b V % fprintf(' Uniform rand-deviate Monte Carlo:\n'); n = 10000; srtn = sqrt(n); a=-3; b = +2; V = b-a; fprintf('\nn=%i; a=%6.4f; b=%6.4f; V=%6.4f;\n',n,a,b,V); % erfc(x) = 2\sqrt(pi)*int(exp(-t^2),t=x..infty); % erf(x) = 2\sqrt(pi)*int(exp(-t^2),t=0..x); exact = 0.5*(erfc(a/sqrt(2))-erfc(b/sqrt(2))); sig2unif_exact = 2.5/sqrt(pi)*(erf(b)-erf(a))-exact^2; fprintf('\nexact integral = %10.8f;',exact); fprintf('\nsig2unifexact = %9.4e; sigunifexact = %9.4e;'... ,sig2unif_exact,sqrt(sig2unif_exact)); rand('twister',3); U = a+V*rand(1,n); fuv=fu(U); % Monte Carlo estimators: mu_unif = mean(fuv); fprintf('\nmu_unif=%10.8f;',mu_unif); fprintf('\nmu_unif_abserror=%9.4e;',mu_unif-exact); fprintf('\nmu_unif_relerror=%9.4e%%;\n',100*(mu_unif/exact-1)); % Monte Carlo variance estimators: sig2unif = var(fuv); % MATLAB var(x); gives unbiased variance fprintf('\nsig2_unif=%9.4e;',sig2unif); fprintf('\nsig2_unif_abserror=%9.4e;',sig2unif-sig2unif_exact); fprintf('\nsig2unif_relerror=%9.4e%%;\n'... ,100*(sig2unif/sig2unif_exact-1)); % std. errors: se_unif_exact = sqrt(sig2unif_exact)/srtn; se_unif = sqrt(sig2unif)/srtn; fprintf('\nstderr_unif_exact=%9.4e;',se_unif_exact); fprintf('\nstderr_unif=%9.4e;',se_unif); fprintf('\nstderr_unif_diff=%9.4e;\n\n',se_unif-se_unif_exact); % %%%%% function y = fu(x) global V y = V*exp(-x.*x/2)/sqrt(2*pi); %%%%% % % end mcm09unifnorm.m %