tic
broj_cestica=10;
broj_odeljaka=5;
broj_ciklusa=10000;
raspodela=zeros(broj_odeljaka,broj_cestica+1);
broj_u_odeljcima=zeros(1,broj_odeljaka);
for i=1:broj_ciklusa
 for j=1:broj_cestica
   indeks=floor(rand*(broj_odeljaka))+1;
   broj_u_odeljcima(indeks)=broj_u_odeljcima(indeks)+1;
 end
 for k=1:broj_odeljaka
  raspodela(k,broj_u_odeljcima(k)+1)=raspodela(k,broj_u_odeljcima(k)+1)+1;
 end
 broj_u_odeljcima=zeros(1,broj_odeljaka);
end
raspodela=raspodela./broj_ciklusa;
xosa=0:broj_cestica;
%plot(xosa,sum(raspodela)/broj_odeljaka,'o-')
figure
plot(xosa,sum(raspodela)/broj_odeljaka,'o-',xosa,binopdf(xosa,broj_cestica,1/broj_odeljaka),'*r-')
xlabel('Broj cestica');
ylabel('Verovatnoca');
legend ('Raspodela cestica','Binomna raspodela')
toc