[SciPy-Dev] scipy.optimize.cobyla not consistant in Windows

Casper Skovby casperskovby@gmail....
Tue Jan 3 16:05:51 CST 2012


>
> Hello.


I have realized that scipy.optimize.cobyla.fmin_cobyla does not always give
the same result even though the input is exactly the same when running in
Windows. When running in Linux I do not experience this. Below I have added
an example. In this case as you can see in the example I am trying to
represent a curve with a cubic spline. As input to the spline I have 12
points – two end points and 10 inner points. I am trying to optimize these
10 inner points (both in x and y direction) in order to represent the curve
as good as possible with this spline. But as said it does not always ends
up with the same result. Sometimes you can run the code several times with
the same result but suddenly the result differ.

I realized this bug when I used OpenOpt (openopt.org) because it gave
similar problems. I have a conversation with the Developer of OpenOpt
Dmitrey (see http://forum.openopt.org/viewtopic.php?id=499). When using
OpenOpt in Linux some of the solvers do also give inconsistant results.
Dmitreys conclusion is: Since scipy_cobyla works different in Windows and
Linux, probably something with f2py is wrong (or something in libraries it
involves).

He suggested me to contact this group. Have you any ideas where the bug can
be located?

Kind Regards,
Casper

from matplotlib.pylab import *
import numpy as np
from numpy import linspace
from scipy import interpolate
from scipy.optimize import cobyla

def residual(p, r, y):
tck = interpolate.splrep(np.concatenate(([r[0]], p[0:len(p)/2], [r[-1]])),
np.concatenate(([y[0]], p[len(p)/2::], [y[-1]])))
yFit = interpolate.splev(r, tck)

resid = sum((y-yFit)**2)
  return resid

def constraint(x):
return 1


def FitDistribution(r, y):
rP = linspace(r[0], r[-1], 12)
yP = np.interp(rP, r, y)
 p0 = np.concatenate((rP[1:-1], yP[1:-1]))
        # form box-bound constraints lb <= x <= ub
lb = np.concatenate((np.zeros(len(rP)-2), np.ones(len(rP)-2)*(-np.inf))) #
lower bound
ub = np.concatenate((np.ones(len(rP)-2)*r[-1],
np.ones(len(rP)-2)*(np.inf))) # upper bound
 ftol = 10e-5
 f = lambda p: residual(p, r, y)

pvec = cobyla.fmin_cobyla(f, p0, constraint, iprint=1, maxfun=100000)
#def objective(x):
# return x[0]*x[1]
#
#def constr1(x):
# return 1 - (x[0]**2 + x[1]**2)
#
#def constr2(x):
# return x[1]
#
#x=cobyla.fmin_cobyla(objective, [0.0, 0.1], [constr1, constr2],
rhoend=1e-7, iprint=0)


tck = interpolate.splrep(np.concatenate(([r[0]], pvec[0:len(pvec)/2],
[r[-1]])), np.concatenate(([y[0]], pvec[len(pvec)/2::], [y[-1]])))
yFit = interpolate.splev(r, tck)
 return yFit

r = linspace(1,100)
y = r**2

yFit = FitDistribution(r, y)

plot(r, y, label = 'func')
plot(r, yFit, label = 'fit')
legend(loc=0)
grid()
show()
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.scipy.org/pipermail/scipy-dev/attachments/20120103/399eb317/attachment.html 


More information about the SciPy-Dev mailing list