Slightly more Sage-ified version of the above very nice solution:
import scipy.integrate
a=1.0
b=2.0
def fun(t):
if t<=-b:
return -a
elif f<b:
return t*a/b
else:
return a
g=lambda t:fun(t)
N=100
time_step=0.1
time_end=10.0
t0=0.0
x0=[[0.5*k,0.5*k] for k in range(-10,10)]
def f(x,t):
return [x[1],-x[0]+g(x[0])]
time_range=[t0..time_end, step=time_step]
x0=[[0.5*k,0.5*k] for k in range(-10,10)]
sol_lines = Graphics()
for n in range(10):
sol = scipy.integrate.odeint(f,x0[n],time_range)
sol_lines += line(sol,rgbcolor=hue(.3+n/15.0))
x0,x1=var('x0 x1')
p=plot_vector_field ((x1,-x0+g(x0)),(x0,-9,9),(x1,-7,7))
show(sol_lines + p, figsize = [9,7])
--
To post to this group, send email to [email protected]
To unsubscribe from this group, send email to
[email protected]
For more options, visit this group at
http://groups.google.com/group/sage-support
URL: http://www.sagemath.org