CODEWITHSUNDEEP
pythonmathematical-optimizationodedifferential-equationsgekko
TexasEngineer
TexasEngineer

Reputation: 734

Change ODE initial condition to meet final target

I'm trying to reach a final condition subject to many Differential and Algebraic Equations. There is the option to change parameters in the model or adjust the initial condition to meet the final target. It is easy to do this with an adjustable parameter, but how can the final condition be met with an adjustable initial condition in Python Gekko?

Optimize parameter

import numpy as np
from gekko import GEKKO
import matplotlib.pyplot as plt

# simulate
m = GEKKO()    # create GEKKO model
k = m.Param(0.1)  # constant
y = m.Var(5.0) # create GEKKO variable
m.Equation(y.dt()==-k*y) # create GEKKO equation
m.time = np.linspace(0,20) # time points
m.options.IMODE = 4
m.solve(disp=False)
plt.plot(m.time,y,'r-',linewidth=2,label='k=0.1')

# adjust k to meet target
m = GEKKO()    # create GEKKO model
k = m.FV(0.1); k.STATUS=1  # adjustable parameter
y = m.Var(5.0)
m.Equation(y.dt()==-k*y)
m.time = np.linspace(0,20,50)
m.options.IMODE = 6
p = np.zeros(50); p[-1]=1
final = m.Param(p)
m.Minimize(final*(y-1)**2)
m.solve(disp=False)
kv = k.value[0]
plt.plot(m.time,y,'k:',linewidth=2,label='k='+str(kv))

plt.plot([20],[1],'bo',markersize=5,label='target')
plt.xlabel('time')
plt.ylabel('y(t)')
plt.legend()
plt.show()

Upvotes: 0

Views: 179

Answers (1)

John Hedengren
John Hedengren

Reputation: 14346

Set fixed_initial=False to calculate the initial condition as shown in the Model Building Gekko documentation.

y = m.Var(5.0,fixed_initial=False)

Additional details are in the Gekko Documentation about m.free_initial(y) as a function if the variable needs to be fixed or calculated after it is defined.

Calculate initial condition

import numpy as np
from gekko import GEKKO
import matplotlib.pyplot as plt

# simulate
m = GEKKO()    # create GEKKO model
k = m.Param(0.1)  # constant
y = m.Var(5.0) # create GEKKO variable
m.Equation(y.dt()==-k*y) # create GEKKO equation
m.time = np.linspace(0,20) # time points
m.options.IMODE = 4
m.solve(disp=False)
plt.plot(m.time,y,'r-',linewidth=2,label='IC=5')

# adjust initial condition to meet target
m = GEKKO()    # create GEKKO model
y = m.Var(5.0,fixed_initial=False)
k = 0.1
m.Equation(y.dt()==-k*y)
m.time = np.linspace(0,20,50)
m.options.IMODE = 6
p = np.zeros(50); p[-1]=1
final = m.Param(p)
m.Minimize(final*(y-1)**2)
m.solve(disp=False)
plt.plot(m.time,y,'k:',linewidth=2,label='IC='+str(y.value[0]))

plt.plot([20],[1],'bo',markersize=5,label='target')
plt.xlabel('time')
plt.ylabel('y(t)')
plt.legend()
plt.show()

Upvotes: 1

Related Questions