Modular Optimization and Root Finding

This page contains discussions about a modularized optimization and root finding package. Please make copious changes and add lots of comments and discussions. We should refactor once it becomes large. Once ideas are clearly formed, they should be posted to the SciPy-dev mailing list for comments and suggestions.

Contents

(Why does this macro not work to generate the table of contents?)

Module Structure

For now I am putting everything in an moptimize module for "Modular Optimization". Hopefully this will be able to replace the scipy.optimize module at somepoint. Note that coding standards, etc. can be found at the following sites:

• PEP 8 -- Generic python coding standards.
• TestingGuidelines -- NumPy/SciPy testing guidelines.
• SciPy/DocstringStandard -- NumPy/SciPy docstring standards.

moptimize 
    |-- __init__.py : Initialization of package (each subdir should have one)
    |
    |-- tests : Testing (each subdir should have one.)
    |
    |-- function : Function wrappers, jacobian computations etc.
    |
    |-- methods : General optimization methods
    |
   ...

Interface

There should be a very simple interface with default options and methods specified to get code working. One should then be able to refine the method to suit the particular problem. Ultimately, it might be nice to provide tools that could try different methods and parameter values to determine the "optimum" solver for a given problem.

Here is an example of a paradigm where a class is used to return a new functor which returns the zeros of a function

from moptimize import Zeros

def f(x,n):
    """Return x**2 - n.  Zero is at sqrt(n)."""
    return x*x - n

zeros = Zeros(f=f,x0=0)

sqrt4 = zeros(4.0)

The idea here is that the Zeros class will perform any initializations (that could be somewhat costly) and checks, then return an optimized function (functor) zeros that actually performs the calculation. The argument f could be a simple function or a more complex functor with instructions, for example, about how to compute the Jacobian. The initializer for the Zeros class would query f to see its capabilities.

If you want to use a specific algorithm, then you could specify it at this point. For example:

from moptimize.methods import Broyden
from moptimize import Zeros
...
method = Broyden(ftol=1e-10,xtol=1e-12)
zeros = Zeros(f=f,x0=0,method=method)

Broyden would be a convenience wrapper of a modularized optimization method specifying the step type, line search method, etc. One could also build ones own custom method using the components. The Broyden class should be able to provide a list of options that can easily be selected to customize the method (either by the user or by an automatic performance tuner) as a starting point for a custom optimizer for example. One could subclass:

from moptimize.method import Broyden
from moptimize.method.criteria import RelativeValueCriterion

class MyMethod(Broyden):
    criterion = RelativeValueCriterion

method = MyMethod

or one could build from options:

from moptimize.method import Broyden
from moptimize.method.criteria import RelativeValueCriterion
method = Broyden(criterion=RelativeValueCriterion())

It would also be nice if there was a unified method for passing keyword arguments in at various stages, so one could change individual paremeters at any stage, i.e. the tolerance could be specified when constructing the method, or when calling the actual optimizer:

zeros = Zeros(f=f,x0=x0,ftol=0.001)
approx2 = zeros(4.0)
better2 = zeros(4.0,ftol=1e-16)

Calling Conventions

It would be nice if the calling conventions for each method were as compatible as possible so that it is easy to swap different methods in an out.

Arguments

f or func

Every driver routine will need to specify a function. Should this argument be called f or func. In SciPy, both are used about the same amount.

x0

When an initial argument is supplied, the convention seems to be to use x0.

ftol, xtol, etc.

These seem to be the convention for specifying tolerances. The only problem I have with this is that I would like the option of specifying both absolute and relative tolerances. Perhaps something like xtol and xtol_rel could be used (absolute tolerances are usually easier to deal with in general).

Internal Structure

As per some arguments expressed on the mailing list, the internal modules will expect vectors as row vectors. Thus, for a root-finder, the current step x and result f(x) will be row vectors, and the Jacobian will be structured such that f(x+dx) = f(x) + dx*J + ... where * represents matrix multiplication.

Performance Issues

The modular approach, with great flexibility in choosing and tuning the optimizer works well for extremely costly functions where the overhead of computing the function is the dominant cost and the goal of the optimizer is to minimize the number of function calls. It may not be so suitable for functions that are fast: the overhead of python may become significant and the user might be better off using a routine that is coded directly in C or Fortran (for example, if the optimization occurs in the core of a loop). Hard-coded optimizers, such as already exist in SciPy should be able to be used with a very similar syntax, however, they will not benefit from the modularity. This is something to keep in mind.