function beta=betarv(m,v,n);
% betarv(m,v,n) generates a row vector of length n, the elements of which 
% are Beta random variables with mean m and variance v;
% see G.S. Fishman, 1973, Concepts and Methods in Discrete Event
% Digital Simulation, Wiley, New York, pp. 204-208.
% REQUIRES THE FILE gammarv.m

if v==0
   beta=m*ones(1,n);
elseif v>=m*(1-m)
    disp('Variance of Beta too large given the mean');
    pause;
else
    a=m*(m*(1-m)/v - 1);
	b=(1-m)*(m*(1-m)/v - 1);
	beta=[];
	for i=1:n
        k1=floor(a);
        k2=floor(b);
        if k1==0 & k2==0
            Y=1;
            Z=1;
            while Y+Z>1
                Y=rand^(1/a);
                Z=rand^(1/b);
            end;
        else
            Y=gammarv(a,1,1);
            Z=gammarv(b,1,1);
        end;        
        beta=[beta Y/(Y+Z)];
    end; 
end;