k=.033;g=50;
t1=clock;
a=k+(1-k)/g;n2=840;c=zeros(1,n2);p1=zeros(1,n2);x=1:n2;
for j=1:n2
	t(j)=j*k/(1-k);
	c(j)=gamma(j/(1-k)+1)/gamma(j+1)*gamma(t(j)+2);
	p1(j)=a^(t(j)+1)*(1-a)^(t(j)-1)*c(j);
end
loglog(x,p1)
%plot(x,p1)
duratam=etime(clock,t1)
%stem(p1)
%p1-resd,