Translating a Sage Example for EMT

The nice Sage Math project is something I am very interested in. It is a rather complete package for mathematical computations based on Python. If you are interested in open source math this is probably one good way to start. I must add however that I prefer open source packages directly based on Python with a command line (like IPython) instead of a web based interface that covers everything with another layer. The advantage of the latter is, of course, that it may provide additional services, e.g., for cooperation.

On Windows you will have to run Sage Math in a Virtual Box. You can then open the web interface in your normal Windows browser and work there.

I will certainly come back in this blog to Sage Math. For a start, let me translate an example of Sage Math into Euler Math Toolbox (EMT). It is an example for the famous Lotka-Volterra equations, representing the relation between predator and prey. As far as I can tell, Python uses a fourth order Runge-Kutta method to solve the differential equation numerically. The equations themselves are written in symbolic form. The package contains features for symbolic computations.

### Simulates Lotka-Volterra population dynamics, producing a vector field
### Arrows indicate the direction of population evolution
### The blue contour is one potential stable loop
var('t R F R0 F0 a b c d')
## CONSTANTS
a = 0.04         # rabbit birthrate
b = 0.0005       # probability of predation per encounter
c = 0.1*0.0005   # rabbit -> fox conversion efficiency
d = 0.2          # death rate of foxes
## INITIAL CONDITIONS
R0 = 5000     #initial number of rabbits
F0 = 200      #initial number of foxes

## DIFFERENTIAL EQUATIONS
de_R = (diff(R,t) == a*R - b*R*F)
de_F = (diff(F,t) == c*R*F - d*F)

## CALCULATION PARAMETERS
end_points = 500
stepsize = 1.0
steps = end_points/stepsize
ics = [0,R0,F0]
des = [de_R.rhs(), de_F.rhs()]

## NUMERICAL SOLUTION OF DIFFERENTIAL EQUATIONS
sol = desolve_system_rk4(des,[R,F],ics,ivar=t,
   end_points=end_points,step=stepsize)

## Clean up to graph
sol_t=[]; sol_R=[]; sol_F=[]
for i in range(steps):
    sol_t.append(sol[i][0])
    sol_R.append(sol[i][1])
    sol_F.append(sol[i][2])

a = plot([],figsize=[6,4])
a += line(zip(sol_R,sol_F))
a += plot_vector_field((des[0], des[1]), (R,0,7000), (F,0,225),
  xmin=1500,color='orange')
a.axes_labels(['Rabbits','Foxes']); show(a)

I entered that as is into the Sage Math browser and the following plot was produced. I had to copy the image out of the browser, because I have no idea on how to save the image with the browser interface. Let me also note that I do not understand every command above, though many of them are self explanatory.

1

Let us translate this too Euler Math Toolbox. We use the numerical side of EMT only, although we could have defined some functions in symbolic form. But it would not help much to do so.

>a=0.04; b=0.0005; c=0.1*b; d=0.2;
>R0=5000; F0=200;
>function dR ([R,F]) := a*R-b*R*F;
>function dF ([R,F]) := c*R*F-d*F;
>function f(x,y) := [dR(y),dF(y)];
>t=0:1:500;
>s=ode("f",t,[R0,F0]); 
>{x,y}=vectorfield2("dR","dF",min(s[1]),max(s[1]),min(s[2]),max(s[2]), ...
>   scale=0.005,<plot); 
>plot2d(s[1],s[2],color=blue, ...
>   xl="Rabbits",yl="Foxes");
>plot2d(x,y,>add,style="->",color=orange):

The result is the following, very similar, plot. However, the anti-aliasing is better, but that could be due to my lack of knowledge in Sage Math.

2

2 Gedanken zu „Translating a Sage Example for EMT

  1. Royi

    Hi,

    Could you collaborate with PortableApps to offer this as Portable program?

    It will increase the exposure of Euler.

    Thank You.

    Antworten

Schreibe einen Kommentar

Deine E-Mail-Adresse wird nicht veröffentlicht.

Diese Website verwendet Akismet, um Spam zu reduzieren. Erfahre mehr darüber, wie deine Kommentardaten verarbeitet werden.