Блог

Я попытался решить дифференциальное уравнение, имеющее логическое условие в Гекко. Я знаю,что Гекко не любит эти вещи, но я предполагал, что простая функция if3() для переключения между двумя заданными выражениями (R1, R2) даст какое-то решение. Вот простой пример кода, который не удался — решение не найдено.

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

A=1
B=1e-5
C=2
D=0.01 
G=1
y0=120

m = GEKKO()    # create GEKKO model
nt = 101
m.time = np.linspace(0,100,nt) # time points

y = m.Var(y0) # create GEKKO variable

R1 = m.Intermediate(-(y A-A/(y/C 1)**(B/D))/G*(D*y D*C))
R2 = m.Intermediate(-0.1*y)

z = m.Var()              # This way - Solution not found
z = m.if3(y-60,R1,R2)    #
m.Equation(y.dt()== z)   # 

#m.Equation(y.dt()== R1)
#m.Equation(y.dt()== R2)

# solve ODE
m.options.IMODE = 4
m.options.NODES = 5
m.solve(disp=False)

# plot results
plt.plot(m.time,y)
plt.xlabel('time')
plt.ylabel('y(t)')
 

К сожалению, кажется, я еще не понял, как преодолеть эти проблемы, прочитав статью о логических условиях в оптимизации.

с уважением,

Радован

Опубликованный сценарий выдает ошибку:

 Exception:  @error: Solution Not Found
 

Чтобы просмотреть дополнительные сведения об ошибке, включите дисплей с disp=True помощью . Это дает дополнительную информацию, указывающую на то, что проблема неосуществима.

  Number of state variables:    2800
 Number of total equations: -  2000
 Number of slack variables: -  800
 ---------------------------------------
 Degrees of freedom       :    0
 
 ----------------------------------------------
 Dynamic Simulation with APOPT Solver
 ----------------------------------------------
Iter:     1 I: -1 Tm:     13.07 NLPi:    9 Dpth:    0 Lvs:    0 Obj:  0.00E 00 Gap:       NaN
 Warning: no more possible trial points and no integer solution
 Maximum iterations
 
 ---------------------------------------------------
 Solver         :  APOPT (v1.0)
 Solution time  :  13.116200000000001 sec
 Objective      :  0.
 Unsuccessful with error code  0
 ---------------------------------------------------
 
 Creating file: infeasibilities.txt
 Use command apm_get(server,app,'infeasibilities.txt') to retrieve file
 @error: Solution Not Found
 

Пара проблем с моделью:

• z = m.Var() не требуется, потому z = if3() что функция определяет переменную.
• IMODE=6 это необходимо, потому что в задаче есть слабые переменные, и для решения проблемы необходим режим оптимизации.
• Попробуй if2() вместо if3() этого . Кажется, что это лучше работает с формой MPCC, чем со смешанной целочисленной формой.
• Решатель APOPT лучше справляется с этой проблемой. Попробуй m.options.SOLVER=1 .

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

A=1; B=1e-5; C=2; D=0.01; G=1
y0=120

m = GEKKO(remote=False)    # create GEKKO model
nt = 101
m.time = np.linspace(0,50,nt) # time points

y = m.Var(y0) # create GEKKO variable

R1 = m.Intermediate(-(y A-A/(y/C 1)**(B/D))/G*(D*y D*C))
R2 = m.Intermediate(-0.1*y)

z = m.if2(y-60,0,1) 
m.Equation(y.dt()== (1-z)*R1 z*R2)

# solve ODE
m.options.IMODE  = 6
m.options.SOLVER = 1
m.options.NODES  = 5
m.solve(disp=True)

# plot results
plt.subplot(2,1,1)
plt.plot(m.time,y)
plt.plot([0,50],[60,60],'r--')
plt.ylabel('y(t)'); plt.grid()
plt.subplot(2,1,2)
plt.plot(m.time,z)
plt.ylabel('z(t)'); plt.grid()
plt.xlabel('time')
plt.show()
 

Альтернативой использованию оператора if в модели является интеграция в цикл и проверка условия y>60 для переключения на другую модель.