Wrapping it up

August 7, 2008

The final deadline of GSoC is approaching quickly. Fortunately, I’m almost finished with my project; the main part of the code was pushed in a few days ago, and in the coming days I only expect to fix bugs and write documentation.

Here are some examples of things that should now work with the hg version of SymPy:


>>> from sympy import *
>>> var('n x')
(n, x)
>>> (log(1+I)/pi).evalf(50)
0.11031780007632579669822821605899884549134487436483 + 0.25*I
>>> (atan(exp(1000))-pi/2).evalf(15, maxprec=1000)
-5.07595889754946e-435
>>> Integral(exp(-x**2), (x, -oo, oo)).evalf(50)
1.7724538509055160272981674833411451827975494561224
>>> Integral(sin(1/x), (x, 0, 1)).transform(x, 1/x).evalf(25, quad='osc')
0.5040670619069283719898561
>>> Sum(1/(9*n**2+4*n+6), (n, 0, oo)).evalf(50)
0.27865709939940301230611638967893239422996461876558

Besides standard formula evaluation, numerical integration is supported, as is summation of infinite hypergeometric series.

Numerical differentiation and summation of generic infinite series is also on its way. The following works with my local copy:


>>> Sum(1/n-log(1+1/n), (n, 1, oo)).evalf(50)
0.57721566490153286060651209008240243104215933593992
>>> Derivative(-gamma(x), x).evalf(50, subs={x:1})
0.57721566490153286060651209008240243104215933593992

The series is summed using the Euler-Maclaurin formula; in fact, the high-order derivatives are computed symbolically behind the scenes. (Because of the need for high-order derivatives, the Euler-Maclaurin formula makes much more sense in the context of a CAS than purely numerically.) This algorithm does have its drawbacks, so at some points pure extrapolation methods (which are already used for slowly convergent hypergeometric series) will be supported as well.

The use of an adaptive algorithm for evalf turns out to be especially useful for numerical differentiation, as one can just directly use a finite difference formula with a tiny step size, without worrying about cancellation error.

There is a lot more to be said about the details, but I’m leaving that for the to-be-written documentation.

fredrikj.net | Blog index | RSS feed | Follow me on Mastodon | Become a sponsor