GEKKO RTO против MP C MODES - PullRequest
Новая
       18

GEKKO RTO против MP C MODES

2 голосов
/ 20 марта 2020

Это вопрос, полученный из этого one . После публикации моего вопроса я нашел решение (больше похожее на патч, заставляющий оптимизатор оптимизировать). Есть что-то, что сбивает меня с толку. Джон Хеденгрен правильно указывает, что b=1.0 в ODE приводит к недостижимому решению с IMODE=6. Тем не менее, в моей работе с IMODE=3 я нашел решение.

Я пытаюсь понять, что здесь происходит, читая документацию GEKKO для IMODE=3 и 6, но это мне не понятно

IMODE=3

RTO Оптимизация в реальном времени (RTO) - это стационарный режим, который позволяет принимать переменные решения (типы FV или MV с STATUS = 1) ) или дополнительные переменные, превышающие количество уравнений. Целевая функция направляет выбор дополнительных переменных для выбора оптимального выполнимого решения. RTO является режимом по умолчанию для Gekko, если m.options.IMODE не указано.

IMODE=6

MP C Модель прогнозирующего управления (MP C ) реализовано с IMODE = 6 в качестве одновременного решения или с IMODE = 9 в качестве метода последовательной съемки.

Почему b = 1. работает в одном режиме, но не в другом?

Это моя нестабильная работа с IMODE=3 и b=1.0: 

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

m = GEKKO(remote=False)
m.time = np.linspace(0,23,24)

#initialize variables
T_e = [50.,50.,50.,50.,45.,45.,45.,60.,60.,63.,\
       64.,45.,45.,50.,52.,53.,53.,54.,54.,53.,52.,51.,50.,45.]
temp_low = [55.,55.,55.,55.,55.,55.,55.,68.,68.,68.,68.,55.,55.,68.,\
            68.,68.,68.,55.,55.,55.,55.,55.,55.,55.]
temp_upper = [75.,75.,75.,75.,75.,75.,75.,70.,70.,70.,70.,75.,75.,\
              70.,70.,70.,70.,75.,75.,75.,75.,75.,75.,75.]
TOU = [0.05,0.05,0.05,0.05,0.05,0.05,0.05,200.,200.,\
       200.,200.,200.,200.,200.,200.,200.,200.,200.,\
       200.,200.,200.,0.05,0.05,0.05]

b = m.Param(value=1.)
k = m.Param(value=0.05)

u = [m.MV(0.,lb=0.,ub=1.) for i in range(24)]

# Controlled Variable
T = [m.SV(60.,lb=temp_low[i],ub=temp_upper[i]) for i in range(24)]

for i in range(24):
    u[i].STATUS = 1

for i in range(23):
    m.Equation( T[i+1]-T[i]-k*(T_e[i]-T[i])-b*u[i]==0.0 )

m.Obj(np.dot(TOU,u))

m.options.IMODE = 3
m.solve(debug=True)
myu =[u[0:][i][0] for i in range(24)]
print myu
myt =[T[0:][i][0] for i in range(24)]
plt.plot(myt)
plt.plot(temp_low)
plt.plot(temp_upper)
plt.show()
fig, ax1 = plt.subplots()
ax2 = ax1.twinx()
ax1.plot(myu,color='b')
ax2.plot(TOU,color='k')
plt.show()

Результаты:

enter image description here

Temperature

Input and Hourly cost

1 Ответ

2 голосов
/ 20 марта 2020

Разница между невозможным IMODE=6 и возможным IMODE=3 заключается в том, что корпус IMODE=3 позволяет оптимизатору отрегулировать начальное условие температуры. Оптимизатор распознает, что начальное условие может быть изменено, и поэтому он изменяет его до 75, чтобы оба оставались выполнимыми, а также минимизировали будущее потребление энергии.

Feasible with hard constraints

from gekko import GEKKO
import numpy as np

m = GEKKO(remote=False)
m.time = np.linspace(0,23,24)

#initialize variables
T_external = [50.,50.,50.,50.,45.,45.,45.,60.,60.,63.,\
              64.,45.,45.,50.,52.,53.,53.,54.,54.,\
              53.,52.,51.,50.,45.]
temp_low = [55.,55.,55.,55.,55.,55.,55.,68.,68.,68.,68.,\
            55.,55.,68.,68.,68.,68.,55.,55.,55.,55.,55.,55.,55.]
temp_upper = [75.,75.,75.,75.,75.,75.,75.,70.,70.,70.,70.,75.,\
              75.,70.,70.,70.,70.,75.,75.,75.,75.,75.,75.,75.]
TOU_v = [0.05,0.05,0.05,0.05,0.05,0.05,0.05,200.,200.,200.,200.,\
         200.,200.,200.,200.,200.,200.,200.,200.,200.,200.,0.05,\
         0.05,0.05]

b = m.Param(value=1.)
k = m.Param(value=0.05)
T_e = m.Param(value=T_external)
TL = m.Param(value=temp_low)
TH = m.Param(value=temp_upper)
TOU = m.Param(value=TOU_v)

u = m.MV(lb=0, ub=1)
u.STATUS = 1  # allow optimizer to change

# Controlled Variable
T = m.SV(value=75)

m.Equations([T>=TL,T<=TH])
m.Equation(T.dt() == k*(T_e-T) + b*u)

m.Minimize(TOU*u)

m.options.IMODE = 6
m.solve(disp=True,debug=True)

import matplotlib.pyplot as plt
plt.subplot(2,1,1)
plt.plot(m.time,temp_low,'k--')
plt.plot(m.time,temp_upper,'k--')
plt.plot(m.time,T.value,'r-')
plt.ylabel('Temperature')
plt.subplot(2,1,2)
plt.step(m.time,u.value,'b:')
plt.ylabel('Heater')
plt.xlabel('Time (hr)')
plt.show()

Если вы уйдете на другой день (48 часов), вы, вероятно, увидите, что проблема в конечном итоге станет неосуществимой, поскольку меньший нагреватель b=1 не сможет соответствовать нижнему ограничению температуры.

Одним из преимуществ использования IMODE=6 является то, что вы можете написать дифференциальное уравнение вместо того, чтобы выполнять дискретизацию самостоятельно. С IMODE=3 вы используете метод Эйлера для дифференциального уравнения. Дискретизация по умолчанию для IMODE>=4 равна NODES=2, что эквивалентно вашему методу конечных разностей Эйлера. Установка NODES=3-6 повышает точность при ортогональной коллокации на конечных элементах .

...