# Code conversion from GAMS to JuMP (Ipopt Solver)

#1

Is there any better way to convert constraints with multiple indices?
Here’s what I’ve done.

[GAMS]
SETS
i       index of buses         /1*33/
t       index of time slots    /1*24/
s      reduced scenarios      /1*10/
;

* Power flow equations
* Bus: 2~N_Bus
NLECONS1(t,i,s)\$(ord(i) > 1) .. p_g(t,i,s)+p_bes(t,i,s)-p_l(t,i,s)-v(t,i,s)*SUM(j, v(t,j,s)*(GmatDATA(i,j)*cos((2*pi/360)*(theta(t,i,s)-theta(t,j,s)))+BmatDATA(i,j)*sin((2*pi/360)*(theta(t,i,s)-theta(t,j,s)))))=e=0   ;
NLECONS2(t,i,s)\$(ord(i) > 1) .. q_g(t,i,s)+q_bes(t,i,s)-q_l(t,i,s)-v(t,i,s)*SUM(j, v(t,j,s)*(GmatDATA(i,j)*sin((2*pi/360)*(theta(t,i,s)-theta(t,j,s)))-BmatDATA(i,j)*cos((2*pi/360)*(theta(t,i,s)-theta(t,j,s)))))=e=0   ;
* Slack Bus
NLECONS3(t,i,s)\$(ord(i) = 1) .. p_tie(t,s)+p_g(t,i,s)+p_bes(t,i,s)-p_l(t,i,s)-v(t,i,s)*SUM(j, v(t,j,s)*(GmatDATA(i,j)*cos((2*pi/360)*(theta(t,i,s)-theta(t,j,s)))+BmatDATA(i,j)*sin((2*pi/360)*(theta(t,i,s)-theta(t,j,s)))))=e=0  ;
NLECONS4(t,i,s)\$(ord(i) = 1) .. q_tie(t,s)+q_g(t,i,s)+q_bes(t,i,s)-q_l(t,i,s)-v(t,i,s)*SUM(j, v(t,j,s)*(GmatDATA(i,j)*sin((2*pi/360)*(theta(t,i,s)-theta(t,j,s)))-BmatDATA(i,j)*cos((2*pi/360)*(theta(t,i,s)-theta(t,j,s)))))=e=0  ;
[JuMP]
# Bus: 2~N_Bus
ConstraintIndex4 = Array{Tuple{Int64,Int64,Int64}}(N_Scen,N_TimeSlot,N_Bus-1);
for i in 1:N_Scen
for j in 1:N_TimeSlot
for k in 2:N_Bus
ConstraintIndex4[i,j,k-1] = (i,j,k)
end
end
end
@NLconstraint( m, POWERFLOWconstraint1[ (s,t,i) in ConstraintIndex4], Pgen[s,t,i] +Pbess[s,t,i] -Pload[s,t,i] -v[s,t,i] * sum(v[s,t,j]*(data["Gmat"][i,j]*cos((2*pi/360)*(theta[s,t,i]-theta[s,t,j]))+data["Bmat"][i,j]*sin((2*pi/360)*(theta[s,t,i]-theta[s,t,j]))) for j in 1:N_Bus)==0)
@NLconstraint( m, POWERFLOWconstraint2[ (s,t,i) in ConstraintIndex4], Qgen[s,t,i] +Qbess[s,t,i] -Qload[s,t,i] -v[s,t,i] * sum(v[s,t,j]*(data["Gmat"][i,j]*sin((2*pi/360)*(theta[s,t,i]-theta[s,t,j]))-data["Bmat"][i,j]*cos((2*pi/360)*(theta[s,t,i]-theta[s,t,j]))) for j in 1:N_Bus)==0)

# Slack BUS
ConstraintIndex = Array{Tuple{Int64,Int64}}(N_Scen,N_TimeSlot);
for i in 1:N_Scen
for j in 1:N_TimeSlot
ConstraintIndex[i,j] = (i,j)
end
end
@NLconstraint( m, POWERFLOWconstraint3[ (s,t) in ConstraintIndex], Ptie[s,t]+Pgen[s,t,1] +Pbess[s,t,1] -Pload[s,t,1] -v[s,t,1] * sum(v[s,t,j]*(data["Gmat"][1,j]*cos((2*pi/360)*(theta[s,t,1]-theta[s,t,j]))+data["Bmat"][1,j]*sin((2*pi/360)*(theta[s,t,1]-theta[s,t,j]))) for j in 1:N_Bus)==0)
@NLconstraint( m, POWERFLOWconstraint4[ (s,t) in ConstraintIndex], Qtie[s,t]+Qgen[s,t,1] +Qbess[s,t,1] -Qload[s,t,1] -v[s,t,1] * sum(v[s,t,j]*(data["Gmat"][1,j]*sin((2*pi/360)*(theta[s,t,1]-theta[s,t,j]))-data["Bmat"][1,j]*cos((2*pi/360)*(theta[s,t,1]-theta[s,t,j]))) for j in 1:N_Bus)==0)

Ipopt Solvers on GAMS and JuMP show different numbers of nonzeros in equality constraint Jacobian even though total numbers of equality constraints are same.

[GAMS]
Number of nonzeros in equality constraint Jacobian...:   303010
Total number of variables............................:    66864
Total number of equality constraints.................:    39990
Total number of inequality constraints...............:    32640

[JuMP]
Number of nonzeros in equality constraint Jacobian...:  1245790
Total number of variables............................:    66864
Total number of equality constraints.................:    39990
Total number of inequality constraints...............:    32640

GAMS solves the problem in 1~2hrs, but JuMP couldn’t finish the work in days.
I don’t know where this comes from.

[GAMS]

COIN-OR Interior Point Optimizer (Ipopt Library 3.8)
written by A. Waechter

******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
Ipopt is released as open source code under the Common Public License (CPL).
******************************************************************************

NOTE: You are using Ipopt by default with the MUMPS linear solver.
Other linear solvers might be more efficient (see Ipopt documentation).

This is Ipopt version 3.8stable, running with linear solver mumps.

Number of nonzeros in equality constraint Jacobian...:   303010
Number of nonzeros in inequality constraint Jacobian.:    47100
Number of nonzeros in Lagrangian Hessian.............:    87794

Total number of variables............................:    66864
variables with only lower bounds:     8480
variables with lower and upper bounds:    34100
variables with only upper bounds:        0
Total number of equality constraints.................:    39990
Total number of inequality constraints...............:    32640
inequality constraints with only lower bounds:     8640
inequality constraints with lower and upper bounds:        0
inequality constraints with only upper bounds:    24000

[JuMP]

******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
Ipopt is released as open source code under the Eclipse Public License (EPL).
******************************************************************************

This is Ipopt version 3.12.4, running with linear solver mumps.
NOTE: Other linear solvers might be more efficient (see Ipopt documentation).

Number of nonzeros in equality constraint Jacobian...:  1245790
Number of nonzeros in inequality constraint Jacobian.:    72560
Number of nonzeros in Lagrangian Hessian.............: 35581600

Total number of variables............................:    66864
variables with only lower bounds:     8480
variables with lower and upper bounds:    34100
variables with only upper bounds:        0
Total number of equality constraints.................:    39990
Total number of inequality constraints...............:    32640
inequality constraints with only lower bounds:     8640
inequality constraints with lower and upper bounds:        0
inequality constraints with only upper bounds:    24000

#2

JuMP doesn’t do much preprocessing to combine duplicate elements, it makes the solver do some of the work there. Looks like the Hessian in particular is handled better in GAMS than JuMP. You could try the AmplNLWriter package.

#3

Thanks for the comment.
Doesn’t AmplNLWriter calls Ipopt Solver for NLP?
As far as I recall, BonminSolver does.

#4

It’s the same optimization solver but uses a C library for evaluating the Jacobian, Hessian, etc instead of the Julia code used by Ipopt.jl

#5