Time optimal control example GEKKO

Question

I'm trying to implement a time optimal control problem in GEKKO. In particular, I copied this short code snippet. Reported also here for practicality:

from gekko import GEKKO
import matplotlib.pyplot as plt
import numpy as np
 
# set up the gekko model
m = GEKKO()
 
# set up the time (minimize the time with time scaling)
m.time = np.linspace(0, 1, 100)
 
# set up the variables
POSITION = m.Var(value=0, ub=330, lb=0)
VELOCITY = m.Var(value=0, ub=33, lb=0)
m.fix_final(VELOCITY, 0)
m.fix_final(POSITION, 300)
 
# set up the value we modify over the horizon
tf = m.FV(value=500, lb=0.1)
tf.STATUS = 1
 
# set up the MV
u = m.MV(integer=True, lb=-2, ub=1)
u.STATUS = 1
 
# set up the equations
m.Equation(POSITION.dt() / tf == VELOCITY)
m.Equation(VELOCITY.dt() / tf == u)
 
# set the objective
m.Obj(tf)
 
# set up the options
m.options.IMODE = 6     # optimal control
m.options.SOLVER = 3    # IPOPT
 
# solve
m.solve(disp=False)
 
# print the time
print("Total time taken: " + str(tf.NEWVAL))
 
# plot the results
plt.figure()
plt.subplot(211)
plt.plot(np.linspace(0,1,100)*tf.NEWVAL, POSITION, label='Position')
plt.plot(np.linspace(0,1,100)*tf.NEWVAL, VELOCITY, label='Velocity')
plt.ylabel('Z')
plt.legend()
plt.subplot(212)
plt.plot(np.linspace(0,1,100)*tf.NEWVAL, u, label=r'$u$')
plt.ylabel('u')
plt.xlabel('Time')
plt.legend()
plt.show()

As is, it works fine, but when I want to remove the constraint on the final value of the VELOCITY.

If I comment the line m.fix_final(VELOCITY, 0) the result does not change. It seems like it assumes that the final velocity should be zero, regardless. Furthermore, if I change the final velocity from zero to any other number, I get an error from GEKKO: Exception: @error: Solution Not Found.

The solution should be trivial to find, in particular, if no constraint is imposed on the final veclocity, the optimal control would be to keep accelerating the whole time ().

Any help would be much appreciated! :)

Answer 1

Change the final constraints from m.fix_final(VELOCITY, 0) and m.fix_final(POSITION, 300) to:

p = np.zeros(100); p[-1] = 1
last = m.Param(p)
m.Equation(last*(POSITION-300)>=0)

This applies an inequality constraint at the last node so that POSITION>=300, but it could also be an equality constraint. We also sometimes use a soft constraint such as m.Minimize(last*(POSITION-300)**2) if an infeasible solution prevents the solver from achieving the final condition. Instead, it will try to get the solution as close to the final constraint as possible. When fixing the final value with m.fix_final(), the derivative is also fixed to zero because that variable is no longer calculated. This is a known limitation of Gekko as reported here.