A useful link: scipy#6326, which details a similar ODE refactor.

I think we're almost ready to start writing. The outline is looking pretty complete.

I think we should highlight that minimize is (effectively) trying to be a class. I fleshed that out more and brought it to the top. Your most recent commit made me realize this is an issue.

I think we should follow the format of Python PEPs, https://github.com/python/peps. We can leave the outline in Markdown.

Some are literal black boxes and implemented in Fortran/C.

I don't think we should say that we're going to extract the Fortran/C code back into Python. That's really well tested over many years, and would be resisted by the community.

LBFGSB has Fortran, but it returns every step, meaning we can have an Optimizer object. If there are minimizers/solvers that run to completion (leastsq?) within C/Fortran, then I'm not sure that we could translate those easily to Optimizer.

Yes, leastsq is totally Fortran contained - little chance of making an Optimizer out of that, unless it's possible to call Fortran step by step, but that would lead to overhead (possibly unacceptable for many).

BTW I did a little bit of profiling for NelderMead. minimize rosen for starting [1, 2, 3, 4]
fmin is around 24ms
NelderMead is around 29ms.
So the class based approach is around 20% slower. Not sure where that comes from. cython might recover that.

• minimize is trying to be a class
  • method: should be subclasses

Playing devils advocate - the original solvers were fmin, f_min_bfgs, etc. Then they were collected together in minimize to introduce uniformity, which occurred to a degree. Now we're proposing that they separate again to individual classes. You can see how people can be suspicious of yet another re-write.

they were collected together in minimize to introduce uniformity, which occurred to a degree

Classes introduce another level of uniformity, correct? Subclasses can inherit from the base class, and with that comes all the parameters and functions. And overridden functions should return values of the same type/shape/etc (though there's no way to really check).

Here are issues that have arisen because of special cases:

• some optimize.minimize methods modify the parameter vector in-place scipy/scipy#8419 (comment): "some optimize.minimize methods modify the parameter vector in-place", fixed by calling callback with x.copy() in optimimze.minimize(method='L-BFGS-B', options={'disp': True}) includes "dead code" print statements scipy/scipy#8417
• scipy.optimize has broken my trust. scipy/scipy#8373 "scipy.optimize has broken my trust." which doesn't have the most civil discussion
  • From issue filer: "This has cost me dozens of hours of debugging time, only to learn that something is wrong with L-BFGS-B that causes the output of options={'disp': True} to not be the the cost function at a given parameter vector."
  • From scipy.optimize has broken my trust. scipy/scipy#8373 (comment) "It'll have to be fixed individually in each minimizer implementation, though, sometimes in Python, sometimes in C, depending on the minimizer."
• Bug in L-BFGS-B minimize? scipy/scipy#5832 looks like an issue with Fortran/Python
• minimize_scalar doesn't seem to honor "disp" option scipy/scipy#6019 "minimize_scalar doesn't honor disp option", classes could help with this
• BUG: L-BFGS-B does two more iterations than specified in maxiter scipy/scipy#7854: "BUG: L-BFGS-B does two more iterations than specified in maxiter"
• return value of scipy.optimize.minimize not consistent for 1D scipy/scipy#6673, "return value of scipy.optimize.minimize not consistent for 1D", result inconsistent for different algs
• optimize - check that maxiter is not exceeded. scipy/scipy#8375, "optimize - check that maxiter is not exceeded."

Can you give the proposal a read over? I added a section on goals and want to make sure we're on the same page.

One issue is the point

hides all details. Some are literal black boxes and implemented in Fortran/C.

I don't think we want to be pulling code out of Fortran and C. As it's written it gives the impression that we're going to do all that calculation in Python. That's an enormous task. The core of the proposal is that pretty much all optimisers have a loop, and that loop can be naturally run as an iterator . So for all the optimizers that have their loop in Python, the loop is extracted and put into __next__. For example LBFGSB has a loop, in that loop there is a call to Fortran. That Fortran call will stay as is, but the loop that calls it will now form the core of the __next__ method.

However, for something like leastsq the entire thing is in external code (minpack), calling back into python for func evaluation. At this stage we don't want to be delving into minpack and changing it. However, with the design we have it will allow someone at a future stage to do that.

It's not obvious what is meant by "black box". Replacing current loops with an iterator is not going to cure that.

I have another idea, which might sway the argument for a proportion of the community. If we put Optimizer, Function and the Optimizer subclass into Cython it might allow the use of a LowLevelCallable. I know that this provides a big speedup for things like integration where the entire operation then happens in C, it may be that it we could get big speed improvements if everything can happen in Cython, without having to call back into Python.
However, I'm not super familiar with Lowlevelcallable, it'll need some investigating.

bug reports, personal experience and anecdotal evidence.

We need hard evidence here. @jaimefrio, do you reckon you'd be able to chip in here?

What is meant by a "minimal class interface"?

Not sure that multiprocessing is immediately introducible as a result of this PR

Thanks for your comments, this is starting to get fleshed out. I've made some edits; could you take another look?

Not sure that multiprocessing is immediately introducible as a result of this PR

Maybe we should mention multiprocessing in a "possible future work" section or "possible third party work this enables"? I can point to certain async optimization algorithms that require multiprocessing. I'm fine leaving this out.

Optimizer, Function and the Optimizer subclass into Cython

You seem to be on top of this; I'll let you handle it.

bug reports, personal experience and anecdotal evidence. ... We need hard evidence here

We should also ask the managers of other libraries this question. Let's finish a draft of the report then ask them. I've added a note at the top of the draft about what to ask them about.

It's not obvious what is meant by "black box". Replacing current loops with an iterator is not going to cure that.

What if I don't want my optimizer to do a line search, or what if I want to change some of the hyperparameters in the line search? I think the methods that use a line search should have a function Optimizer._line_search. This would allow the advanced user to modify it or any hyperparameters.

What is meant by a "minimal class interface"?

One that wasn't specific to minimize or even SciPy I guess. I think I was imagining some class Optimizer with 4 subclasses, one for each of minimize, minimize_scalar, root and linprog. I'd imagine this class has a couple basic functions (__next__, results, converged, ...) that all raise a NotImplementedError.

Then, each of the methods in those functions could subclass from their parents (e.g., Optimize is a parent of Minimize is a parent of NealderMead). This could be used in solve_ivp I guess too.

I'm not saying this proposal should implement all 4 subclasses of Optimizer, but I think we should lay the groundwork for any future work.

I need to rework the outline. In the second section, we list briefly list the features of this rewrite. The rest of the paper should be expand upon each of those features.

I've been pondering a bit more. I can understand the reluctance to have 'another API', and whilst that's probably the main issue you're concerned with, I think for me by far the main attraction is that the changes will make maintenance and introducing new features a lot easier. I think we should give the latter more prominence.

One of my other thoughts, and a lesson I learnt with writing DifferentialEvolutionSolver, is that it could be a good idea to keep the Optimizer/Function classes private for a couple of releases. They would still be the core functionality for the current fmin, etc, but this approach has the advantage of:

1. being able to tweak implementation of the classes before they become public and it's too late to change anything.
2. being able to tweak implementation if we discover a problem whilst implementing other scalar minimizers.
3. Give us time to convert the other scalar minimizers before a big release.

keep the Optimizer/Function classes private for a couple of releases

I also think that's a good idea. I hadn't considered it but am on board.

I think we should give the latter more prominence.

I've moved the "other API" section to the further down.

I think we should follow the outline of an actual PEP. This is what I plan to do next.

• PEP 465, the matmul operator: https://www.python.org/dev/peps/pep-0465/
• PEP 1, Purpose and Guidelines: https://www.python.org/dev/peps/pep-0001/
• NumPy's NEP format: https://github.com/numpy/numpy/blob/master/doc/neps/nep-template.rst, though it doesn't look to be extensively reviewed

As with anything I'd like to see the evidence for that before concurring. For example, I thought Cython would speed up this PR, but there were no significant gain when I experimented with that. I hope that the real gain would be to use LowLevelCallable from Cython to reduce overhead from Python function calls.

@rgommers, I made the formatting changes you suggested, such as the header table. I think we're happy with the PEP to be placed in a scipy/* repo at this point.

@rgommers I'm also satisfied with the proposal (though it's not perfect and could use improvement).

I think SEP or SPEP are the best acronyms. Maybe SciPEP?

Thanks, looks good. I'll have a look at creating that repo now.

Okay either https://github.com/numpy/neps is still in flux, or we are keeping the sources for NEPs within the numpy repo itself after all, I can't remember. If the answer is the latter, let's put SciPEPs at doc/source/scipeps in the SciPy repo.