Of course N=100 in "my" code and repeated x0=[[0.5*k,0.5*k] for k in range(-10,10)] in Marshall one are superfluous :)
Andrzej Chrzeszczyk On 21 Mar, 13:58, kcrisman <[email protected]> wrote: > On Mar 20, 9:55 pm, Marshall Hampton <[email protected]> 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 [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
