ENH:optimize: Rewrite nnls in Python #18570
Conversation
|
Good work here! One question from me: Are there benchmarks for small matrices? |
|
Ah good one, let me run some benchmarks |
| """ | ||
| m, n = A.shape | ||
|
|
||
| AtA = A.T @ A |
There was a problem hiding this comment.
Condition number? Looks like an SVD use case.
There was a problem hiding this comment.
Didn't understand your comment. Where do you want to use SVD?
There was a problem hiding this comment.
You are constructing A.T @ A, which has the condition number a square of the condition number of A. And you're using it as a l.h.s. of a linear solve routine.
There was a problem hiding this comment.
True, albeit partial solves. The algorithm goes through submatrices of AtA to only deal with active constraints hence the A[P, P] dance.
|
For small arrays indeed there is some performance loss with this PR. The quickest way is to cythonize it obviously. But not sure if it is worth it since I don't know how much this is being used.. |
|
Also it'd be useful IMO to compare to https://github.com/scipy/scipy/blob/v1.10.1/scipy/optimize/_lsq/bvls.py, which is user-facing through |
This is impressive. For small problems there is a minor performance loss but I don't think it must be a blocker here. Cythonization can be done in a follow up. |
|
What do the optimize residents think about this rewrite? Should we go full Cython or get this in? |
|
Sorry, not familiar with this area. Id say if the python code is clean, and there's a performance benefit, then a merge should occur as is |
|
OK let's get it in. Since we still have some time for v.1.12 we will have the chance to remedy should something arise. Thanks everyone. |
Reference issue
Part of #18566, as one of the low-hanging fruits.
What does this implement/fix?
atolto control the inner iteration behavior. Since the algorithm errors out, the user needs something to relax the problem. I usedatoljust in case in the future we needrtol, so we don't need to deprecatetolkeyword.Since I am not super up-to-date with
optimizeetiquette, please let me know if I'm going against the grain.This PR now opens up the possibility to solve multiple RHS through vBenthem, Keenan, https://doi.org/10.1002/cem.889 but since there are too many things happening, I'd like to leave it to later as an enhancement. Also now code is a bit more readable (or at least familiar syntax) so probably if needed can be Cythonizable.