How to represent an if function controlled by three variables？
回答済みI am reproducing a paper on optimizing the configuration of power systems.During this period, I encounter a problem.I don't know how to express the formula in the figure below.M1,N1,N2 are all variables.Mw,PG,P1,P2 are all known.I tried to write the following code.But I don't know whether it's right or wrong.But the entire program reported an error：Model is infeasible.So I need some help and I earnestly hope for your answer.
eps = 0.0001
M_0 = 10000
M_1 = 10000
M_2 = 5
M_3 = 10000
for i in range(8760):
m.addConstr(b1[i] == a5[i] + a1[i])
m.addConstr(b2[i] == a5[i] + a2[i])
m.addConstr(b3[i] == a6[i] + a3[i])
m.addConstr(b4[i] == a6[i] + a4[i])
m.addConstrs(b1[i] >= 2*c1[i]-eps for i in range(8760))
m.addConstrs(b1[i] <= M_2*c1[i]+2-eps for i in range(8760))
m.addConstrs(b2[i] >= 2*c2[i]-eps for i in range(8760))
m.addConstrs(b2[i] <= M_2*c2[i]+2-eps for i in range(8760))
m.addConstrs(b3[i] >= 2*c3[i]-eps for i in range(8760))
m.addConstrs(b3[i] <= M_2*c3[i]+2-eps for i in range(8760))
m.addConstrs(b4[i] >= 2*c4[i]-eps for i in range(8760))
m.addConstrs(b4[i] <= M_2*c4[i]+2-eps for i in range(8760))
m.addConstrs(c1[i] + c2[i] + c3[i] + c4[i] == 1 for i in range(8760))
m.addConstrs(((P_G[t]/P1)*a2[t]-eps <= N1 for t in range(8760)) , "Power11")
m.addConstrs((N1 <= (P_G[t]/P1)*a1[t]+(M_1/P1)*a2[t]-2*eps for t in range(8760)) , "Power12")
m.addConstrs(((P_G[t]/P2)*a4[t]-eps <= N2 for t in range(8760)) , "Power21")
m.addConstrs((N2 <= (P_G[t]/P2)*a3[t]+(M_1/P2)*a4[t]-2*eps for t in range(8760)) , "Power22")
m.addConstrs((M_w[d]*a6[d*24+h]-eps <= M1_t[d,h] for d in range(365) for h in range(24)) , "Power31")
m.addConstrs((M1_t[d,h] <= M_w[d]*a5[d*24+h] + M_0*a6[d*24+h]-2*eps for d in range(365) for h in range(24)) , "Power32")
m.addConstrs((a1[i]+a2[i] == 1 for i in range(8760)) , "indicator_constr1" )
m.addConstrs((a3[i]+a4[i] == 1 for i in range(8760)) , "indicator_constr2" )
m.addConstrs((a5[i]+a6[i] == 1 for i in range(8760)) , "indicator_constr3" )
#情况2：可再生能源净功率大于制淡负荷功率上限
m.addConstrs(((c9[d*24+h] == 1) >> (n1_t[d,h] == N1) for d in range(365) for h in range(24)) , "c1n1")
m.addConstrs(((c9[d*24+h] == 1) >> (p1_t[d,h] == N1*P1) for d in range(365) for h in range(24)) , "c1p1")
m.addConstrs(((c9[d*24+h] == 1) >> (p2_t[d,h] == P_G[d*24+h]-p1_t[d,h]) for d in range(365) for h in range(24)) , "c1p2")
m.addConstrs(((c9[d*24+h] == 1) >> (n2_t[d,h] == p2_t[d,h]/P2) for d in range(365) for h in range(24)) , "c1n2")
#情况1：可再生能源净功率在制淡负荷功率的上下限之间
m.addConstrs(((c10[d*24+h] == 1) >> (n1_t[d,h] == P_G[d*24+h]/P1) for d in range(365) for h in range(24)) , "c2n1")
m.addConstrs(((c10[d*24+h] == 1) >> (p1_t[d,h] == P_G[d*24+h]) for d in range(365) for h in range(24)) , "c2p1")
m.addConstrs(((c10[d*24+h] == 1) >> (n2_t[d,h] == 0) for d in range(365) for h in range(24)) , "c2n2")
m.addConstrs(((c10[d*24+h] == 1) >> (p2_t[d,h] == 0) for d in range(365) for h in range(24)) , "c2p2")
#情况4：完成淡水需求订单后，可再生能源净功率大于蓄电池充电上限
m.addConstrs(((c3[d*24+h] == 1) >> (n2_t[d,h] == N2) for d in range(365) for h in range(24)) , "c3n2")
m.addConstrs(((c3[d*24+h] == 1) >> (p2_t[d,h] == N2*P2) for d in range(365) for h in range(24)) , "c3p2")
m.addConstrs(((c3[d*24+h] == 1) >> (p1_t[d,h] == P_G[d*24+h]-p2_t[d,h]) for d in range(365) for h in range(24)) , "c3p1")
m.addConstrs(((c3[d*24+h] == 1) >> (n1_t[d,h] == p1_t[d,h]/P1) for d in range(365) for h in range(24)) , "c3n1")
#情况3：完成淡水需求订单后，可再生能源净功率在蓄电池充电的上下限之间
m.addConstrs(((c4[d*24+h] == 1) >> (n2_t[d,h] == P_G[d*24+h]/P2) for d in range(365) for h in range(24)) , "c4n2")
m.addConstrs(((c4[d*24+h] == 1) >> (p2_t[d,h] == P_G[d*24+h]) for d in range(365) for h in range(24)) , "c4p2")
m.addConstrs(((c4[d*24+h] == 1) >> (n1_t[d,h] == 0) for d in range(365) for h in range(24)) , "c4n1")
m.addConstrs(((c4[d*24+h] == 1) >> (p1_t[d,h] == 0) for d in range(365) for h in range(24)) , "c4p1")
-
Sorry, I need to make some additional modifications.c9 should be modified to c1，and c10 should be modified to c2.
I have made the following constraints regarding M1：
for d in range(365):
for h in range(25):
if h == 0:
m.addConstr((Q1_t[d,0] == 1))
m.addConstr((M1_t[d,0] == 1))
m.addConstr((Q2_t[d,0] == 1))
else:
m.addConstr((Q1_t[d,h] == gp.quicksum(p1_t[d,j] for j in range(0,h))))
m.addConstr((M1_t[d,h] == 0.2*Q1_t[d,h]))
m.addConstr((Q2_t[d,h] == gp.quicksum(p2_t[d,j] for j in range(0,h))))0 -
-
-
Gurobi Staff
We discuss how to model conditional constraints in the Knowledge Base article How do I model conditional statements in Gurobi? It looks like your formulation already follows the modeling discussed in the Knowledge Base article, but it might have too many uses of \(\epsilon\).
It is hard to follow your code, because the variable names in your code do not match the parameters and variables in your hand-written formula.
In order to investigate the source of infeasibility, please refer to How do I determine why my model is infeasible?
Best regards,
Jaromił0 -
-
Thank you for your reply.It has helped me a lot
0 -
サインインしてコメントを残してください。
コメント
3件のコメント