On Mar 20, 9:55 pm, Marshall Hampton <hampto...@gmail.com> wrote:
> 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])

Marshall, is this something that could be wrapped easily into a Sage
function (perhaps along with plot_slope_field) for general phase
portraits, or is it only going to work for certain types of 'nice'
DEs?

- kcrisman

-- 
To post to this group, send email to sage-support@googlegroups.com
To unsubscribe from this group, send email to 
sage-support+unsubscr...@googlegroups.com
For more options, visit this group at 
http://groups.google.com/group/sage-support
URL: http://www.sagemath.org

Reply via email to