[SciPy-dev] SciPy improvements

Ondrej Certik ondrej@certik...
Wed Apr 11 12:02:42 CDT 2007


I am studying theoretical physics and I have collected a lot of useful
python code, that I believe could go to SciPy. So first I'll review
what I have and if you find it interesting, I would like to discuss
how I could implement it in SciPy.

1) Optimization


I have a Differential Evolution optimizer, Simplex optimizer, mcmc
(not well tested yet), I took a code from someone else, but adapted
the interface to the SciPy's one:

def fmin_de(f,x0,callback=None,iter=None):

Those are unconstrained optimizers. Then I have a constrains code,
that applies a logistic function to the fitting variable and allows me
to do constrained optimization. For example the L-BFGS with my
constrains converges 7x faster, than the original L-BFGS-B on my

2) Nonlinear solvers

I have written these nonlinear solvers for the problem R(x) = 0, where
x and R has a dimension "n":

    broyden1 - Broyden's first method - is a quasi-Newton-Raphson method for
        updating an approximate Jacobian and then inverting it
    broyden2 - Broyden's second method - the same as broyden1, but updates the
        inverse Jacobian directly
    broyden3 - Broyden's second method - the same as broyden2, but instead of
        directly computing the inverse Jacobian, it remembers how to construct
        it using vectors, and when computing inv(J)*F, it uses those vectors to
        compute this product, thus avoding the expensive NxN matrix
    broyden_generalized - Generalized Broyden's method, the same as broyden2,
        but instead of approximating the full NxN Jacobian, it construct it at
        every iteration in a way that avoids the NxN matrix multiplication.
        This is not as precise as broyden3.
    anderson - extended Anderson method, the same as the broyden_generalized,
        but added w_0^2*I to before taking inversion to improve the stability
    anderson2 - the Anderson method, the same as anderson, but formulated

I use them in the self-consistent cycle of the Density Functional
Theory (so I use a terminology of DFT literature in the names of the

Also I am writing a BFGS solver with linesearch, that should behave
even better than the Broyden scheme.

Of course I am trying to use SciPy's code (like linesearch) wherever possible.

3) PETSC bindings

I found these nice petsc bindings:


I believe this could also be an optional package in SciPy. Because if
SciPy has some sparse matrices code, then it should definitely also
has this.

4) Finite element code?

I have my own code, that uses libmesh:


and calls tetgen and parses input from gmsh etc. Can convert the mesh,
can refine it, can solve it, etc. Webpages are here:


I am not sure here, if it should belong to SciPy. Probably not.

5) Symbolic manipulation in Python


We'll have some google summer of code students working on this and
also I am not sure if it belongs to SciPy. However, this project looks
very promising.

So that's it.

I have some comments to SciPy:

1) Documentation

Virtually none, I just use the source code to understand, what SciPy
can do and how. But the docstrings are good though. I would suggest to
update the


more often (for example I didn't find there the new l-bfgs code).

2) What is the official NumPy's page?

I believe it should be


however, it points to a sourceforge page. I believe this homepage
should contain all the relevant information about numpy and contain
links to the fee based documentation and possibly some tutorials.
The documentation of NumPy is scattered across many pages and I find
it confusing.

I know you have some list here:


But I am quite confused from the whole SciPy page. I think less but
unconfusing information is better, but that's just my opinion. I think
the front page of both SciPy and NumPy should be clean and simple with
a clear link to a documentation. Like


However, the documentation:


is confusing, because except the fee based guide to NumPy, it's not
clear, what is the official SciPy doc and what the best way of
learning SciPy is.

So I am interested in your opinions and then I'll just integrate my
code into scipy.optimize and send you a patch to this list?


More information about the Scipy-dev mailing list