GEKKO 多相优化:无解,-1 自由度

GEKKO Multiphase Optimization: No solution, -1 degree of freedom

我现在正在学习 GEKKO,并尝试编写一个简单的一维控制示例,其中一个块具有动态特性

dx/dt = v
dv/dt = 你

有:-5 <= u <= 5

其中 x 是位置,v 是速度,u 是应用的驱动必须在 1 个时间单位内从位置 0 移动到 1,初始速度和最终速度为 0。objective 是:

最小绝对值(u)

它与单相一起工作,结果是 this Bang-Off-Bang 控制。

我现在想尝试将其实现为具有 3 个链接阶段的多阶段轨迹:

  1. u = 5
  2. u = 0
  3. u = -5

其中每个阶段的时间是控制变量。

问题是求解器没有找到解。自由度为 -1,其中给定的 objective 理论上它们应该为 1(总的最终时间),如果 objective 被等式约束替换(这导致求解器说明自由度更小,大约-100)。

这是我的代码:

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

m = GEKKO(remote=True)  # Model

n = 3  # Number of phases: Bang-Off-Bang
max_u = 5.0  # Control can either be max_u, 0 or -max_u

# Options
m.options.NODES = 4  # Two intermediate nodes between collocation points
m.options.SOLVER = 1  # APOPT
m.options.IMODE = 6  # MPC Direct Collocation
m.options.MAX_ITER = 500  # To prevent endless loops
m.options.MV_TYPE = 0  # MVs are constant between endpoints
m.options.DIAGLEVEL = 1  # Show some diagnostics

# Time from 0 to 1 in 31 collocation points for every phase each
m.time = np.linspace(0, 1, 31)

tf = [m.FV(value=0.3) for i in range(n)]
for i in range(n):
    tf[i].STATUS = 1  # make final time controllable

# Collocation variables x and v
x = [m.Var(value=0, fixed_initial=False) for i in range(n)]
v = [m.Var(value=0, fixed_initial=False) for i in range(n)]
u = [max_u, 0, -max_u]  # For the three phases

# Fix start- and endpoint (Inbetween points can be fixed with m.fix())
m.fix_initial(x[0], val=0)
m.fix_initial(v[0], val=0)
m.fix_final(x[n-1], val=1)
m.fix_final(v[n-1], val=0)

# Differential equations describing the system scaled by tf
for i in range(n):
    m.Equation(x[i].dt() == v[i]*tf[i])
    m.Equation(v[i].dt() == u[i]*tf[i])

# Connect phases at endpoints
for i in range(n-1):
    m.Connection(x[i+1], x[i], 1, 'end', 1, 'end')
    m.Connection(x[i+1],'calculated', pos1=1, node1=1)
    m.Connection(v[i+1], v[i], 1, 'end', 1, 'end')
    m.Connection(v[i+1],'calculated', pos1=1, node1=1)

# Make final time equal to 1
m.Minimize((m.sum(tf)-1)**2)
#m.Equation(m.sum(tf) == 1)  # Could be used instead of objective

m.open_folder()

# Run optimization with Diagnostics
m.solve(disp=True)

希望有人能帮帮我,先谢谢了!

软终端约束通常更适合收敛到解决方案。这是硬终端约束的替换:

f = np.zeros(31); f[-1]=1; final=m.Param(f)
m.Minimize(final*(x[n-1]-1)**2)
m.Minimize(final*v[n-1]**2)

#m.fix_final(x[n-1], val=1)
#m.fix_final(v[n-1], val=0)

这给出了三个链接阶段的最佳解决方案。

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

m = GEKKO(remote=True)  # Model

n = 3  # Number of phases: Bang-Off-Bang
max_u = 5.0  # Control can either be max_u, 0 or -max_u

# Options
m.options.NODES = 4  # Two intermediate nodes between collocation points
m.options.SOLVER = 1  # APOPT
m.options.IMODE = 6  # MPC Direct Collocation
m.options.MAX_ITER = 500  # To prevent endless loops
m.options.MV_TYPE = 0  # MVs are constant between endpoints
m.options.DIAGLEVEL = 0  # Show some diagnostics

# Time from 0 to 1 in 31 collocation points for every phase each
m.time = np.linspace(0, 1, 31)

tf = [m.FV(value=0.3) for i in range(n)]
for i in range(n):
    tf[i].STATUS = 1  # make final time controllable

# Collocation variables x and v
x = [m.Var(value=0, fixed_initial=False) for i in range(n)]
v = [m.Var(value=0, fixed_initial=False) for i in range(n)]
u = [max_u, 0, -max_u]  # For the three phases

# Fix start- and endpoint (Inbetween points can be fixed with m.fix())
m.fix_initial(x[0], val=0)
m.fix_initial(v[0], val=0)

f = np.zeros(31); f[-1]=1; final=m.Param(f)
m.Minimize(final*(x[n-1]-1)**2)
m.Minimize(final*v[n-1]**2)

#m.fix_final(x[n-1], val=1)
#m.fix_final(v[n-1], val=0)

# Differential equations describing the system scaled by tf
for i in range(n):
    m.Equation(x[i].dt() == v[i]*tf[i])
    m.Equation(v[i].dt() == u[i]*tf[i])

# Connect phases at endpoints
for i in range(n-1):
    m.Connection(x[i+1], x[i], 1, 'end', 1, 'end')
    m.Connection(x[i+1],'calculated', pos1=1, node1=1)
    m.Connection(v[i+1], v[i], 1, 'end', 1, 'end')
    m.Connection(v[i+1],'calculated', pos1=1, node1=1)

# Make final time equal to 1
m.Minimize((m.sum(tf)-1)**2)
#m.Equation(m.sum(tf) == 1)  # Could be used instead of objective

#m.open_folder()

# Run optimization with Diagnostics
m.solve(disp=True)

# Generate plot
t = [m.time*tf[i].value[0] for i in range(3)]
t[1] += t[0][-1]
t[2] += t[1][-1]
for i in range(3):
    plt.plot(t[i],x[i].value)
    plt.plot(t[i],v[i].value)
plt.xlabel('Time')
plt.show()    

缺少2个自由度(-1 DOF)是因为在使用m.fix_final(x[n-1], val=1)m.fix_final(v[n-1], val=0)时端点处的导数也是固定的。这是 Gekko 的一个已知问题以及变量的定义方式。如果确实需要硬约束,还有另一种方法(创建连接到端点的固定 FV)。