GEKKO NL solver failed with initial guess of 0

Viewed 32

I'm new to GEKKO and test it with a simple non linear optimization problem, but it fires an error exception of "Solution Not Found" when initial guess is 0. I have two issues with it:

(1) Why ?

(2) How can I know which initial guess is ok for the solver?

My code was written according to the Nonlinear Regression example given here: https://github.com/BYU-PRISM/GEKKO/blob/master/docs/examples.rst

Measure points file can be downloaded from here: https://www.dropbox.com/s/ai5vbd2u5gx8r17/measure.npy?dl=0

# Model:
#            x
# y = -------------
#      1- (a1*x)^2
#
# Finding a1 for fitting model as much as possible for the given measures points


# loading given measure points (file exist at https://www.dropbox.com/s/ai5vbd2u5gx8r17/measure.npy?dl=0)
measure_data = np.load("C:/measure.npy")
xm = measure_data[0]
ym = measure_data[1]

# GEKKO model
m = GEKKO()

# parameters
x = m.Param(value=xm)
a = [m.FV(value=0.000) for i in range(1)]
for par in a:
    par.STATUS=1

# variables
y = m.CV(value=ym)
y.FSTATUS=1

m.Equation(y==x/(1-(a[0]*x)**2))
    
# regression mode
m.options.IMODE = 2

# optimize
m.solve(disp=False)


p = [par.value[0] for par in a] 
optimized_y = xm/(1-(p[0]*xm)**2)

plt.figure(1)

plt.plot(xm,ym,'k', label = "measurements")
plt.plot(xm,optimized_y,'r', label = "optimized_y")

plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.show()

The solver gives a solution if I set the initial guess to 0.001, but it was just playing with value with no real understanding what i'm doing.

enter image description here

1 Answers

An initial guess of 0 shouldn't be a problem, but it depends on how the solver searches for the solution and the data. If the solver iterates to a[0]*x==1 then it could result in divide-by-zero.

m.Equation(y==x/(1-(a[0]*x)**2))

Try reformulating the equation as:

m.Equation(y*(1-(a[0]*x)**2)==x)

If the equation can't be rearranged to avoid the singularity, another thing to try is setting a constraint z>=1e-5 to avoid divide-by zero such as:

z = m.Var(lb=1e-5)
m.Equation(z==(1-(a[0]*x)**2))
m.Equation(y==x/z)

For future questions, please post a complete problem with sample that reproduces the issue and can be tested. This question is missing measure.npy.

Related