function r = binornd(n,p,dims)
% generates a matrix of binomial random numbers from the 
% binomial distribution with parameters n and p
% INPUTS
%  n:  number of simulated Bernoulli trials (e.g. coinflips)
%  p:  probability of a 1
%  dims = list of dimensions for the output array r
% OUTPUT
%  r:  array of sample counts of the Binomial distribution
r = sum(rand([dims,n]) < p, length(dims)+1);