I see that I did make that simplification. You can get some improvement by cutting the stepsize in the RK4 algorithm down. This can by done:

```
y=function('y',x)
f(x,y)=1-2*x*y
soln=desolve_rk4(f(x,y),y,ics=[0,0],ivar=x,end_points=16,stepsize=0.001)
p=line(soln)
p+=plot_vector_field([1/sqrt(1+f(x,y)^2),f(x,y)/sqrt(1+f(x,y)^2)],(x,0,20),(y,0,0.8))
show(p)
```

You can see that near $x=16$ the approximation starts falling apart. With this sensitivity, I recommend a more robust solver. To use `desolve_odeint`

, you need to have a system. So, below I use $dx/dt=1, dy/dt=1-2xy$ as the system.

```
var('x,y')
g(x,y)=1-2*x*y
f=[SR(1),g(x,y)]
soln=desolve_odeint(f,[0,0],srange(0,20,0.05),[x,y])
p=line(soln)
p+=plot_vector_field([1/sqrt(1+g(x,y)^2),g(x,y)/sqrt(1+g(x,y)^2)],(x,0,20),(y,0,0.8))
show(p)
```

(The `SR(1)`

ensures that Sage knows `1`

is a symbolic function. `desolve_odeint`

appears to need this.)

To get your plot to start at $x=-10$, you want to start at (0,0) and run back to $x=-10$. Then, combine this with your other plot.

```
var('x,y')
g(x,y)=1-2*x*y
f=[SR(1),g(x,y)]
soln1=desolve_odeint(f,[0,0],srange(0,-10,-0.05),[x,y])
soln2=desolve_odeint(f,[0,0],srange(0,20,0.05),[x,y])
p=line(soln1)+line(soln2)
p+=plot_vector_field([1/sqrt(1+g(x,y)^2),g(x,y)/sqrt(1+g(x,y)^2)],(x,-10,20),(y,-0.8,0.8))
show(p)
```

Also, for fun, you might like playing with `pplane`

at http://math.rice.edu/%7Edfield/dfpp.html. This has a very good solver and is handy for some uses. And, it will solve forwards and backwards. (The original pplane program runs in matlab, but this link has a version that runs in java in your browser.)