英文:
Dyn.Opt problem solution doesn't converge on GEKKO on varying no. of time pts or values of maniplatd var or colc. nodes
问题
I'm solving a Dynamic Optimization problem on gekko but the solution doesn't converge if I vary the number of time points or initial/lb/ub values for manipulated variable or change the collocation nodes. What could be the issue?
我正在解决一个关于gekko的动态优化问题，但如果我改变时间点的数量、操作变量的初始值/下限/上限值，或者改变配点节点，解决方案就无法收敛。可能的问题是什么？
Just changed the values mentioned above, and on some combinations, the solution doesn't converge at all, and the following error message appears:
刚刚更改了上述提到的值，在某些组合下，解决方案根本无法收敛，出现以下错误消息：
raise Exception(response)
Exception: @error: Solution Not Found
Code:
from gekko import GEKKO
import numpy as np
import matplotlib.pyplot as plt
m = GEKKO()
nt = 101 # no. of time steps
m.time = np.linspace(0,147,nt)
# ... (其他代码部分未提供翻译)
代码：
from gekko import GEKKO
import numpy as np
import matplotlib.pyplot as plt
m = GEKKO()
nt = 101 # 时间步数
m.time = np.linspace(0,147,nt)
# ... (其他代码部分未提供翻译)
英文:
I'm solving a Dynamic Optimization problem on gekko but the solution doesn't converge if I vary the number of time points or initial/lb/ub values for manipulated variable or change the collocation nodes. What could be the issue?
Just changed the values mentioned above, and on some combinations, the solution doesn't converge at all, and the following error message appears:
raise Exception(response)
Exception: @error: Solution Not Found
Code:
from gekko import GEKKO
import numpy as np
import matplotlib.pyplot as plt
m = GEKKO()
nt = 101 # no. of time steps
m.time = np.linspace(0,147,nt)
k17 = 0.0508
k26 = 1
k27 = 0.0577
k28 = 0.0104
k32 = 2
k33 = 2
k37 = 0.0016
k38 = 0.0107
k47 = 0.006
k48 = 0.072
k53 = 2
k57 = 0.0201
k58 = 0.082
k61 = 2
k62 = 2
k63 = 2
k72 = 2
k77 = 0.0133
k78 = 0.011
k87 = 0.081
k88 = 0.0148
kd = 0.06
w = 0.043
s1 = m.Var(value=0,lb=0)
s2 = m.Var(value=3,lb=0)
s3 = m.Var(value=0.01,lb=0)
s4 = m.Var(value=0.46,lb=0)
s5 = m.Var(value=0.27,lb=0)
p1 = m.Var(value=0.51,lb=0)
p2 = m.Var(value=0.33,lb=0)
p3 = m.Var(value=0.3,lb=0)
p4 = m.Var(value=0.34e-4,lb=0)
p5 = m.Var(value=2.01,lb=0)
p6 = m.Var(value=0.05,lb=0)
X = m.Var(value=0.09,lb=0)
Xd = m.Var(value=0.02,lb=0)
V = m.Var(value=5,lb=0,ub=10)
u = m.MV(value=0.05,lb=0,ub=0.1) # manipulated variable
u.STATUS = 1
u.DCOST = 0
f1 = (8.443e-4)*X*s1/(8.989e5 + s1)
f2 = (2.481e6)*X*s1*s3/((6.495e4 + s1)*(7.076e2 + s3))
f3 = (3.968e5)*X*s1*s3/((3.723e4 + s1)*(2.782e3 + s3))
f4 = (1.09e2)*X*s3/(0.019+s3)
f5 = 7.283*X*s4/(1.92e3 + s4)
f6 = (3.337e5)*X*s2*s5/((2.719e4 + s2)*(4.488e4 + s5))
f7 = (3.977e3)*X*s2/(9.324e3 + s2)
f8 = (6.697e-6)*X*s2/(0.537+s2)
f9 = (3.261e4)*X*s2/(6.683e5 + s2)
P = np.zeros(nt)
P[-1] = 1.0
final = m.Param(value=P)
# Equations
m.Equation(V.dt()==u)
m.Equation(s1.dt()==-f1-f2-f3-k17*f7+(2-s1)*u/V)
m.Equation(s2.dt()==-f6-k27*f7-k28*f8-f9-s2*u/V)
m.Equation(s3.dt()==-k32*f2-k33*f3-f4+f6-k37*f7-k38*f8+f9-s3*u/V)
m.Equation(s4.dt()==-f5+f6-k47*f7-k48*f8-s4*u/V)
m.Equation(s5.dt()==k53*f3+f5-f6-k57*f7-k58*f8-s5*u/V)
m.Equation(p1.dt()==k61*f1+k62*f2+k63*f3-p1*u/V)
m.Equation(p2.dt()==k72*f2-k77*f7-k78*f8-p2*u/V)
m.Equation(p3.dt()==f4-k87*f7-k88*f8-p3*u/V)
m.Equation(p4.dt()==f8-p4*u/V)
m.Equation(p5.dt()==f7-p5*u/V)
m.Equation(p6.dt()==f5+f9-p6*u/V)
m.Equation(X.dt()==w*X-kd*X*Xd-X*u/V)
m.Equation(Xd.dt()==kd*X*Xd-Xd*u/V)
m.Obj(-final*V*p4)
m.options.IMODE = 6
m.options.NODES = 4
#m.options.COLDSTART=2
#m.options.MAX_ITER=1000
m.solve(disp=True)
p4_ = np.multiply(p4.value,1000)
plt.figure(1)
plt.subplot(2,1,1)
plt.plot(m.time,u.value,'r-')
plt.subplot(2,1,2)
plt.plot(m.time,p4_,'b--')
答案1
得分: 0
以下是代码部分不需要翻译的内容：
There are several reasons that the solver fails to reach a solution with `Solution Not Found`. Here are some of the reasons and things to try:
1. If the solver reports `No Feasible Solution` then the bounds are too restrictive. Widen the variable bounds until a solution is obtained.
2. If the solver reaches the maximum iterations, try changing the solver with `m.options.SOLVER=1` or letting the solver continue for more iterations with `m.options.MAX_ITER=500`.
3. If the solver stops with the message `Unbounded Solution` then try setting bounds on the `MV` values.
4. Try initialization with no degrees of freedom by setting `{FV}.STATUS=0` and `{MV}.STATUS=0`). The simulation should have the same number of equations and variables as a square problem. Switching to `m.options.IMODE=7` is a sequential simulation and can also help with the initialization. Setting `m.options.COLDSTART=1` or `m.options.COLDSTART=2` can help find an initial solution that helps the optimization. Additional initialization resources are available in the [Dynamic Optimization course][1].
Please post a copy of your code for more specific suggestions. There are several [Dynamic Optimization Benchmark problems][2] posted to the Dynamic Optimization Course.
**Response to Edit**
Try setting a smaller time scale until the problem converges. Here is a successful solution with a final time of `0.08`.
```python
m.time = np.linspace(0,0.08,11)
The solver does not find a solution after this time point. It reports an infeasible solution. One cause of this may be variable p4 that makes it infeasible because of the lower bound of zero.
希望这些内容对你有所帮助。如果你有任何其他问题，欢迎提出。
<details>
<summary>英文:</summary>
There are several reasons that the solver fails to reach a solution with `Solution Not Found`. Here are some of the reasons and things to try:
1. If the solver reports `No Feasible Solution` then the bounds are too restrictive. Widen the variable bounds until a solution is obtained.
2. If the solver reaches the maximum iterations, try changing the solver with `m.options.SOLVER=1` or letting the solver continue for more iterations with `m.options.MAX_ITER=500`.
3. If the solver stops with the message `Unbounded Solution` then try setting bounds on the `MV` values.
4. Try initialization with no degrees of freedom by setting `{FV}.STATUS=0` and `{MV}.STATUS=0`). The simulation should have the same number of equations and variables as a square problem. Switching to `m.options.IMODE=7` is a sequential simulation and can also help with the initialization. Setting `m.options.COLDSTART=1` or `m.options.COLDSTART=2` can help find an initial solution that helps the optimization. Additional initialization resources are available in the [Dynamic Optimization course][1].
Please post a copy of your code for more specific suggestions. There are several [Dynamic Optimization Benchmark problems][2] posted to the Dynamic Optimization Course.
**Response to Edit**
Try setting a smaller time scale until the problem converges. Here is a successful solution with a final time of `0.08`.
```python
m.time = np.linspace(0,0.08,11)
The solver does not find a solution after this time point. It reports an infeasible solution. One cause of this may be variable p4 that makes it infeasible because of the lower bound of zero.
Here is the diagnostic code with the plots:
from gekko import GEKKO
import numpy as np
import matplotlib.pyplot as plt
m = GEKKO()
m.time = np.linspace(0,0.08,11)
nt = len(m.time)
k17 = 0.0508
k26 = 1
k27 = 0.0577
k28 = 0.0104
k32 = 2
k33 = 2
k37 = 0.0016
k38 = 0.0107
k47 = 0.006
k48 = 0.072
k53 = 2
k57 = 0.0201
k58 = 0.082
k61 = 2
k62 = 2
k63 = 2
k72 = 2
k77 = 0.0133
k78 = 0.011
k87 = 0.081
k88 = 0.0148
kd = 0.06
w = 0.043
s1 = m.Var(value=0,lb=0)
s2 = m.Var(value=3,lb=0)
s3 = m.Var(value=0.01,lb=0)
s4 = m.Var(value=0.46,lb=0)
s5 = m.Var(value=0.27,lb=0)
p1 = m.Var(value=0.51,lb=0)
p2 = m.Var(value=0.33,lb=0)
p3 = m.Var(value=0.3,lb=0)
p4 = m.Var(value=0.34e-4,lb=0)
p5 = m.Var(value=2.01,lb=0)
p6 = m.Var(value=0.05,lb=0)
X = m.Var(value=0.09,lb=0)
Xd = m.Var(value=0.02,lb=0)
V = m.Var(value=5,lb=0.01)
u = m.MV(value=0.05,lb=0,ub=0.1) # manipulated variable
u.STATUS = 1
u.DCOST = 0
f1 = m.Intermediate((8.443e-4)*X*s1/(8.989e5 + s1))
f2 = m.Intermediate((2.481e6)*X*s1*s3/((6.495e4 + s1)*(7.076e2 + s3)))
f3 = m.Intermediate((3.968e5)*X*s1*s3/((3.723e4 + s1)*(2.782e3 + s3)))
f4 = m.Intermediate((1.09e2)*X*s3/(0.019+s3))
f5 = m.Intermediate(7.283*X*s4/(1.92e3 + s4))
f6 = m.Intermediate((3.337e5)*X*s2*s5/((2.719e4 + s2)*(4.488e4 + s5)))
f7 = m.Intermediate((3.977e3)*X*s2/(9.324e3 + s2))
f8 = m.Intermediate((6.697e-6)*X*s2/(0.537+s2))
f9 = m.Intermediate((3.261e4)*X*s2/(6.683e5 + s2))
P = np.zeros(nt)
P[-1] = 1.0
final = m.Param(value=P)
# Equations
m.Equation(V.dt()==u)
m.Equation(s1.dt()==-f1-f2-f3-k17*f7+(2-s1)*u/V)
m.Equation(s2.dt()==-f6-k27*f7-k28*f8-f9-s2*u/V)
m.Equation(s3.dt()==-k32*f2-k33*f3-f4+f6-k37*f7-k38*f8+f9-s3*u/V)
m.Equation(s4.dt()==-f5+f6-k47*f7-k48*f8-s4*u/V)
m.Equation(s5.dt()==k53*f3+f5-f6-k57*f7-k58*f8-s5*u/V)
m.Equation(p1.dt()==k61*f1+k62*f2+k63*f3-p1*u/V)
m.Equation(p2.dt()==k72*f2-k77*f7-k78*f8-p2*u/V)
m.Equation(p3.dt()==f4-k87*f7-k88*f8-p3*u/V)
m.Equation(p4.dt()==f8-p4*u/V)
m.Equation(p5.dt()==f7-p5*u/V)
m.Equation(p6.dt()==f5+f9-p6*u/V)
m.Equation(X.dt()==w*X-kd*X*Xd-X*u/V)
m.Equation(Xd.dt()==kd*X*Xd-Xd*u/V)
m.Maximize(final*V*p4)
m.options.IMODE = 6
m.options.NODES = 3
m.options.SOLVER = 1
m.options.COLDSTART=0
#m.options.MAX_ITER=1000
m.solve(disp=True)
p4_ = np.multiply(p4.value,1000)
plt.figure(1)
plt.subplot(2,1,1)
plt.plot(m.time,u.value,'r-')
plt.subplot(2,1,2)
plt.plot(m.time,p4_,'b--')
plt.figure(2)
plt.plot(m.time,s1,label='s1')
plt.plot(m.time,s2,label='s2')
plt.plot(m.time,s3,label='s3')
plt.plot(m.time,s4,label='s4')
plt.plot(m.time,s5,label='s5')
plt.legend()
plt.figure(3)
plt.plot(m.time,p1,label='p1')
plt.plot(m.time,p2,label='p2')
plt.plot(m.time,p3,label='p3')
plt.plot(m.time,p4,label='p4')
plt.plot(m.time,p5,label='p5')
plt.plot(m.time,p6,label='p6')
plt.legend()
plt.figure(4)
plt.subplot(3,1,1)
plt.plot(m.time,s1,label='s1')
plt.legend()
plt.subplot(3,1,2)
plt.plot(m.time,s3,label='s3')
plt.legend()
plt.subplot(3,1,3)
plt.plot(m.time,p4,label='p4')
plt.legend()
plt.show()
通过集体智慧和协作来改善编程学习和解决问题的方式。致力于成为全球开发者共同参与的知识库，让每个人都能够通过互相帮助和分享经验来进步。
评论