Hello,
I’m trying to run a model for reconfiguring electrical power grids that considers non-technical losses separately. To do this, I’m adding parameters (PDNT and QDNT) and a variable (Isqr_NT) to an existing model (MISOCP) which uses CPLEX. However, it returns this message (CPLEX 12.6.3.0: QP Hessian is not positive semi-definite).
I searched on this forum and found some commands (like option cplex_options ‘reqconvex=3’;), but I was unsuccessful.
Can anybody help me?
Thanks.
#Sets
set N; #Set of nodes (buses)
set B within N cross N; #Set of branches
#Parameters
param Tb{N}; #Node (Bus) type (Load= 0, Subestation= 1)
param PD{N}; #Active demand
param QD{N}; #Reactive demand
param PDNT{N}; #Active demand of nontechnical losses
param QDNT{N}; #Reactive demand of nontechnical losses
param R{B}; #Resistance
param X{B}; #Reactance
param Z2{B}; #Square of the impedance magnitude
param Imax{B}; #Current limit for main branches
param Sbase; #Base Power
param Vnom; #Nominal voltage
param Vbase; #Base voltage
param Vmax; #Maximum voltage magnitude
param Vmin; #Minimum voltage magnitude
var Isqr{B} >= 0; #Square of the current magnitude
var Isqr_NT{B} >= 0; #Square of the current magnitude for nontechnical losses
var P{B}; #Active power flow
var Q{B}; #Reactive power flow
var b{B}; #Slack variable
var y{B} >= 0, <= 1; #Status of the switch
var beta{N cross N} binary; #For radiality constraint
#Objective function
minimize FO: #(Per unit)
sum {(i,j) in B}(R[i,j]*Isqr[i,j]) + sum {(i,j) in B}(R[i,j]*Isqr_NT[i,j]);
#Constraints
subject to active_balance {i in N}:
sum {(k,i) in B}(P[k,i]) + PS[i] = sum {(i,j) in B}(P[i,j]) + sum {(i,j) in B}(R[i,j]*Isqr[i,j] + R[i,j]*Isqr_NT[i,j]) + PD[i] + PDNT[i];
subject to reactive_balance {i in N}:
sum {(k,i) in B}(Q[k,i]) + QS[i] = sum {(i,j) in B}(Q[i,j]) + sum {(i,j) in B}(X[i,j]*Isqr[i,j] + X[i,j]*Isqr_NT[i,j]) + QD[i] + QDNT[i];
subject to voltage_drop {(i,j) in B}:
Vsqr[i] = Vsqr[j] + b[i,j] - 2*(R[i,j]*P[i,j] + X[i,j]*Q[i,j]) + Z2[i,j]*Isqr[i,j]+ Z2[i,j]*Isqr_NT[i,j];
subject to current_calculation {(i,j) in B}:
Vsqr[j]*Isqr[i,j] + Vsqr[j]*Isqr_NT[i,j] >= P[i,j]^2 + Q[i,j]^2;
subject to slack_variable_1 {(i,j) in B}:
b[i,j] >= -(Vmax^2 - Vmin^2)*(1 - y[i,j]);
subject to slack_variable_2 {(i,j) in B}:
b[i,j] <= (Vmax^2 - Vmin^2)*(1 - y[i,j]);
subject to current_max {(i,j) in B}:
Isqr[i,j] + Isqr_NT[i,j] <= Imax[i,j]^2*y[i,j];
subject to voltage_limits {i in N}:
Vmin^2 <= Vsqr[i] <= Vmax^2;
#Additional constraints to decrease the run time
subject to P_max1 {(i,j) in B}:
P[i,j] <= Vmax*Imax[i,j]*y[i,j];
subject to P_max2 {(i,j) in B}:
-Vmax*Imax[i,j]*y[i,j] <= P[i,j];
subject to Q_max1 {(i,j) in B}:
Q[i,j] <= Vmax*Imax[i,j]*y[i,j];
subject to Q_max2 {(i,j) in B}:
-Vmax*Imax[i,j]*y[i,j] <= Q[i,j];
subject to radiality1 {(i,j) in B}:
beta[i,j] + beta[j,i] = y[i,j];
subject to radiality2 {(i,j) in B: Tb[i]==1}:
beta[i,j] = 0;
subject to radiality3 {(i,j) in B: Tb[j]==1}:
beta[j,i] = 0;
subject to radiality4 {i in N: Tb[i]==0}:
sum {j in N: (i,j)in B || (j,i) in B}(beta[i,j])=1;