# [SciPy-dev] SciPy improvements

Robert Kern robert.kern@gmail....
Wed Apr 11 15:09:13 CDT 2007

```Ondrej Certik wrote:
> Hi,
>
> 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.

Excellent! Thank you! One thing to be aware of is that scipy uses the BSD
permission from the others who have contributed to the code you are submitting.

> 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):

Well, fmin() is already an implementation of the simplex algorithm. How does
yours compare? We can't include the MCMC optimizer until we have an
implementation of Metropolis-Hastings in scipy itself; we're not going to depend
on an external PyMC. As for the differential evolution code, with all respect to
you and Jame Phillips, it's not the best way to implement that algorithm in
Python. It's a straight translation of the C++ code so it doesn't make use of
numpy at all. I have an implementation that does:

http://svn.scipy.org/svn/scipy/trunk/Lib/sandbox/rkern/diffev.py

It was written for pre-numpy scipy, so it may need some sprucing-up before it works.

> 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
> problem.

Interesting. Let's toss it up on the Cookbook first and pound on it a bit. I
have qualms about applying such transformations to the domains of target
functions and then using derivative-based optimizers on them, but those qualms
might be baseless. Still, I'd rather experiment first before putting them into
scipy.

> 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
>         multiplication.
>     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
>         differently
>     linear_mixing
>     exciting_mixing
>
> 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
> methods).
>
> 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.

That's fantastic. I'd love to see them. Are they in chemev? I don't see them.

> 3) PETSC bindings
>
> I found these nice petsc bindings:
>
> http://cheeseshop.python.org/pypi/petsc4py/0.7.2
>
> 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.

I don't know. It has a fine (and probably better) existence separate from scipy.

> 4) Finite element code?
>
> I have my own code, that uses libmesh:
>
> http://libmesh.sourceforge.net/
>
> 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.

I think you are right. You can't really get around the licenses of your
dependencies, here.

> 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.

Again, I think it has a fine existence separate from scipy. A reason to bring it
into scipy would be such that other parts of scipy would use it to implement
their own stuff. Otherwise, I don't think there is much point.

> -----------------
> 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
>
> http://www.scipy.org/doc/api_docs/
>
> more often (for example I didn't find there the new l-bfgs code).
>
> 2) What is the official NumPy's page?

http://numpy.scipy.org

> I believe it should be
>
> http://numpy.org/
>
> however, it points to a sourceforge page.

Correct. I don't know who still owns that domain.

> 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:
>
> http://www.scipy.org/MigratingFromPlone
>
> 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
>
> http://www.pytables.org/moin/PyTables
> http://matplotlib.sourceforge.net/

page. Because the front page is special, I'd recommend submitting your
modifications as a ticket to our Trac (see below) instead of just editing it.
The other pages on the wiki, please modify them as you see fit.

> However, the documentation:
>
> http://www.scipy.org/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.

There really is no official scipy doc at this time. That's part of the problem.

> 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?

Register an account with the scipy Trac (click "Register" in the upper-right
corner):

http://projects.scipy.org/scipy/scipy

Then make a new ticket and attach your patch to that. Submit enough patches, and
we'll just give you SVN access.

--
Robert Kern

"I have come to believe that the whole world is an enigma, a harmless enigma