v[1..20]'=w[j]*(r*v[j]-p*x[j]*v[j]-q*z*v[j]+eps)
x[1..20]'=c*v[j]-b*x[j]-u*v*x[j]
z'=k*v-b*z-u*v*z
par r=2.2,p=2,c=1,b=.01
par u=1,q=2.4,k=1
par pp=.1,eps=.001
w0=1
v=sum(0,19)of(shift(v1,i'))
aux vload=v
aux div=sum(0,19)of(shift(v1,i')^2/max(v,.00001)^2)
markov w[1..20] 2 0
{0} {pp*w[j-1]}
{0} {0}
@ total=200,meth=euler,bound=100000
@ xhi=200,yp=vload,ylo=0,yhi=5
done