Я работаю над задачей динамической оптимизации, очень похожей на эту: сокращение модели с промежуточными звеньями в GEKKO. Я хочу найти лучший температурный профиль, который максимизирует концентрацию биодизеля в последний момент времени (x4(tf)), и сравнить его с полиномиальным температурным профилем.
Я хотел бы сделать его плавным и дифференцируемый. Используя только T = m.MV(value=298,ub=398,lb=298) и T.STATUS = 1, это нереально для переменной состояния. Я пробовал играть со всевозможными значениями для T.DCOST, T.DMAX, MV_STEP_HOR,... но ничего не получилось, даже не похоже на пример выше.
Есть ли еще один параметр настройки, который мне здесь не хватает?
С уважением и заранее большое спасибо!
Вот мой код:
import numpy as np
import matplotlib.pyplot as plt
from gekko import GEKKO
m = GEKKO()
tf = 100 # Tempo final # Final time
m.time = np.linspace(0, tf, tf) # Intervalo de tempo # Time interval
# Temperatura
T = m.MV(value=298, ub=338, lb=298)
T.STATUS = 1
T.DCOST = 1e-5
T.DMAX = 0.7
#print(T.DCOST)
# Variáveis de decisão
# Decision variables
x1 = m.Var(value=0.3226)
x2 = m.Var(value=0)
x3 = m.Var(value=0)
x4 = m.Var(value=0)
x5 = m.Var(value=1.9356)
x6 = m.Var(value=0)
# Constantes da Equação de Arrhenius para velocidade de reação
# Arrhenius Equation
A1, b1 = 3.92e7, 6614.83
A2, b2 = 5.77e5, 4997.98
A3, b3 = 5.88e12, 9993.96
A4, b4 = 0.98e10, 7366.64
A5, b5 = 5.35e3, 3231.18
A6, b6 = 2.15e4, 4824.87
k1 = m.Intermediate(A1*m.exp(-b1/T))
k2 = m.Intermediate(A2*m.exp(-b2/T))
k3 = m.Intermediate(A3*m.exp(-b3/T))
k4 = m.Intermediate(A4*m.exp(-b4/T))
k5 = m.Intermediate(A5*m.exp(-b5/T))
k6 = m.Intermediate(A6*m.exp(-b6/T))
# Modelo dinâmico
# Dynamic model
m.Equation(x1.dt() == -k1*x1*x5 + k2*x2*x4)
m.Equation(x2.dt() == k1*x1*x5 - k2*x2*x4 - k3*x2*x5 + k4*x3*x4)
m.Equation(x3.dt() == k3*x2*x5 - k4*x3*x4 - k5*x3*x5 + k6*x6*x4)
m.Equation(x4.dt() == k1*x1*x5 - k2*x2*x4 + k3*x2*x5 - k4*x3*x4 + k5*x3*x5 - k6*x6*x4)
m.Equation(x5.dt() == -(k1*x1*x5 - k2*x2*x4 + k3*x2*x5 - k4*x3*x4 + k5*x3*x5 - k6*x6*x4))
m.Equation(x6.dt() == k5*x3*x5 - k6*x6*x4)
p = np.zeros(tf)
p[-1] = 1.0
final = m.Param(value=p)
m.Maximize(x4*final)
m.options.IMODE = 6
m.solve(disp=False)
print(f'x4(tf): {x4.value[-1]:0.4f}')
print(T.value[-1])
# Plot de concentração e perfil de temperatura
# Concentration and temperature plot
plt.figure(dpi=300)
plt.plot(m.time, x4.value, label='Biodiesel')
plt.xlabel('Tempo (min)')
plt.ylabel('Concentração (mol/L)')
plt.legend(loc='best')
plt.grid(True)
plt.figure(dpi=300)
plt.plot(m.time, T.value, label=f'GEKKO {T.DCOST}', c='r')
# plt.axhline(y=338, c='y', linestyle='--')
# plt.axhline(y=298, c='y', linestyle='--')
plt.xlabel('Tempo (min)')
plt.ylabel('Temperatura (K)')
plt.legend(loc='best')
plt.grid(True)
Подробнее здесь: https://stackoverflow.com/questions/786 ... n-in-gekko