Joaquim Puig's profile - activity

If you had your solution already coded for "desolve_system_rk4" you could try "desolve_system_odeint" which has the same syntax and uses an implicit method by default (you can change the method, of course). It should be much faster than the rk4.

Hi,

In a public SAGE server, it is possible to "publish" it so that non-registered users can view it, edit it after logging and downloading it. It would be nice to be able to publish it in wordpress, blogger or other blogging sites in a more or less simple way. Which is the best way to do it? Thanks in advance

In my opinion the easiest way to do the fit is to use the find_fit function, which is standard in SAGE. For example:

sage: x = [1, 3, 4]   
sage: y = [1.2, 4.5, 3.6]
sage: var('a,b,t')
sage: model(t)=a*t+b
sage: data=zip(x,y)
sage: fit=find_fit(data,model,solution_dict=True) # Dictionary is optional but handy
sage: model.subs(fit) # Show the model
t |--> 0.92142857142840007*t + 0.64285714285636364
sage: plot(model.subs(fit),(t,0,5))+points(data,size=20,color='red') # Nice plot

When you have a huge number of points and you want just a polynomial fit, I found that it is (numerically) better to use the polyfit function from numpy:

sage: import numpy as np
sage: a,b=np.polyfit(x,y,1)
sage: a,b
(0.92142857142857137, 0.64285714285714401)

Hi,

We recently found that our SAGE 4.7 in Ubuntu 11.04 has the following problem with the notebook. We use firefox 6 and everything works fine except that we cannot insert cells, even if the blue line apperas. If we press ALT+ENTER, then we can evaluate the input and produce a new cell, but we would like to be able to insert cells as usual.

Thanks in advance,

Joaquim

I have now noticed that after version 4.7.1 (maybe before) it is not necessary to recompile SAGE.

Hi,

I don't know your problem, but sometimes it is possible to introduce the differential equation for A(t) in the set of equations that you pass to the ode solver? Then you don't have to interpolate the solutions.

By the way, you can also try desolve_odeint, which is faster and can deal with compiled code?

Good luck!

Joaquim

You can use your C/C++ code directly, using the Ctypes library, which is included in SAGE. You can find a nicely documented example on how to call your libraries on Numerical Sage documentation:

http://www.sagemath.org/doc/numerical_sage/ctypes.html

First of all you have to install the package tseriesChaos, which is a R package and not a regular Sage package.

As far as I know, to install R packages you need to have compiled version of sage. You go to Sage website, download source, follow the instructions, wait for a while (depending on your hardware) and that's all. Then type:

sage: r.install_packages('tseriesChaos')

and the package will install downloading and compiling packages automatically. The final message will be

Please restart Sage in order to use 'tseriesChaos'.

so restart Sage. After that load the library

sage: r.library('tseriesChaos')

and it should load without problems.

In my opinion, if you just want to use the 'tseriesChaos' package it might be easier to use R directly. With the steps above you, after doing "sage -R" from the shell you can do the following example from tseriesChaos package:

annaquim@aljibe:~$ sage -R

> require(tseriesChaos)
Loading required package: tseriesChaos
Loading required package: odesolve
> output <-lyap_k(lorenz.ts, m=3, d=2, s=200, t=40, ref=1700, k=2, eps=4)
Finding nearests
Keeping  1700  reference points
Following points
> plot(output)
> lyap(output, 0.73, 2.47)
(Intercept)      lambda 
  0.1185032   0.7123131

You can also consider reprogramming the routines using sage or the included numpy/scipy.

In my opinion, the funcion desolve_odeint (in versions >=4.6) provides what you want. Basically it takes a symbolic expression as the vector field, it compiles it to a fast float expression and then integrates it using the scipy.odeint solver. I have never used with hundreds of variables, but let us know your experiences...

In your case I would do the following:

sage: var('x, y')
sage: mu = 10.0
sage: dy = (y, -x + mu*y*(1-x**2))
sage: dy_fc = [expr for expr in dy]
sage: xx=srange(0,50.05,0.001) # I would plot more points
sage: ics=[1.0,0.0]
sage: sol=desolve_odeint(dy_fc,ics,xx,[x,y])

and if you want a nice plot

sage: line(sol)+plot_vector_field(dy_fc,(x,-2,2),(y,-15,15))

Hi,

Another option is to use the function desolve_odeint (from version 4.6 on, see the docs on desolve_odeint?), which uses the odeint solver internally, but takes as an input symbolic funcions . For instance, to integrate the Lorenz attractor, you would do

 sage: x,y,z=var('x,y,z') # Declare the variables 
 sage: lorenz=[10.0*(y-x),x*(28.0-z)-y,x*y-(8.0/3.0)*z]
 sage: times=srange(0,50.05,0.05) # Integration Time 
 sage: ics=[0,1,1] # and initial conditions
 sage: sol=desolve_odeint(lorenz,ics,times,[x,y,z]) #integrate

to plot the attractor, simply:

sage: line3d(sol)

I do not know what you mean exactly by the double pendulum, but here it is something similar:

sage: theta1,theta2,x1,x2=var('theta1,theta2,x1,x2')
sage: dpendulum=[x1,x2,-sin(theta1)-0.1*sin(theta2),-sin(theta2)+0.1*sin(theta1)]
sage: # Time and initial conditions
sage: times=srange(0,50.05,0.05) 
sage: ics=[1,1,0,0] 
sage: sol=desolve_odeint(dpendulum,ics,times,[theta1,theta2,x1,x2])

and a nice plot is here

sage: graphics_array([line(zip(sol[:,0],sol[:,2])),line(zip(sol[:,1],sol[:,3]))])

Joaquim Puig