function mcm1test09
% Figs. 9.8a, 9.8b mcm1test:  Monte Carlo Method (2009), 
%     nx = 1 dim, uniform dist, 
%     I = int(h(x),x=a..b), h(x) = sqrt(1-x^2), -1 <= a < b <= +1;
%     technically, g(x) = (b-a)h(x) = (b-a)*sqrt(1-x^2) to account for 
%     uniform density f(x) = 1/(b-a) on [a,b], so I = meang.
%
clc; clear
%
fprintf('Monte Carlo Test of 1-dim Uniform Dist. on (a,b)');
fprintf('\n   with h(x)=sqrt(1-x^2) and g(x) = (b-a)h(x):\n');
a = 0; b = +1; % -1 <= a < b <= +1;
% integral of h(x) = sqrt(1-x^2); on [a,b]:
IntExact = 0.5*(asin(b)-asin(a))+0.5*(b*sqrt(1-b^2)-a*sqrt(1-a^2)); 
MugExact = IntExact;
Sigg = sqrt((b-a)^2*(1-(b^2+a*b+a^2)/3)-MugExact^2);
fprintf('\nk        n muhatn mugExact   sighatn     Sigg    stderrn  AbsErrorg\n');
kmax = 7; 
n = zeros(1,kmax); meang = zeros(1,kmax); sigg = zeros(1,kmax); 
sigdrn = zeros(1,kmax); error = zeros(1,kmax);
for k = 1:kmax
   rand('state',0); % set state or seed
   n(k) = 10^k; % sample size, k = log10(n(k)) ;
   x = a+(b-a)*rand(n(k),1); % get n(k) X 1 random sample on (a,b);
   g = (b-a)*sqrt(1-x.^2); % vectorized g;
   meang(k) = mean(g); % E[g(X)];
   sigg(k) = std(g); % sqrt(sigmag^2), sigmag^2 = unbiased variance of g;
   sigdrn(k) = sigg(k)/sqrt(n(k));
   error(k) = abs(meang(k)-MugExact);
   fprintf('%1i %8i %6.4f   %6.4f %9.3e %9.3e  %9.3e %9.3e\n'...
       ,k,n(k),meang(k),MugExact,sigg(k),Sigg,sigdrn(k),error(k))
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
scrsize = get(0,'ScreenSize');
ss = [3.0,2.8,2.6,2.4,2.2,2.0];
%
nfig = 1;
figure(nfig);
kv = 1:kmax;
plot(kv,meang,'k-o','linewidth',3,'MarkerSize',12);
hold on
plot(kv,sigg,'k-x','linewidth',3,'MarkerSize',12);
axis([min(kv) max(kv) 0 1]);
title('Monte Carlo Results, Uniform Dist., g(x) = (b-a)*sqrt(1-x^2)'...
    ,'Fontsize',36,'FontWeight','Bold');
xlabel('log(n), Log_{10} Sample Size'...
    ,'Fontsize',36,'FontWeight','Bold');
ylabel('g-Moments \mu_n, \sigma_n'...
    ,'Fontsize',36,'FontWeight','Bold');
legend('\mu_n, Mean-est.','\sigma_n, StdDev-est.','Location','Best');
set(gca,'Fontsize',32,'FontWeight','Bold','linewidth',3);
set(gcf,'Color','White','Position' ...
    ,[scrsize(3)/ss(nfig) 70 scrsize(3)*0.60 scrsize(4)*0.80]);
hold off
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
nfig = nfig+1;
figure(nfig);
kv = 1:kmax;
plot(kv,log10(sigdrn),'k-o','linewidth',3,'MarkerSize',12); 
hold on
plot(kv,log10(error),'k-x','linewidth',3,'MarkerSize',12); 
% ymin = min(min(log10(sigdrn)),min(log10(error))); 
% ymax = max(max(log10(sigdrn)),max(log10(error)));
axis tight; %axis([min(kv) max(kv) ymin ymax]);
title('Monte Carlo Errors, Uniform Dist., g(x) = (b-a)*sqrt(1-x^2)'...
    ,'Fontsize',36,'FontWeight','Bold');
xlabel('log(n), Log_{10} Sample Size'...
    ,'Fontsize',36,'FontWeight','Bold');
ylabel('g-Errors log(StdError_n), log(AbsError_n)'...
    ,'Fontsize',36,'FontWeight','Bold');
legend('log_{10}(StdError_n)','log_{10}(AbsError_n)','Location','Best');
set(gca,'Fontsize',32,'FontWeight','Bold','linewidth',3);
set(gcf,'Color','White','Position' ...
    ,[scrsize(3)/ss(nfig) 70 scrsize(3)*0.60 scrsize(4)*0.80]); 
hold off
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% end mcm1test.m
%