Numerical Optimization Forum - Support

Re: problem with interalg

dummy@example.com (cytochrome) — Mon, 15 Aug 2011 20:03:46 +0000

Dmitrey wrote:

hi,
I've fixed it 2 minutes later after committing constraints handling for interalg:
http://trac.openopt.org/openopt/changes … w_path=%2F
(I forgot to commit another one file baseProblem.py)
try updating from svn repository

Hi Dmitrey,

Thank you for your amazing responsiveness!

-c

Re: problem with interalg

dummy@example.com (Dmitrey) — Mon, 15 Aug 2011 20:00:12 +0000

hi,
I've fixed it 2 minutes later after committing constraints handling for interalg:
http://trac.openopt.org/openopt/changes … w_path=%2F
(I forgot to commit another one file baseProblem.py)

try updating from svn repository

problem with interalg

dummy@example.com (cytochrome) — Mon, 15 Aug 2011 19:54:48 +0000

Hi,

I am experiencing an error trying to use the interalg functionality.

System:

OpenOpt 0.34
Windows 7
Enthought Python 7.1-1

I ran exactGlobalNLP.py and got the following traceback:

------------------------- OpenOpt 0.34 -------------------------
solver: interalg_0.21   problem: unnamed    type: GLP
 iter   objFunVal   
    0  1.163e+02 
Traceback (most recent call last):
  File "C:\Users\breisfel\Documents\Programming\python\exactGlobalNLP.py", line 34, in 
    r = p.minimize(solver)
  File "C:\Python27\lib\site-packages\openopt-0.34-py2.7.egg\openopt\kernel\baseProblem.py", line 196, in 
    self.minimize = lambda *args, **kwargs: minimize(self, *args, **kwargs)
  File "C:\Python27\lib\site-packages\openopt-0.34-py2.7.egg\openopt\kernel\baseProblem.py", line 774, in minimize
    return runProbSolver(p, *args, **kwargs)
  File "C:\Python27\lib\site-packages\openopt-0.34-py2.7.egg\openopt\kernel\runProbSolver.py", line 238, in runProbSolver
    solver(p)
  File "C:\Python27\lib\site-packages\openopt-0.34-py2.7.egg\openopt\solvers\UkrOpt\interalg_oo.py", line 225, in __solver__
    C = p._FD.nonBoxConsWithTolShift
AttributeError: EmptyClass instance has no attribute 'nonBoxConsWithTolShift'

Re: OpenOpt for non linear problems

dummy@example.com (Dmitrey) — Fri, 12 Aug 2011 09:05:37 +0000

however, using stats.norm.pdf as exp(...) has the following serious drawback for interalg: it has several occurences of argument and thus decreases quality of interval analysis, hence interalg may work essentially worse.

Re: OpenOpt for non linear problems

dummy@example.com (Dmitrey) — Fri, 12 Aug 2011 07:30:25 +0000

gbenoit, your knowledge of python language seems to be absolutely inappropriate for solving the problems like this one.
In my last message I've already provided all the info you're need to solve it.
If you intend to use openopt w/o funcdesigner, you should use
def Psi(x):
alpha, a0, a1, b1 = x
return math.log(stats.norm.pdf(ThetaSto[j],...)
and x0 = [0.2,2,0.2,0.1]

All Python-scientific folk use numpy.exp, numpy.sin etc instead of math.exp, math.sin etc.

In correct form FD code should look like this (you should change values of j, ThetaSto and maybe chose more appropriate implicitBounds):

from numpy import *
from openopt import *
from FuncDesigner import *

alpha, a0, a1, b1 = oovars('alpha a0 a1 b1')

ThetaSto = array([1,2,3]) 
j = 1    

Psi = (log(sqrt(2*pi)*(a0 + a1*cos(j*2*pi/365) + b1*sin(j*2*pi/365))*sqrt((1-exp(-2*alpha))/(2*alpha))\
                                     )-ThetaSto[j]**2*(1-exp(-alpha))**2/(2*((a0 + a1*cos(j*2*pi/365) + b1*sin(j*2*pi/365))*\
                                     sqrt((1-exp(-2*alpha))/(2*alpha)))**2))
                                     
# parameters estimation
StartPoint = {alpha:0.2, a0:2, a1:0.2, b1:0.1}

p = NLP(Psi, StartPoint, goal='max') 

#r = p.solve('ralg') - I'm not sure it will work with all those multiextrema funcs like cos, sin, nonconvex log
r = p.solve('interalg', implicitBounds = (-1000, 1000), fTol = 1e-10)

print(r(alpha, a0, a1, b1))

however, you have sqrt((1-exp(-2*alpha))/(2*alpha)), that is not defined over whole R^nVars, and with interalg you may encounter something like

unhandled ValueError
"numpy.nanargmin raises on a.shape[axis]==0; Bottleneck too."
File: nanargmin.pyx, Line: 495

(or another message if you don't have bottleneck installed). Bugfix for this will be released in few days along with new constrained interalg version. As for some other solvers like ralg or de, you should add constraint (1-exp(-2*alpha))/(2*alpha) > 0.

Re: OpenOpt for non linear problems

dummy@example.com (gbenoit) — Fri, 12 Aug 2011 06:47:30 +0000

Dmitrey wrote:

So what about other issues? What do you mean with log(stats.norm.pdf(3 parameters)).sum()?
If it was log(stats.norm.pdf(parameter_1)) + log(stats.norm.pdf(parameter_3)) + log(stats.norm.pdf(parameter_3)), here's my version by FD (giveas auto derivatives and interval analysis for interalg):
from numpy import *
from openopt import *
from FuncDesigner import *

alpha, a0, a1, b1 = oovars('alpha a0 a1 b1')

ThetaSto = array([1,2,3])
j = 1    # actually ThetaSto is an existing array correctly defined (12775x1)

#pdf = lambda x: exp(-x**2/2)/sqrt(2*pi)

logpdf = lambda x: -x**2 / 2.0 - 0.5*log(2*pi)
# however, that fixed constant 0.5*log(2*pi) doesn't matter for optimal point searching


Psi = logpdf(ThetaSto[j]) + logpdf(ThetaSto[j]*exp(-alpha)) + logpdf((a0 + a1*cos(2*pi*j/365) + b1*sin(2*pi*j/365))*\
                                      sqrt((1-exp(-2*alpha))/(2*alpha)))

# parameters estimation
StartPoint = {alpha:0.2, a0:2, a1:0.2, b1:0.1}

p = NLP(Psi, StartPoint, goal='max') 

#r = p.solve('ralg') - I'm not sure it will work with all those multiextrema funcs like cos, sin, nonconvex log
r = p.solve('interalg', implicitBounds = (-1000, 1000), fTol = 1e-10)

print(r(alpha, a0, a1, b1))
#[500.0, 0.0, 0.0, 0.0]
I intend to add some stats funcs like erf, pdf, but I have some more urgent things for now.
In our country we don't have neither internal study nor internship, could you explain what openopt is used for (which task?)
Regards, D.

Thank you for your answer. I have been away for a few days that's why I did not answer your last message.
I have "simplified" the formula (applying the log on the previous one so I don't need stats.norm.pdf anymore) :

Psi = lambda alpha, a0, a1, b1: (-math.log(math.sqrt(2*Pi)*(a0 + a1*math.cos(j*2*pi/365) + b1*math.sin(j*2*pi/365))*math.sqrt((1-math.exp(-2*alpha))/(2*alpha))\
)-ThetaSto[j]**2*(1-math.exp(-alpha))**2/(2*((a0 + a1*math.cos(j*2*pi/365) + b1*math.sin(j*2*pi/365))*\
math.sqrt((1-math.exp(-2*alpha))/(2*alpha)))**2)).sum()

I have been using oovars and 'interalg' instead of 'ralg' but I still get the same error :

File "C:\Python26\lib\site-packages\openopt-0.34-py2.6.egg\openopt\kernel\ooMisc.py", line 202, in setNonLinFuncsNumber
setattr(p, 'n'+userFunctionType, asfarray(FV(*(X, ) + args, **kwargs)).size)
TypeError: () takes exactly 4 arguments (1 given)

Concerning the kind of study carried out, this is a study on a temperature model.

Regards,

Re: Problem with problem, galileo and de

dummy@example.com (Dmitrey) — Thu, 11 Aug 2011 05:23:26 +0000

This is not mess. oovars yields ooarray, that is derived class of numpy.ndarray and has all its attributes. All solvers except of interalg can handle both oovars() and oovar(), and interalg cannot due to interval analysis, that is currently performed for single-sized oovar only (or ooarray of them). Maybe in long-term future it will be implemented for oovar with size > 1.

Re: Problem with problem, galileo and de

dummy@example.com (baryluk) — Wed, 10 Aug 2011 20:31:43 +0000

Thanks.

Are there some plans to consolidate and unify oovars / oovar mess?

Re: Problem with problem, galileo and de

dummy@example.com (Dmitrey) — Tue, 09 Aug 2011 22:26:04 +0000

You could use
Theta_DNA_YY1<[1,2,3], Theta_DNA_YY1> [-1,-20,-300]
(or numpy arrays).
If you use oovars(n) (for interalg) instead of oovar(size=n), then you could use

Theta_DNA_YY1[i] < ub[i], Theta_DNA_YY1[i] > lb[i]

(however, those Theta_DNA_YY1<[1,2,3], Theta_DNA_YY1> [-1,-20,-300] should work as well)

Re: Problem with problem, galileo and de

dummy@example.com (baryluk) — Tue, 09 Aug 2011 22:06:42 +0000

What if I want different bound on components of vector Theta_DNA_YY1. I want to use interalg, so would like stay in given bounds, and some of them will be different in single vector. Should I manually rescale problem by multipling Theta_DNA_YY1 by some numerical constants?

Thanks, I changed code as you suggested and it works. I used vector form previously but in actuall full code I have different constrains on components, and also first (vector) form had some other problems previously (I do not remember what problems).

Re: Problem with problem, galileo and de

dummy@example.com (Dmitrey) — Tue, 09 Aug 2011 20:49:10 +0000

Correct:

1) this code should be removed:
for i in xrange(Zmiennych):
constraints.extend([-1.0e6 <= Theta_DNA_YY1_i, Theta_DNA_YY1_i <= 1.0e6])

2) use
constraints = [0 < K, K < 1000,Theta_DNA_YY1>-1e6, Theta_DNA_YY1 < 1e6]
instead. Also, you may use implicitBounds = (-1e6, 1e6) or mere implicitBounds = 1e6 instead, it doesn't makes problem constrained for those solvers that don't require finite box bounds.

Problem with problem, galileo and de

dummy@example.com (baryluk) — Tue, 09 Aug 2011 20:40:44 +0000

Hi, I'm trying to use openopt, but find lots of bugs, and when I resolve one, or find some kind of workaround, then next I just find next one.

Here is simple model, in which I try to fit some parameters to experimental data, essentially using least-square model + additional parameter, making it non-linear. It is here presented as one in 3 variables (original have about Zmiennych=200 + 1 for K), but error is the same.

from openopt import GLP, oosolver
from FuncDesigner import *

from numpy import zeros, ones

#K_ = oovars(1)
#K = K_[0]
K = oovar(size=1)
A0 = 2.0e-6
B0 = 5.0e-6

Zmiennych = 2

Theta_DNA = zeros(Zmiennych) # , int/float)
Theta_YY1 = zeros(Zmiennych)
Theta_Total = zeros(Zmiennych)
Theta_DNA_YY1 = oovar(size=Zmiennych)

for i in range(Zmiennych):
    Theta_YY1[i] = 4.5
    Theta_DNA[i] = 8.7
    Theta_Total[i] = 4.1


S_AB = 0.5*(A0 + B0 + K*1.0e-6 - sqrt( (A0 + B0 + K*1.0e-6)**2 - 4.0*A0*B0 ) )
S_A = (A0 - S_AB)
S_B = (B0 - S_AB)

Roznice = 1.0e-3*(S_A*Theta_DNA + S_B*Theta_YY1 + S_AB*Theta_DNA_YY1) - Theta_Total

f = Roznice[0]**2 + Roznice[1]**2
#f = sum(Roznice**2)

constraints = [0 < K, K < 1000]     # Was trying to use <=, but it doesn't help
for i in xrange(Zmiennych):
    constraints.extend([-1.0e6 <= Theta_DNA_YY1[i], Theta_DNA_YY1[i] <= 1.0e6])

startPoint = {K: 0.32, Theta_DNA_YY1: ones(Zmiennych)}

p = GLP(f, startPoint, constraints=constraints)

#solver = oosolver('de')
solver = oosolver('galileo')

r = p.solve(solver, plot=1)

Using de, I got

FuncDesigner warning: Probably scipy installation could speed up running the code involved

------------------------- OpenOpt 0.34 -------------------------
solver: de   problem: unnamed    type: GLP
 iter   objFunVal   log10(maxResidual)   
    0  3.362e+01            -100.00 
OpenOpt Error: this solver requires finite lb, ub: lb <= x <= ub
Traceback (most recent call last):
  File "./bug2.py", line 47, in 
    r = p.solve(solver, plot=1)
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/kernel/baseProblem.py", line 235, in solve
    return runProbSolver(self, *args, **kwargs)
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/kernel/runProbSolver.py", line 238, in runProbSolver
    solver(p)
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/solvers/UkrOpt/de_oo.py", line 74, in __solver__
    if not p.__isFiniteBoxBounded__(): p.err('this solver requires finite lb, ub: lb <= x <= ub')
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/kernel/oologfcn.py", line 16, in ooerr
    raise OpenOptException(msg)
oologfcn.OpenOptException: this solver requires finite lb, ub: lb <= x <= ub

and one using galileo

FuncDesigner warning: Probably scipy installation could speed up running the code involved

------------------------- OpenOpt 0.34 -------------------------
solver: galileo   problem: unnamed    type: GLP
OpenOpt Error: the solver galileo cannot handle 'A' data
Traceback (most recent call last):
  File "./bug2.py", line 47, in 
    r = p.solve(solver, plot=1)
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/kernel/baseProblem.py", line 235, in solve
    return runProbSolver(self, *args, **kwargs)
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/kernel/runProbSolver.py", line 219, in runProbSolver
    nErr = check(p)
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/kernel/check.py", line 19, in check
    p.err('the solver ' + p.solver.__name__ + ' cannot handle ' + "'" + fn + "' data")
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/kernel/oologfcn.py", line 16, in ooerr
    raise OpenOptException(msg)
oologfcn.OpenOptException: the solver galileo cannot handle 'A' data

With interalg,

FuncDesigner warning: Probably scipy installation could speed up running the code involved

------------------------- OpenOpt 0.34 -------------------------
solver: interalg_0.21   problem: unnamed    type: GLP
 iter   objFunVal   log10(maxResidual)   
    0  3.362e+01            -100.00 
OpenOpt Error: 
            solver interalg_0.21 requires finite lb, ub: 
            lb <= x <= ub 
            (you can use "implicitBoounds")
            
Traceback (most recent call last):
  File "./bug2.py", line 48, in 
    r = p.solve(solver, plot=1)
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/kernel/baseProblem.py", line 235, in solve
    return runProbSolver(self, *args, **kwargs)
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/kernel/runProbSolver.py", line 238, in runProbSolver
    solver(p)
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/solvers/UkrOpt/interalg_oo.py", line 46, in __solver__
    ''' % self.__name__)
  File "/usr/local/lib/python2.6/dist-packages/openopt-0.34-py2.6.egg/openopt/kernel/oologfcn.py", line 16, in ooerr
    raise OpenOptException(msg)
oologfcn.OpenOptException: 
            solver interalg_0.21 requires finite lb, ub: 
            lb <= x <= ub 
            (you can use "implicitBoounds")

I really do not understand what is going on in both cases. Name (even automatically created) for which data, variable is a problem would be helpful. In first case, I given all constrains, being finite. In second case, I have no idea what solver is asking me about 'A' data.

I'm fighting with this for more than a week, without success. Was trying to use oovar or oovars, in multiple combinations, but then other problems arises. Help will be appreciated. Thanks.

Re: OpenOpt for non linear problems

dummy@example.com (Dmitrey) — Thu, 04 Aug 2011 16:04:56 +0000

So what about other issues? What do you mean with log(stats.norm.pdf(3 parameters)).sum()?

If it was log(stats.norm.pdf(parameter_1)) + log(stats.norm.pdf(parameter_3)) + log(stats.norm.pdf(parameter_3)), here's my version by FD (giveas auto derivatives and interval analysis for interalg):

from numpy import *
from openopt import *
from FuncDesigner import *

alpha, a0, a1, b1 = oovars('alpha a0 a1 b1')

ThetaSto = array([1,2,3])
j = 1    # actually ThetaSto is an existing array correctly defined (12775x1)

#pdf = lambda x: exp(-x**2/2)/sqrt(2*pi)

logpdf = lambda x: -x**2 / 2.0 - 0.5*log(2*pi)
# however, that fixed constant 0.5*log(2*pi) doesn't matter for optimal point searching


Psi = logpdf(ThetaSto[j]) + logpdf(ThetaSto[j]*exp(-alpha)) + logpdf((a0 + a1*cos(2*pi*j/365) + b1*sin(2*pi*j/365))*\
                                      sqrt((1-exp(-2*alpha))/(2*alpha)))

# parameters estimation
StartPoint = {alpha:0.2, a0:2, a1:0.2, b1:0.1}

p = NLP(Psi, StartPoint, goal='max') 

#r = p.solve('ralg') - I'm not sure it will work with all those multiextrema funcs like cos, sin, nonconvex log
r = p.solve('interalg', implicitBounds = (-1000, 1000), fTol = 1e-10)

print(r(alpha, a0, a1, b1))
#[500.0, 0.0, 0.0, 0.0]

I intend to add some stats funcs like erf, pdf, but I have some more urgent things for now.
In our country we don't have neither internal study nor internship, could you explain what openopt is used for (which task?)

Regards, D.

Re: OpenOpt for non linear problems

dummy@example.com (gbenoit) — Thu, 04 Aug 2011 15:47:32 +0000

Dmitrey wrote:

I see something like some bugs in your code.
stats.norm.pdf(param1, param2, param3, param4) seems to work incorrectly, maybe you meant stats.norm.pdf([param1, param2, param3, param4]) instead?
--> I tried in the shell and at least this part seems correct.
--> I don't want the probability for a N(0,1) but the probability for a N(m, sigma) : my second parameter is the mean (m), my third parameter the standard deviation (sigma)

If you want to maximize Psi function, you should use p = NLP(Psi, StartPoint, goal = 'max'), elseware minimization will be used (by default).
There are other issues I'll comment after your answer.
What is EDF Research and Development (url) and what is "internal studies"? Education?

--> EDF = Electricité de France (french electricity supplier) url: www.edf.com
--> internal study = internship

--> thank you anyway for your quick answer

Re: OpenOpt for non linear problems

dummy@example.com (Dmitrey) — Thu, 04 Aug 2011 15:10:27 +0000

I see something like some bugs in your code.

stats.norm.pdf(param1, param2, param3) seems to work incorrectly, maybe you meant stats.norm.pdf([param1, param2, param3]) instead?

If you want to maximize Psi function, you should use p = NLP(Psi, StartPoint, goal = 'max'), elseware minimization will be used (by default).

There are other issues I'll comment after your answer.

What is EDF Research and Development (url) and what is "internal studies"? Education?