Free Python optimization framework

Friday, August 31, 2007

ALGENCAN news

I continue trying to connect ALGENCAN to openopt.
Yestorday I informed ALGENCAN developers about minor bug in toyprob.py, and they had fix that one.

One more bug, more serious, is that they continue using obsolete Numeric module in one of files (at least, in pywrapper.c) instead of numpy, despite we agreed long tima ago they will use numpy everywhere. I informed them today and I hope they will fix it soon.
As for me I have temporary solution - replacing x by asfarray(x) in all python funcs related to ALGENCAN. So it works now, at least for UC problems, I will check constrained today. However, there are still lots of problems, and 2-nd order problems (that use d2f, d2c, d2h) are still unavailable in openopt (at all, not only for ALGENCAN).

2 more steps

2 more steps were implemented in optimizers module:
- Goldfeld
- Goldstein-Price.
Matthieu and I have some difficalties in other 2-nd order steps related to negative curvature, because available Python libraries (including numpy) hasn't modified Cholesky factorization, that is required. As we agreed no other than Python code is included to openopt, it means no fortran, etc, only python and numpy, so those fortran-written routines are inappropriate. The proposed article contains 2-pages algorithm w/o vectorization, so Python implementation will be too slow and would take much time for implementation, unittests, bugfixes.
So we decided to switch to other tasks.

Wednesday, August 22, 2007

Barzilai-Borwein algorithms

both Barzilai-Borwein algs have been written and tested.
Matthieu informed me he intends to rename some variables, especially containing word "step" because it yields a big ambiguity (alpha_step or x_step, and some more). previously direction also had "step" in variable name.

However, I have some problems with default sigma1, sigma2, alpha1, alpha2, rho, sigma values (from non-monotone Barzilai-Borwein alg). My initial guesses were too bad, as I noticed. Currently it works so-so - but maybe only for my tests.

Tuesday, August 21, 2007

UC solvers

Cubic Interpolation Search has been connected to Matthieu's code.
BarzilaiBorwein alg is written but not connected yet, I still need some more digging in Matthieu's code to connect it correctly.

Thursday, August 16, 2007

l_bfgs_b from scipy.optimize have been connected to openopt

One more box-bounded NLP solver have been connected to openopt as scipy_lbfgsb.
Is solves the example from previous post much faster: ~0.1 sec.

I tried to connect cobyla, but some problems have been encountered. The main one is: all constraints should be passed as sequence of funcs. As far as I learned, even linear constraints Ax<=b or x>=lb should be separated into funcs, each one returning single value. At least, current python-cobyla binding requires this approach. Of course, I can implement the trick, but it can slow calculations very much; maybe, some changes in Python-cobyla binding will turn out to be much more effective solution. But I have no time for now, I'm busy with other things.

One more example have been added to /examples directory: nlp1.py. It uses lincher, but any other solver, capable of handling all constraints, can be use instead.

Wednesday, August 15, 2007

tnc have been connected to openopt

(tnc is scipy.optimize NLP solver for box-bounded problems, requires scipy installed)

example of tnc usage:
----------------------------
from scikits.openopt import NLP
from numpy import *

N = 15

ff = lambda x: ((x-15)**4).sum() + 64*((x-1)**2).sum()
p = NLP(ff, cos(arange(N)))
p.lb = arange(N)
p.ub = 5 + arange(N)
r = p.solve('scipy_tnc')
----------------------------
So it should yield something like
istop: 4 (|| F[k] - F[k-1] || < funtol)
Solver: Time elapsed = 0.12 CPU Time Elapsed = 0.12
objFunValue: 10370.1398294 (feasible)

r.xf should be equal to
array([ 5. , 6. , 7. , 8. ,
8.72421634, 8.72421837, 8.72421093, 8.72420843,
8.72420744, 9. , 10. , 11. ,
12. , 13. , 14. ])


if N = 150

lincher:
itn 37 : Fk= 887272370.281 maxResidual= 5.24576648786e-09
Solver: Time elapsed = 4.16 CPU Time Elapsed = 4.01
objFunValue: 887272370.281 (feasible)

scipy_tnc:
istop: 4 (|| F[k] - F[k-1] || < funtol)
Solver: Time elapsed = 3.02 CPU Time Elapsed = 2.93
objFunValue: 887273034.135 (feasible)
maxResidual = 0

So, despite tnc has a little bit worst objFun value, it solve the example much quicker and is more suitable for box-bounded problems - of course, because it's intended for that ones and only for that ones, while lincher handles all constraints.
On the other hand, I'm not fond of CVXOPT QP solver, using better one could significantly increase lincher speed.

However, this is the only one example I had tested.
Maybe tomorrow I'll connect and publish lbfgsb and cobyla results.

----------------------------

+ today some examples to lincher have been written, I'll commit them tomorrow after minor changes.

Tuesday, August 14, 2007

news from ALGENCAN

I received another one letter from ALGENCAN developers. They informed me about some changes they made according to my suggestions (like using numpy.array instead of Python list etc).

So I don't think all is implemented very nice for now, especially the way they handle linear constraints in input data. But anyway it's already much more better and could be connected, have I enough time - now I'm a little bit out of GSoC schedule, so, maybe, I will write ALGENCAN binding in September or so. It still requires much more work than lbfgsb or tnc or cobyla (from scipy), because it has much more input parameters and non-standard (vs scipy) interface. Maybe I will write bindings to the ones before GSoC finish.

Also, ALGENCAN can use 2nd derivatives (as well as HESSmult func), but 2-nd order problems (all those d2f, d2c, d2h) are not implemented in openopt yet.

So, I sent them some more suggestions, like using iterfcn for to enable real-time openopt graphical output and user-defined stop criteria.

openopt tickets 33 and 34

I have been out of internet connection for some days and now I continue my work.
So I was informed by Nils Wagner about the ticket33 and ticket34
I didn't encountered 33 (as Nils wrote it's already had been fixed); fix to ticket 34 have been committed to svn, as well as some other changes. One of them is stand-alone line-search optimizer based on Armijo rule, so it didn't require scipy or Matthieu code. On the other hand, the latter can be very useful when no constraints are present (however, UC solvers like bfgs etc usually will be better in the case), or if the user's problem has just 1 constraint (in the case max(constraints) = the constraint = smooth func).

I tried to close the tickets but I have some problems with register/login, I hope Jeff Strunk will fix them soon.

Tuesday, August 7, 2007

line-search: still not ready

I have tried 3 line-search routines: scipy.optimize.line_search, scipy.optimize.linesearch.line_search and Matthieu's. Every one has it's own drawbacks, as I informed scipy dev mail list. So, at last, Matthieu have installed cvxopt (because lincer requires qp solver from that one; Matthieu had some problems with cvxopt compile), and begin to debug his routines. So he suspects the problem is in the algorithm he had took from a book, and he wants to implement other one (during 1-2 days).

Friday, August 3, 2007

lincher: Problems with line-search subproblem

Today I began converting openopt docstrings to new format + connecting Matthieu's line-search optimizer to lincher (instead of scipy.optimize.line_search that seems to be obsolete and no maintained any more, noone can (or want) answer about old_fval and old_old_fval parameters)
However, I encountered some difficulties, and some problems were not solved in the Matthieu's code, now he is looking at the ones.
Currently I try to use StrongWolfePowellRule, 1st iteration is performed ok but 2nd iteration returns the same point (I found that line-search solver returns same point):
>>> starting solver lincher (BSD license) with problem unnamed
itn 0: Fk= 8596.39550577 maxResidual= 804.031293334
itn 1 : Fk= 291.063506679 maxResidual= 0.686929320414
N= 672.327620905
itn 2 : Fk= 291.063506679 maxResidual= 0.686929320414
N= 2040.4800449
solver lincher has finished solving the problem unnamed
istop: 4 (|| F[k] - F[k-1] || < elapsed =" 0.53" elapsed =" 0.5" residual =" 0.686929320414)">>> lsF(x0)
1692.7290772406229
>>> lsF(x0+1e-4*direction)
1691.0958067992992#grow down
>>> lsF(x0+1e-3*direction)
1681.525182413252
>>> lsF(x0+1e-2*direction)
1681.5725957967709# starting to grow up again
>>> lsF(x0+1e-1*direction)
1684.4384883880082

I informed Matthieu and I hope to obtain an answer tomorrow.

Wednesday, August 1, 2007

more info on ticket 464

some email messages to/from Nils Wagner + some hours were elapsed for to found one more bug related to ticket 464.
So it turns out that asfarray(matrix([[0.3]])) returns array in numpy 1.0.1 and matrix in numpy 1.0.4.dev3937.
I think 1st way (array) is more correct.
Nils have promised to organize a ticket.

tickets: ready to be closed

all tickets assigned to me (+some more from Nils Wagner) now are ready to be closed.
The only one remaining is 399 (connecting PSwarm), but maybe it will be assigned to other person, because it requires connecting Python and C.
Now I'm trying to connect some Matthieu's solvers to openopt, for my own usage (in solvers lincher and ralg) as well.