# SchnakenRD.ode
# Reaction diffusion system

f(u,v)=gamma*(a-u+u^2*v)
g(u,v)=gamma*(b-u^2*v)

#********* Periodic Boundaries ********
u1'=f(u1,v1)+du*(u2-2*u1+u99)/h^2
v1'=g(u1,v1)+dv*(v2-2*v1+v99)/h^2


#Neumann (no flux) Boundaries
#u1'=f(u1,v1)+du*(u2-u1)/h^2
#v1'=g(u1,v1)+dv*(v2-v1)/h^2
%[2..99]
u[j]'=f(u[j],v[j])+du*(u[j+1]+u[j-1]-2*u[j])/h^2
v[j]'=g(u[j],v[j])+dv*(v[j+1]+v[j-1]-2*v[j])/h^2
%

#********* Periodic Boundaries ********
u100=u1
v100=v1

#Neumann (no flux) Boundaries
#u100'=f(u100,v100)+du*(u99-u100)/h^2
#v100'=g(u100,v100)+dv*(v99-v100)/h^2

par a=.2,b=2.0,gamma=0.2,du=0.01,dv=1,h=.1
@ meth=cvode,bandup=2,bandlo=2,total=100
done