3238.sagews

%auto
# make sure to evaluate this cell to load the functions
URL='http://sage.math.washington.edu/home/jason/sage-field-plots/'
files=['lic_internal.pyx', 'plotfields.py','streamline.py']
for f in files:
    load(URL+f)

Compiling /home/sagenbws/.sage//temp/mod.math.washington.edu/549//tmp_0.pyx...

var('y')
F=[-1 - x^2 + y, 1 + x - y^2]
v=plot_vector_field(F,(x,-3,3),(y,-3,3),aspect_ratio=1)
v

q=plot_streamlines(F, (x,-3,3), (y,-3,3))
q

p=plot_field_texture(F, (x,-3,3), (y,-3,3),plot_points=400,alpha=0.9)
p

p+q

%auto
# make sure to evaluate this cell to load the functions
    
# based on code from oddrobot, http://sagenb.org/home/pub/1532/
from sage.calculus.desolvers import desolve_system_rk4
x,y,t=var('x y t')
class DESolution:
    def __init__(self,system,time_range,initial,stepsize=0.05):
        
        self.tvar=time_range[0]
        self._times=srange(time_range[1],time_range[2],stepsize)
        self.vars=[v for v,_ in initial]
        self.dim=len(self.vars)
        self._soln=desolve_odeint(system, ics=[v for _,v in initial], times=self._times, dvars=self.vars, ivar=self.tvar)
        
    def phase_plot(self,vars=None,color='blue',**kwargs):
        # find which indices the specified variables are
        if vars is not None:
            vars_index=[self.vars.index(v) for v in vars]
        elif self.dim<=3:
            vars_index=range(self.dim)
        else:
            vars_index=range(2)
            
        p=line(self._soln[:,vars_index],color=color,**kwargs)
        # add an arrow head showing which way we are going around the phase line
        half=int(self._soln.shape[0]/2)
        p+=arrow(self._soln[half,vars_index], self._soln[half+1,vars_index],color=color)
        if len(vars_index)==2:
            p.axes_labels([str(self.vars[v]) for v in vars_index])
        return p
        
    def coordinates(self,colors=None,**kwargs):
        if colors is None:
            colors=rainbow(len(self.vars))
        p=Graphics()
        p+=sum(line(zip(self._times,self._soln[:,i]), color=colors[i], legend_label=str(self.vars[i]),**kwargs) for i in range(self.dim))
        return p

var('x,y')
F=[x-y^2,x-y]
solution=DESolution(F,[t,0,4],[(x,1),(y,0.5)])
solution.coordinates(['red','blue'])

solution.phase_plot([x,y])

plot_streamlines(F,(x,0,2),(y,0,2))+solution.phase_plot([x,y],thickness=3)

plot_field_texture(F,(x,0,2),(y,0,2),alpha=0.8)+solution.phase_plot([x,y],thickness=3)