【发布时间】：2021-05-07 09:04:01
【问题描述】：
我正在尝试解决 MINLP 问题，首先使用 IPOPT 求解器获得初始解，然后使用 APOPT 获得混合整数解。但是，调用 APOPT 求解器时出现以下错误：
错误：异常：访问冲突在文件 ./f90/cqp.f90 的第 359 行回溯：不可用，使用 -ftrace=frame 或 -ftrace=full 进行编译错误：找不到“results.json”。检查上面的其他错误详细信息 Traceback（最近一次调用最后一次）：文件“optimisation_stack.py”，第 244 行，在 Optimise_G(t,ob, jofbuses, q, qc, s, oa, k, l, T, G_previous, C , Y, G_previous, G_max, G_min) 文件“optimisation_stack.py”，第 134 行，在 Optimise_G sol = MINLP(xinit, A, B, A_eq, B_eq, LB,UB, t, ob, jofbuses, q, qc, s , oa, k, l, T, G_previous, C, Y, G_previous) 文件“optimisation_stack.py”，第 215 行，在 MINLP m_APOPT.solve(disp = False) 文件“C:\Users\Zineb\AppData\Local\程序\Python\Python37\lib\site-packages\g ekko\gekko.py”，第 2227 行，在解决 self.load_JSON() 文件“C:\Users\Zineb\AppData\Local\Programs\Python\Python37\lib\site-packages\gekko\gk_post_solve.py”中，行13、在load_JSON f = open(os.path.join(self._path,'options.json')) FileNotFoundError: [Errno 2] No such file or directory: 'C:\Users\Zineb\AppData\Local\Temp \tmptdgafg1zgk_model1\options.json'
我的代码如下，我尽量简化：
import numpy as np
from gekko import GEKKO
# Define matrices A,A_eq, and vectors b, b_eq for the optimization
def Optimise_G(t,ob, jofbuses, q, qc, s, oa, k, l, T, G_next, C, Y, G_previous, G_max, G_min):
Mbig_1 = T*C
Mbig_2 = C
nb_phases = len(G_next)
b_max = len(t)
no_lanegroups = len(q)
A_eq = np.zeros(((nb_phases+1)*b_max + 1, (3*nb_phases+3)*b_max+nb_phases))
for i in range(nb_phases):
A_eq[0][i] = 1
#B_eq = np.zeros(((nb_phases+1)*b_max + 1, 1))
B_eq = np.zeros((nb_phases+1)*b_max + 1)
B_eq[0] = C - sum(Y[0:nb_phases])
counter_eq = 0
# G(i)=Ga(i,b)+Gb(i,b)+Gc(i,b)
for b in range(b_max):
for i in range(nb_phases):
counter_eq = counter_eq + 1
A_eq[counter_eq][i] = 1
A_eq[counter_eq][nb_phases*(b+1)+ i] = -1
A_eq[counter_eq][nb_phases*b_max + nb_phases*(b+1) + i] = -1
A_eq[counter_eq][2*nb_phases*b_max + nb_phases*(b+1) + i] = -1
# ya(b)+y(b)+y(c)=1
for b in range(b_max):
counter_eq = counter_eq + 1
A_eq[counter_eq][3*nb_phases*b_max + nb_phases + b] = 1
A_eq[counter_eq][(3*nb_phases+1)*b_max + nb_phases + b] = 1
A_eq[counter_eq][(3*nb_phases+2)*b_max + nb_phases + b] = 1
B_eq[counter_eq] = 1
A = np.zeros((no_lanegroups + (2*3*nb_phases+4)*b_max, (3*nb_phases+3)*b_max+nb_phases))
B = np.zeros(no_lanegroups + (2*3*nb_phases+4)*b_max)
counter = -1
# Sum Gi (i in Ij)>=Gj,min
for j in range(no_lanegroups):
counter = counter + 1
for i in range(k[j], l[j]+1):
A[counter][i-1] = -1
B[counter] = -C*qc[j]/s[j]
# ya(b)G_lb(i)<=Ga(i,b), yb(b)G_lb(i)<=Gb(i,b), yc(b)G_lb(i)<=Gc(i,b)
for b in range(b_max):
for i in range(nb_phases):
counter = counter + 1
A[counter][nb_phases*(b+1)+i] = -1
A[counter][3*nb_phases*b_max + nb_phases + b] = G_min[i]
B[counter] = 0
counter = counter + 1
A[counter][nb_phases*b_max + nb_phases*(b+1) + i] = -1
A[counter][(3*nb_phases+1)*b_max + nb_phases + b] = G_min[i]
B[counter] = 0
counter = counter + 1
A[counter][2*nb_phases*b_max + nb_phases*(b+1) +i] = -1
A[counter][(3*nb_phases+2)*b_max + nb_phases + b] = G_min[i]
B[counter] = 0
# ya(b)Gmax(i)>=Ga(i,b), yb(b)Gmax(i)>=Gb(i,b), yc(b)Gmax(i)>=Gc(i,b)
for b in range(b_max):
for i in range(nb_phases):
counter = counter + 1
A[counter][nb_phases*(b+1) +i] = 1
A[counter][3*nb_phases*b_max + nb_phases + b] = -G_max[i]
B[counter] = 0
counter = counter + 1
A[counter][nb_phases*b_max + nb_phases*(b+1) + i] = 1
A[counter][(3*nb_phases+1)*b_max + nb_phases + b] = -G_max[i]
B[counter] = 0
counter = counter + 1
A[counter][2*nb_phases*b_max + nb_phases*(b+1) +i] = 1
A[counter][(3*nb_phases+2)*b_max + nb_phases + b] = -G_max[i]
B[counter] = 0
# (1-yc(b))t(b)<=(T-1)C+sum(Gi(1:l(jofbuses(b))))+sum(Y(1:l(jofbuses(b))-1))
for b in range(b_max):
counter = counter + 1
A[counter][0:l[jofbuses[b]-1]] = -np.ones((1,l[jofbuses[b]-1]))
A[counter][(3*nb_phases+2)*b_max+nb_phases+b] = -t[b]
B[counter] = -t[b] + (T-1)*C + sum(Y[0:l[jofbuses[b]-1]-1])
# (T-1)C+sum(Gi(1:l(jofbuses(b))))+sum(Y(1:l(jofbuses(b))-1))<=yc(b)t(b)+(1-yc(b))Mbig_1
for b in range(b_max):
counter = counter + 1
A[counter][0:l[jofbuses[b]-1]] = np.ones((1,l[jofbuses[b]-1]))
A[counter][(3*nb_phases+2)*b_max+nb_phases+b] = -t[b] + Mbig_1
B[counter] = Mbig_1 - (T-1)*C - sum(Y[0:l[jofbuses[b]-1]-1])
# -Mbig_2(1-yb(b))<=db(b)=right-hand side of Equation (6)
for b in range(b_max):
counter = counter + 1
constant = q[jofbuses[b]-1]/s[jofbuses[b]-1]*(t[b] - (T-1)*C + sum(G_previous[l[jofbuses[b]-1]:nb_phases]) + sum(Y[l[jofbuses[b]-1] -1:nb_phases]))+ (T-1)*C + sum(Y[0:k[jofbuses[b]-1]-1]) - t[b]
A[counter][0:k[jofbuses[b]-1]-1] = -np.ones((1,k[jofbuses[b]-1]-1))
A[counter][(3*nb_phases+1)*b_max + nb_phases + b] = Mbig_2
B[counter] = constant + Mbig_2
# db(b)<=Mbig_2 yb(b)
for b in range(b_max):
counter = counter + 1
constant = q[jofbuses[b]-1]/s[jofbuses[b]-1]*(t[b] - (T-1)*C +sum(G_previous[l[jofbuses[b]-1]:nb_phases]) + sum(Y[l[jofbuses[b]-1] -1:nb_phases]))+ (T-1)*C + sum(Y[0:k[jofbuses[b]-1]-1]) - t[b]
A[counter][0:k[jofbuses[b]-1]-1] = np.ones((1,k[jofbuses[b]-1]-1))
A[counter][(3*nb_phases+1)*b_max + nb_phases + b] = -Mbig_2
B[counter] = -constant
#Lower Bound LB
LB_zeros = np.zeros(3*b_max*(nb_phases+1))
G_min = np.array(G_min)
LB = np.append(G_min, LB_zeros)
#Upper Bound UB
UB = np.ones(3*b_max)
G_max = np.array(G_max)
for i in range(3*b_max+1):
UB = np.concatenate((G_max,UB))
xinit = np.array([(a+b)/2 for a, b in zip(UB, LB)])
sol = MINLP(xinit, A, B, A_eq, B_eq, LB ,UB, t, ob, jofbuses, q, qc, s, oa, k, l, T, G_previous, C, Y, G_previous)
目标函数：
def objective_fun(x, t, ob, jofbuses, q, qc, s, oa, k, l, T, G_next, C, Y, G_previous):
nb_phases = len(G_next)
b_max = len(t)
no_lanegroups = len(q)
obj = 0
obj_a = 0
obj_b = 0
G = x[0:nb_phases]
for j in range(no_lanegroups):
delay_a = 0.5*q[j]/(1-q[j]/s[j]) * (pow((sum(G_previous[l[j]:nb_phases]) + sum(G[0:k[j]-1]) + sum(Y[l[j]-1:nb_phases]) + sum(Y[0:k[j]-1])),2) + pow(sum(G[l[j]:nb_phases]) + sum(G_next[0:k[j]-1]) + sum(Y[l[j]-1:nb_phases]) + sum(Y[0:k[j]-1]),2))
obj = obj + oa*delay_a
obj_a = obj_a + oa*delay_a
for b in range(b_max):
delay_b1 = x[(3*nb_phases+1)*b_max + nb_phases + b]*(q[jofbuses[b]-1]/s[jofbuses[b]-1] * (t[b] - (T-1)*C + sum(G_previous[l[jofbuses[b]-1]:nb_phases]) + sum(Y[l[jofbuses[b]-1] -1:nb_phases])) + (T-1)*C - t[b] + sum(Y[0:k[jofbuses[b]-1]-1]))
delay_b2 = x[(3*nb_phases+2)*b_max + nb_phases + b-1]*(q[jofbuses[b]-1]/s[jofbuses[b]-1] * (t[b] - (T-1)*C - sum(Y[0:l[jofbuses[b]-1]-1])) + T*C + sum(G_next[0:k[jofbuses[b]-1]-1]) + sum(Y[0:k[jofbuses[b]-1]-1]) - t[b])
delay_b3 = sum(x[nb_phases*b_max + nb_phases*b:nb_phases*b_max + nb_phases*b+k[jofbuses[b]-1]-1]) - q[jofbuses[b]-1]/s[jofbuses[b]-1]*sum(x[2*nb_phases*b_max + nb_phases*b:2*nb_phases*b_max + nb_phases*b +l[jofbuses[b]-1]])
delay_b = delay_b1+delay_b2 +delay_b3
obj = obj + delay_b*ob[b]
obj_b = obj_b + delay_b*ob[b]
return obj
MINLP 求解器：
def MINLP(xinit, A, B, A_eq, B_eq, LB ,UB, t, ob, jofbuses, q, qc, s, oa, k, l, T, G_next, C, Y, G_previous):
nb_phases = len(G_next)
b_max = len(t)
## First Solver: IPOPT to get an initial guess
m_IPOPT = GEKKO(remote = False)
m_IPOPT.options.SOLVER = 3 #(IPOPT)
# Array Variable
rows = nb_phases + 3*b_max*(nb_phases+1)#48
x_initial = np.empty(rows,dtype=object)
for i in range(3*nb_phases*b_max+nb_phases+1):
x_initial[i] = m_IPOPT.Var(value = xinit[i], lb = LB[i], ub = UB[i])#, integer = False)
for i in range(3*nb_phases*b_max+nb_phases+1, (3*nb_phases+3)*b_max+nb_phases):
x_initial[i] = m_IPOPT.Var(value = xinit[i], lb = LB[i], ub = UB[i])#, integer = True)
# Constraints
m_IPOPT.axb(A,B,x_initial,etype = '<=',sparse=False)
m_IPOPT.axb(A_eq,B_eq,x_initial,etype = '=',sparse=False)
# Objective Function
f = objective_fun(x_initial, t, ob, jofbuses, q, qc, s, oa, k, l, T, G_next, C, Y, G_previous)
m_IPOPT.Obj(f)
#Solver
m_IPOPT.solve(disp = True)
####################################################################################################
## Second Solver: APOPT to solve MINLP
m_APOPT = GEKKO(remote = False)
m_APOPT.options.SOLVER = 1 #(APOPT)
x = np.empty(rows,dtype=object)
for i in range(3*nb_phases*b_max+nb_phases+1):
x[i] = m_APOPT.Var(value = x_initial[i], lb = LB[i], ub = UB[i], integer = False)
for i in range(3*nb_phases*b_max+nb_phases+1, (3*nb_phases+3)*b_max+nb_phases):
x[i] = m_APOPT.Var(value = x_initial[i], lb = LB[i], ub = UB[i], integer = True)
# Constraints
m_APOPT.axb(A,B,x,etype = '<=',sparse=False)
m_APOPT.axb(A_eq,B_eq,x,etype = '=',sparse=False)
# Objective Function
f = objective_fun(x, t, ob, jofbuses, q, qc, s, oa, k, l, T, G_next, C, Y, G_previous)
m_APOPT.Obj(f)
#Solver
m_APOPT.solve(disp = False)
return x
定义参数并调用Optimise_G：
#Define Parameters
C = 120
T = 31
b_max = 2
G_base = [12,31,12,11,1,41]
G_max = [106, 106, 106, 106, 106, 106]
G_min = [7,3,7,10,0,7]
G_previous = [7.3333333, 34.16763, 7.1333333, 10.0, 2.2008602e-16, 47.365703]
Y = [2, 2, 3, 2, 2, 3]
jofbuses = [9,3]
k = [3,1,6,4,1,2,5,4,2,3]
l = [3,2,6,5,1,3,6,4,3,4]
nb_phases = 6
oa = 1.25
ob = [39,42]
t = [3600.2, 3603.5]
q = [0.038053758888888886, 0.215206065, 0.11325116416666667, 0.06299876472222223,0.02800455611111111,0.18878488361111112,0.2970903402777778, 0.01876728472222222, 0.2192723663888889, 0.06132227222222222]
qc = [0.04083333333333333, 0.2388888888888889, 0.10555555555555556, 0.0525, 0.030555555555555555, 0.20444444444444446,0.31083333333333335, 0.018333333333333333, 0.12777777777777777, 0.07138888888888889]
s = [1.0, 1.0, 1.0, 1.0, 0.5, 1.0, 1.0, 0.5, 0.5, 0.5]
nb_phases = len(G_base)
G_max = []
for i in range(nb_phases):
G_max.append(C - sum(Y[0:nb_phases]))
Optimise_G(t,ob, jofbuses, q, qc, s, oa, k, l, T, G_previous, C, Y, G_previous, G_max, G_min)
有没有办法解决这个问题？ 非常感谢！
【问题讨论】：
标签： python optimization gekko