[SciPy-dev] Abstract vectors in optimization

[Ben FrantzDale]

David, Robert, et al.

I think you understand what I am talking about. This is just about
optimization (although in general I would argue that functions should take
as general an interface as possible for the same reasons).

Attached is all the code I changed to make fmin_cg work with a
SimpleHilbertVector class.


If people are receptive to this approach, I have a few next questions:
1. What is the right interface to provide? Is the interface I described the
right (Pythonic) way to do it, or is there a better interface? (If it were
Templated C++, I'd do it with overloaded free functions, but I think it
needs to/should be done with class methods.)
2. Which functions should be changed? (fmin_cg, clearly; what else? The full
list in optimize.py is fmin, fmin_powell, fmin_bfgs, fmin_ncg, fmin_cg, and
fminbound.)
3. How would one deal with numpy arrays? Is it best to wrap them in the
interface (as shown in the attached code) or is it better to add methods to
them? If the right way is to wrap them, should all the optimization
functions have that as a first step -- check for numpy arrays and wrap as
needed?


Given an appropriate interface, the changes to optimize.py are essentially
search-and-replace and I would be happy to do it.

Does Travis Oliphant (whose name is on optimize.py) read this list?


―Ben

On Tue, Jan 6, 2009 at 8:28 AM, David Cournapeau <
david at ar.media.kyoto-u.ac.jp> wrote:

> Ben FrantzDale wrote:
> > David,
> >
> > Thanks for the reply. Let me describe the problem I've thought most
> > about to give you a sense of where I am coming from:
> >
> > The optimization problem I have worked most with is atomistic
> > simulation (particularly relaxation) in which you have millions of
> > atoms all interacting with a force function. The force function is
> > sparse but as the atoms move, the sparseness pattern changes (with
> > sparseness structures that can be updated automatically when the state
> > vector is "moved"). Furthermore, the atoms can be distributed via MPI,
> > so while they happen to be a size-3n array of doubles locally, local
> > atom 1,234 at one moment may not be the same atom as local atom 1,234
> > at the next step. (I don't know the details of PyMPI or how garbage
> > collection might interact with MPI, but using SciPy for
> > massively-parallel optimization is another motiviation for this post.)
> >
> > A common solver to use to relax the state of this sort of system is
> > nonlinear CG. However, I've seen many separate implementations of CG
> > hand-written to operate on these data structures. As a software
> > engineer, this is sad to see; really one nonlinear CG solver should be
> > reusable enough that it can be dropped in and used. I tried writing a
> > sufficiently-generic CG solver in C++, but designing the right family
> > of types, particularly with memory-management issues in mind, got out
> > of hand. Given the cost function is the expensive part of the sorts of
> > problems I am interested in, the overhead of Python is tiny, and SciPy
> > has nice optimization libraries, which lead me to consider
> > generalizing SciPy functions.
> >
>
> I think I simply did not understand what you meant by scipy; for me,
> scipy means the whole package - and in that case, slicing and other
> advanced features of numpy certainly are useful. But it looks like you
> had in mind only the optimize package; in that context, Robert's answer
> is certainly better than mine.
>
> >
> >
> > I agree with Robert's comment about non-Euclidean metrics, such as
> > working with SO(3). Similarly, I have thought it would be interesting
> > to optimize function spaces (as in finite-element solvers) in which it
> > should be possible to add two functions over the same domain but with
> > different discretizations, so that __add__ would implicitly remesh one
> > or the other functions rather than having __add__ require the same
> > basis functions on the LHS and RHS so it can just add corresponding
> > scalar coefficients.
>
> That can definitely be useful, and not only for optimization.
> Optimization is not a field I can claim to know a lot about, but those
> more abstract arrays with overridable 'core' aspects certainly is useful
> in other problems as well. Taking an example nearer to my own field, a
> package called lastwave for wavelets has some neat tricks related to
> indexing, including automatic interpolation which can be very useful for
> signal processing, and is a bit similar to your example of 'remeshing'
> if I understand correctly.
>
> cheers,
>
> David
> _______________________________________________
> Scipy-dev mailing list
> Scipy-dev at scipy.org
> http://projects.scipy.org/mailman/listinfo/scipy-dev
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20090106/6059896b/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: Hilbert_cg_example.py
Type: text/x-python
Size: 17465 bytes
Desc: not available
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20090106/6059896b/attachment.py>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: testit.py
Type: text/x-python
Size: 726 bytes
Desc: not available
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20090106/6059896b/attachment-0001.py>