u0'=-mu*u0+r*(birth-u0)
u[1..49]'=-mu*u[j]+r*(u[j-1]-u[j])
u50'=-mu1*u50+r*u49
total=sum(0,50)of(shift(u0,i'))
resource=sum(0,50)of(rr(i')*shift(u0,i'))
rr(x)=1/(1+exp(-(x-xt)/xshp))
par xshp=5,xt=30,cap=20,b=.2,mu=.05,r=.8,mu1=2
birth=max(b*total*(1-resource/cap),0)
init u0=2
@ total=500,dt=.25
done