Good Optimization Modeling Practices

with GAMS
All You Wanted to Know About Practical Optimization but Were Afraid to Ask

Andres Ramos
http://www.iit.comillas.edu/aramos/
Andres.Ramos@comillas.edu

Good Optimization Modeling Practices with GAMS. July 2016 1

“The disciples who received my instructions, and could themselves
comprehend them, were seventy-seven individuals. They were all
scholars of extraordinary ability.” Confucius

Good Optimization Modeling Practices with GAMS. July 2016 2

Do not confuse the ingredients of the recipe
• Mathematical formulation
• LP, MIP, QCP, MCP
• Language
• GAMS
• Solver
• CPLEX, GUROBI, PATH
• Solver algorithm
• Primal simplex, dual simplex, interior
point
• Input/output interfaces
• Text file, CSV, Excel, Matlab, Access
• Operating system
• Windows, Linux, MacOS
• User-developed algorithm
• Benders decomposition, Lagrangian
relaxation, GA
• Stochastic extensions
• EMP

Source:
Good http://www.gams.com/presentations/present_modlang.pdf
Optimization Modeling Practices with GAMS. July 2016 3
Questions to deal with
• What you don’t know how to do it in GAMS?

• What is the most advanced features you know in GAMS?

• What is the most important advantage/disadvantage of GAMS for

you?

• What you would like to do and has not been able to?

Good Optimization Modeling Practices with GAMS. July 2016 4

 Few and practical tips & tricks
• It’s not a systematic approach to teach basic/advanced GAMS features,
just selected features I have used in several models

I have gambas I have chopitos

I have croquetas I have jamón
I have morcillas I have ensalá
I have una hueva mu bien aliña.

• It is optimization for
shepherds (i.e., practitioners)
Good Optimization Modeling Practices with GAMS. July 2016 5
Contents
1. Programming Style

2. GAMS Code

3. Look and Feel

4. Mathematical Formulation

5. Advanced Algorithms

Good Optimization Modeling Practices with GAMS. July 2016 6

 1. Programming Style
2. GAMS Code
3. Look and Feel
4. Mathematical Formulation
5. Advanced Algorithms
1

Programming Style
Programming
• Discipline whose control is basic in many engineering projects
• Science: thinking, discipline, rigorousness and experimentation
• Art: beauty and elegance

• A good design is fundamental

• Before writing any code the optimization problem must be written
algebraically
• Learning by reading
• Coding by gradual refinement, incremental implementation
• Use a mockup for development and verification of the model
• Be careful with the details (“God is in the detail”)

B.W. Kernighan and P.J. Plauger, The Elements of Programming Style, McGraw Hill,
Good Optimization Modeling Practices with GAMS. July 2016
New York, 1978 http://en.wikipedia.org/wiki/The_Elements_of_Programming_Style 8
Which is the most Spanish beautiful landscape?

Aigüestortes
National Park Beach of
Menorca
Mediterranean Guadalquivir
oak wood marshland

Good Optimization Modeling Practices with GAMS. July 2016 9

General recommendations
• Act according to the Pareto principle
• 20 % takes to create the first prototype
• 80 % of code development is devoted to http://www.jetbrains.com/

maintenance and refinement

• MAINTAINABILITY and reusability are crucial
• Code is developed to be read by humans,
not by machines. Write code for
understanding the model, not for
obscuring it.
• Say what you mean, simply and directly.
• Don’t stop with your first draft. Refine it.
Good Optimization Modeling Practices with GAMS. July 2016 10
http://blog.codinghorror.com/when-understanding-means-rewriting/
Procrastination
• Don’t procrastinate when coding

Good Optimization Modeling Practices with GAMS. July 2016 11

Documentation. Comments
• It is a crucial task in code development
• In particular, GAMS was born to explicitly include documentation into the
code.
• Code must be self-documented
• Illustrative comments and well localized
• Make sure comments and code agree.
• Don’t just echo the code with comments - make every
comment count.
• Don’t comment bad code or tricks - rewrite it.
• Don’t patch bad code - rewrite it.
• Don’t over-comment.

Good Optimization Modeling Practices with GAMS. July 2016 12

Code style
• Any project manager ought to define the style before starting up a
multiple participant project (or maybe just for his/her own help)
• Systematic and consistent use of uppercase and lowercase letters
• Use lowercase letters instead of uppercase. We are more used to read
lowercase letters.
• Clean code, take care of the aesthetics when coding
• Aesthetics is as important as the content. Code must be read at a glance.
• Format the code to help the reader understand it.
• Indent to show the logical structure of a program.
• Keep coherence in the coding rules (indent in repetitive sentences)
• Align code to show patterns.
• Make the reading easier (parallelism among consecutive similar sentences,
indent)
• Use meaningful and long names for identifiers. Same use of identifiers
in different parts of the code.

Good Optimization Modeling Practices with GAMS. July 2016 13

Efficiency vs. Clarity
• Make it clear and right before you make it faster
• Keep it simple to make it faster
• Don’t sacrifice clarity for small gains in efficiency

Good Optimization Modeling Practices with GAMS. July 2016 14

 1. Programming Style
2. GAMS Code
3. Look and Feel
4. Mathematical Formulation
5. Advanced Algorithms
2

GAMS Code
Search, compare and if you find something better use it
v. 24.8.4

v. 0.5.1

v. 5.0.1

Good Optimization Modeling Practices with GAMS. July 2016 17

Optimization under
Microsoft Excel
SolverStudio
(http://solverstudio.org/)

OpenSolver (https://opensolver.org/)
Solver from FrontlineSolvers
(http://www.solver.com/premium-
solver-platform)
Good Optimization Modeling Practices with GAMS. July 2016 18
DCLF in Microsoft Excel

Good Optimization Modeling Practices with GAMS. July 2016 19

Pyomo. Optimization modeling in Python
http://www.pyomo.org/documentation

Good Optimization Modeling Practices with GAMS. July 2016 20

Interfaces, languages, solvers Solver
IBM CPLEX
Gurobi
Interface
(graphical) XPRESS

Excel Mathematical Language Algebraic Language GLPK

Access GAMS CBC

SQL AMPL PATH

Matlab AIMMS
Python Pyomo
Julia JuMP
MatLab
Solver Studio

Good Optimization Modeling Practices with GAMS. July 2016 21

Learning by reading first, and then by doing
• GAMS Model Libraries
(https://www.gams.com/modlibs/)

• Decision Support Models in the Electric Power Industry

(http://www.iit.comillas.edu/aramos/Ramos_CV.htm#ModelosAyudaDecision)

Good Optimization Modeling Practices with GAMS. July 2016 22

GAMS (General Algebraic Modeling System)

Good Optimization Modeling Practices

GAMS with GAMS.
birth: 1976July 2016Bank23
World slide
Developing in GAMS
• Development environment gamside

aaa.gpr
• Documentation
• GAMS Documentation Center https://www.gams.com/help/index.jsp
• GAMS Support Wiki http://support.gams.com/
• Bruce McCarl's GAMS Newsletter https://www.gams.com/maillist/newsletter.htm
• Users guide Help > GAMS Users Guide
• Solvers guide Help > Expanded GAMS Guide

• Model: FileName.gms
• Results: FileName.lst
• Process log: FileName.log
Good Optimization Modeling Practices with GAMS. July 2016 24
My first minimalist optimization model
positive variables x1, x2
variable z

equations of, e1, e2, e3 ;

of .. 3*x1 + 5*x2 =e= z ;

e1 .. x1 =l= 4 ;
e2 .. 2*x2 =l= 12 ;
e3 .. 3*x1 + 2*x2 =l= 18 ;

model minimalist / all /

solve minimalist maximizing z using LP

max = 3 +5
,
≤4
2 ≤ 12
3 +2 ≤ 18
, ≥0
Good Optimization Modeling Practices with GAMS. July 2016 25
Blocks in a GAMS model
• Mandatory
variables
equations
model
solve

• Optional
sets: (alias)
• alias (i,j) i and j can be used indistinctly
• Checking of domain indexes
data: scalars, parameters, table

Good Optimization Modeling Practices with GAMS. July 2016 26

 My first transportation model (classical organization)
sets min
I origins / VIGO, ALGECIRAS /
J destinations / MADRID, BARCELONA, VALENCIA /

parameters ≤ ∀#
pA(i) origin capacity
/ VIGO 350
ALGECIRAS 700 /
≥ $ ∀%
pB(j) destination demand
/ MADRID 400 A. Mizielinska y D. Mizielinski Atlas del mundo: Un insólito viaje por las
BARCELONA 450 ≥0 mil curiosidades y maravillas del mundo Ed. Maeva 2015
VALENCIA 150 /

table pC(i,j) per unit transportation cost

MADRID BARCELONA VALENCIA
VIGO 0.06 0.12 0.09
ALGECIRAS 0.05 0.15 0.11

variables
vX(i,j) units transported
vCost transportation cost

positive variable vX

equations
eCost transportation cost
eCapacity(i) maximum capacity of each origin
eDemand (j) demand supply at destination ;

eCost .. sum[(i,j), pC(i,j) * vX(i,j)] =e= vCost ;

eCapacity(i) .. sum[ j , vX(i,j)] =l= pA(i) ;
eDemand (j) .. sum[ i , vX(i,j)] =g= pB(j) ;

model mTransport / all /

solve mTransport using LP minimizing vCost

Good Optimization Modeling Practices with GAMS. July 2016 27

General structure of GAMS sentences
• Commenting
• Lines with * in the first column
• $OnText $OffText to comment many lines

• No distinction between uppercase and lowercase letters

• Parenthesis (), square bracket [] or braces {} can be used
indistinctly to distinguish levels
• Language reserved words appear in bold
• Sentences end with a ;
• Can be suppressed when the following word is a reserved one

Good Optimization Modeling Practices with GAMS. July 2016 28

Basic input/output in text format
• Data input from a text file
$include FileName.txt

• display IdentifierName (shows its content or value)

• Data output to a text file

file InternalName / ExternalName.txt /
put InternalName
put IdentifierName
putclose InternalName

• Specific options to control the output format

• Put Writing Facility

Good Optimization Modeling Practices with GAMS. July 2016 29

Functions and operators
• +, -, *, /, ** or power(x,n)
• abs, arctan, sin, cos, ceil, floor, exp, log, log10, max, min, mod, round, sign, sqr, sqrt,
trunc, normal, uniform
• gyear, gmonth, gday, ghour, gminute, gsecond, gdow, gleap, jdate, jnow, jstart, jtime
• lt <, gt >, eq =, ne <>, le <=, ge >=
• not, and, or, xor
• diag(set_element,set_element)={1,0}
• sameas(set_element,set_element)={T,F}
• ord, card ordinal and cardinal of a set, SetName.pos ordinal of a set
• sum, prod, smax, smin
• inf, eps are valid as data

Model temporal license

abort $[jstart > jdate(2016,7,15)] 'License for this model has expired and it cannot be used any more. Contact the developers'

Good Optimization Modeling Practices with GAMS. July 2016 30

$ Operator in assignments, summations, constraints
• Sets a condition
$(value > 0) $(number1 <> number2)

• On the left of an assignment, it does the assignment ONLY if the

condition is satisfied
if (condition,
DO THE ASSIGNMENT
) ;
• On the right of an assignment, it does the assignment ALWAYS and if
the condition is not satisfied it assigns value 0
if (condition,
DO THE ASSIGNMENT
else
ASSIGNS VALUE 0
) ; Good Optimization Modeling Practices with GAMS. July 2016 31
Dynamic sets
• Efficiency is strongly related to the use of
dynamic sets
• Subsets of static sets whose content may
change by assignments
sets m months / m01 * m12 /
mp(m) even months
display m ;
mp(m) $[mod(ord(m),2) = 0] = yes ;
display mp ;
mp(‘m03') = yes ;
display mp ;
mp(m) $[ord(m) = 4] = no ;
display mp ;
According to legend Roland's
• Fundamental elements in developing GAMS Breach was cut by Count Roland
models with his sword Durendal in an
attempt to destroy that sword,
after being defeated during the
• Must be used systematically to avoid the Battle of Roncesvalles in 778.
formulation of superfluous equations,
variables or assignments
Good Optimization Modeling Practices with GAMS. July 2016 32
Use and abuse of dynamic sets
sets
w weeks / w01 * w52 /
h hours / h001 * h168 /
nd nodes / node01 * node99 /
ln(nd,nd) lines
spring( w) spring weeks / w13 * w25 /
days5 ( h) first 5 days of the week / h001 * h120 /
sprday(w,h) spring working days

alias (nd,ni,nf)

parameters
pDemand(w,h,nd ) demand in each node
pFlow (w,h,nd,nd) flow in each line ;

sprday(w,h) $[spring(w)*days5(h)] = yes ;

* these sentences are equivalent

pDemand(w,h,nd ) $[spring(w)*days5(h)] = uniform(-0.5,0.5) ;

pFlow (w,h,nd,nd ) $[spring(w)*days5(h)] = uniform(-1.0,1.0) ;

pDemand(w,h,nd ) $sprday(w,h) = uniform(-0.5,0.5) ;

pFlow (w,h,nd,nd ) $sprday(w,h) = uniform(-1.0,1.0) ;

pDemand(sprday,nd ) = uniform(-0.5,0.5) ;
pFlow (sprday,nd,nd) = uniform(-1.0,1.0) ;

Good Optimization Modeling Practices with GAMS. July 2016 33

Index shifting. Lag and lead
• t=J,F,M,A,MA,J,JU,AU,S,O,N,D
vReserve(t) + pInflow(t) - vOutflow(t) =e= vReserve(t+1)

• Vector values out of the domain are 0

vReserve(’D’) + pInflow(’D’) - vOutflow(’D’) =e= 0

• Circular sequence of an index (++, --)

t=J,F,M,A,MA,J,JU,AU,S,O,N,D
vReserve(t) + pInflow(t) - vOutflow(t) =e= vReserve(t++1)
vReserve(’D’) + pInflow(’D’) - vOutflow(’D’) =e= vReserve(’J’)

• Inverted order sequence of PP index even though t is traversed in

increasing order
PP(t+[card(t)-2*ord(t)+1])

Good Optimization Modeling Practices with GAMS. July 2016 34

Analytical computation of constraints and variables
• Important to know the
estimated size of the
optimization problem and
dependence with respect
to basic elements
• It can be used for detecting
formulation errors
• Use LimRow/LimCol
• Suitable to know the
constraint matrix structure
(GAMSChk)
option LP=GAMSChk

Good Optimization Modeling Practices with GAMS. July 2016 36

Source: MPES
Number of equations
and variables for MHE

Good Optimization Modeling Practices with GAMS. July 2016 37

How big is a big optimization problem
• Memory requirements for creating the model (GAMS)
• 1 GB for every 1 million rows
• Memory requirements for solving the model (solver)
• Depends on the difficulty in solving the model
• Integrality gap for MIP problems

Good Optimization Modeling Practices with GAMS. July 2016 38

Avoid creation of superfluous constraints and variables

• Or how to achieve a compact formulation (small size of the

constraint matrix or small density)
• Some redundant constraints can introduce a tighter model, see later
• However, introduce logical conditions (with a $ in GAMS) in the
creation of equations or the use of variables to avoid superfluous
ones
• Reduction rules: mathematical reasoning or common sense based
on the problem context
• Flows by non existing connections in a network
• Solvers can detect some of these superfluous equations/variables
but it is more efficient to avoid their creation (pre-processing)
• Profile, ProfileTol

Good Optimization Modeling Practices with GAMS. July 2016 39

Debugging an optimization model
• Grammar error
• Read the error and click in the red line of the log file

• Infeasibility detection
• Soft (elastic) constraints
• Introduce a deficit or surplus variable in each equation and penalize them in
the objective function
• Detect the smallest core of infeasible constraints by the LP solver
(option Irreducible Infeasible Subsets iis)
• Once known, they have to be deleted or modified

Good Optimization Modeling Practices with GAMS. July 2016 40

LP Performance issues and their suggested
resolution

E. Klotz, A.M. Newman Practical guidelines for solving difficult linear programs Surveys in Operations
Research and Management Science 18 (1-2), 1-17, Oct 2013 10.1016/j.sorms.2012.12.001
Good Optimization Modeling Practices with GAMS. July 2016 41
Efficiency in GAMS code usage (loop)

set i / 1*2000 / set i / 1*2000 /

alias (i,ii) alias (i,ii)
parameter X(i,i) parameter X(i,i) ;
loop ((i,ii),
X(i,ii) = 4 ;
) ; X(i,ii) = 4 ;

75.3 s 0.3 s

Good Optimization Modeling Practices with GAMS. July 2016 43

Observer effect
• Changes that the act of observation will make
on a phenomenon being observed

option Profile=10, ProfileTol=0.01

set i / 1*2000 /
set i / 1*2000 / alias (i,ii)
alias (i,ii) parameter X(i,i)
parameter X(i,i)
loop ((i,ii),
X(i,ii) = 4 ;
loop ((i,ii),
X(i,ii) = 4 ; ) ;
) ;

74.2 s 72.9 s

Good Optimization Modeling Practices with GAMS. July 2016 44

Efficiency in GAMS code usage (index order)

option Profile=10, ProfileTol=0.01 option Profile=10, ProfileTol=0.01

set i / 1*200 / set i / 1*200 /

j / 1*200 / j / 1*200 /
k / 1*200 / k / 1*200 /
parameter X(k,j,i), Y(i,j,k) ; parameter X(i,j,k), Y(i,j,k) ;

Y(i,j,k) = 2 ; Y(i,j,k) = 2 ;

X(k,j,i) = Y(i,j,k) X(i,j,k) = Y(i,j,k)

4.5 s 1.3 s

Good Optimization Modeling Practices with GAMS. July 2016 45

These constructs also exist in GAMS
loop (set,
) ;

while (condition,
) ;

repeat
until condition;

if (condition,
else
) ;

for (i=beginning to/downto end by increment,

) ;

Good Optimization Modeling Practices with GAMS. July 2016 46

 $Phantom null

Detection of sets
nd
ndref (nd)
refnd (nd)
nodes
current reference
subset of reference
node
nodes
/
/
/
node01
node01
node01
* node19 /
/
/

isolated nc
nod
ln
(nd)
(nd)
(nd,nd)
subnet(nd,nd,nd)
current connected
subset of connected
lines
subnetworks
nodes
nodes
/
/
null
null
/
/

subnetworks parameters
pAux1 auxiliary / 0 /
pAux2 auxiliary / 1 /

alias (nd,n1,n2,ni,nf)

file out / out.gms / put out ;

* create a naïve network, a chain
ln(ni,nf) $[ni.pos = nf.pos-1] = yes ;

* break these links

ln('node10','node11') = no ;
ln('node15','node16') = no ;

* detection of isolated subnetworks

* for every subnetwork => max number of iterations
loop (n1 $[sum(nod(nd), 1) < card(nd)],
* define the reference node for 2nd+ iterations
ndref(nd) $[n1.pos > 1 and not refnd(nd) and nd.pos = smin(n2 $[not nod(n2)], n2.pos)] = yes ;
* empty the set of connected nodes
nc(nd ) = no ;
* connect the reference node
nc(ndref) = yes ;
pAux2 = 1 ;

* for every node => max number of iterations

loop (n2 $pAux2,
* count the number of connected nodes
pAux1 = sum[nc, 1] ;
* add nodes to the set of already connected nodes
nc(nf) $ sum[nc $[ln(nc,nf) or ln(nf,nc)], 1] = yes ;
* count the new added nodes
pAux2 $[sum[nc, 1] - pAux1 = 0] = 0 ;

if (pAux2 = 0,
* subnetwork of the connected lines to a reference node
subnet(ln(nc,nf),ndref) = yes ;
subnet(ln(nf,nc),ndref) = yes ;
* subnetwork of the connected nodes
nod ( nc ) = yes ;
* subnetwork of the reference nodes
refnd ( ndref) = yes ;
) ;

) ;
display subnet, nod, refnd ;
* disconnect the reference node
ndref(nd) = no ;Good Optimization Modeling Practices with GAMS. July 2016 47
) ;
 Inverting a matrix, e.g., PTDF
set i / i1*i3 /

table a(i,i) matrix to invert

i1 i2 i3
i1 1
i2 3
i3 5

parameter ainv(i,i) inverted matrix

execute_unload ' GDXForInverse.gdx' i a

execute 'invert GDXForInverse.gdx i a GDXFromInverse.gdx ainv'
execute_load ' GDXFromInverse.gdx' ainv
execute 'del GDXForInverse.gdx GDXFromInverse.gdx'

display a, ainv

* computation of the susceptance matrix of the corridors

pYBUS(c2) = - sum[la(c2,cc), 1/pLineX(la)] ;
pYBUS(nf,ni) $c2(ni,nf) = pYBUS(ni,nf) ;
pYBUS(nf,nf) = - sum[ni, pYBUS(ni,nf)] ;
pYBUS(ni,ndref(nd)) = 0 ;
pYBUS(ndref(nd),nf) = 0 ;
pYBUS(ni,nd) $[not nc(ni)] = 0 ;
pYBUS(nd,nf) $[not nc(nf)] = 0 ;

* obtaining the inverse of pYBUS and saving it into pYBUSInv

execute_unload ' GDXForInverse.gdx' noref pYBUS
execute 'invert GDXForInverse.gdx noref pYBUS GDXFromInverse.gdx pYBUSInv'
execute_load ' GDXFromInverse.gdx' pYBUSInv
execute 'del GDXForInverse.gdx GDXFromInverse.gdx' ;

* computation of the PTDF matrix

Good
pPTDF(la(ni,nf,cc),ngd) Optimization Modeling
= [pYBUSInv(ni,ngd) 48 ;
Practices with GAMS. July 2016 + eps
- pYBUSInv(nf,ngd)]/pLineX(la)
 • Split formulation from data.
• Protect formulation confidentiality
Deployment ready • Secure Work Files
• Control the access to symbol names
• Link the model to a specific license

• Set declaration. Initialization.

• Variables
• Equations
• Model

• Include and manipulate input data and parameters. Prepared

to be
deployed
• Bounds and initialization of variables
• Solve the optimization problem
• Output of the results
Good Optimization Modeling Practices with GAMS. July 2016 49
Transportation model (ready for deployment)
sets Formulation.gms
I origins $include Data.gms Remaining.gms
J destinations
solve mTransport using LP minimizing vCost
parameters
pA(i ) origin capacity
pB(j ) destination demand
pC(i,j) per unit transportation cost
sets Data.gms
variables I origins / VIGO, ALGECIRAS /
vX(i,j) units transported J destinations / MADRID, BARCELONA, VALENCIA /
vCost transportation cost
parameters
positive variable vX pA(i) origin capacity
/ VIGO 350
equations
ALGECIRAS 700 /
eCost transportation cost
eCapacity(i) maximum capacity of each origin
eDemand (j) demand supply at destination ; pB(j) destination demand
/ MADRID 400
eCost .. sum[(i,j), pC(i,j) * vX(i,j)] =e= vCost ; BARCELONA 450
eCapacity(i) .. sum[j, vX(i,j)] =l= pA(i) ; VALENCIA 150 /
eDemand (j) .. sum[i, vX(i,j)] =g= pB(j) ;
table pC(i,j) per unit transportation cost
model mTransport / all / MADRID BARCELONA VALENCIA
VIGO 0.06 0.12 0.09
ALGECIRAS 0.05 0.15 0.11

Generate the runtime model Execute runtime model + data

gams Formulation.gms Save=model Send model.g00 gams Remaining.gms Restart=model

Good Optimization Modeling Practices with GAMS. July 2016 50

Model log
• Open console from GAMSIde for logging messages from
the model
• Code specific for Windows, UNIX/Linux/macOS
$set console
$if '%system.filesys%' == 'MSNT' $set console con
$if '%system.filesys%' == 'UNIX' $set console /dev/tty
$if '%console%.' == '.' abort 'console not recognized'

file console / '%console%' /

sets
day day / day01*day10 /
sc scenario / sc01* sc02 /

put console
loop ((day,sc),
putclose 'Day ' day.tl:0 ' Scenario ' sc.tl:0 ' Elapsed Time ' [(jnow-jstart)*86400]:6:3 ' s' sleep(1)
) ;

$ifthen.MSNT '%system.filesys%' == 'MSNT'

execute 'del pp.txt' ;
$endif.MSNT
$ifthen.UNIX '%system.filesys%' == 'UNIX'
execute 'rm pp.txt' ;
$endif.UNIX Good Optimization Modeling Practices with GAMS. July 2016 51
GAMS Code Conventions
• Must be defined in blocks. For example, a set and all its subsets should constitute one block in
the sets section.
• Names are intended to be meaningful. Follow conventions
• Items with the same name represent the same concept in different models
• Units should be used in all definitions
• Parameters are named pParameterName (e.g., pTotalDemand)
• Variables are named vVariableName (e.g., vThermalOutput)
• Equations are named eEquationName (e.g., eLoadBalance)
• Use short set names (one or two letters) for easier reading
• Alias duplicate the final letter (e.g., p, pp)
• In the case of variables, the blocks should be defined by meaning and not by variable type
(Free (default), Positive, Negative, Binary, Integer, SOS1, SOS2, SemiCont, SemiInt)
• Objective function must be a free variable
• Equations are laid out as clearly as possible, using brackets for readability

Good Optimization Modeling Practices with GAMS. July 2016 52

Stochastic Daily Unit Commitment of Thermal and Hydro
Units (SDUC) (i) (http://iit.upcomillas.es/aramos/StarGenLite_SDUC.zip)
$Title StarGen Lite Stochastic Daily Unit Commitment of Thermal and Hydro Units (SDUC)

$OnText

Developed by

Andrés Ramos
Instituto de Investigación Tecnológica
Escuela Técnica Superior de Ingeniería - ICAI
UNIVERSIDAD PONTIFICIA COMILLAS
Alberto Aguilera 23
28015 Madrid, Spain
Andres.Ramos@comillas.edu

March 15, 2016

$OffText

$OnEmpty OnMulti OffListing

* solve the optimization problems until optimality

option OptcR = 0

Good Optimization Modeling Practices with GAMS. July 2016 53

Stochastic Daily Unit Commitment of Thermal and
Hydro Units (SDUC) (ii)
* definitions

sets
n hour
n1(n) first hour of the day
sc scenario

g generating unit
t (g) thermal unit
h (g) hydro plant

alias (n,nn)

parameters
pDemand (n) hourly load [GW]
pOperReserve (n) hourly operating reserve [GW]
pOperReserveUp (n) hourly operating reserve up [GW]
pOperReserveDw (n) hourly operating reserve down [GW]
pIntermGen (n,sc) stochastic IG generation [GW]
pScenProb (sc) probability of scenarios [p.u.]
pCommitt ( g,n) commitment of the unit [0-1]
pProduct (sc,g,n) output of the unit [GW]
pSRMC (sc, n) short run marginal cost [ € per MWh]

pMaxProd (g) maximum output [GW]

pMinProd (g) minimum output [GW]
pMaxCons (g) maximum consumption [GW]
pIniOut (g) initial output > min load [GW]
pIniUC (g) initial commitment [0-1]
pRampUp (g) ramp up [GW per h]
pRampDw (g) ramp down [GW per h]
pMinTU (g) minimum time down [h]
pMinTD (g) minimum time up [h]
pSlopeVarCost (g) slope variable cost [M€ per GWh]
pInterVarCost (g) intercept variable cost [M€ per h]
pEmissionCost (g) emission cost [M€ per GWh]
pStartupCost (g) startup cost [M€]
pShutdownCost (g) shutdown cost [M€]
pMaxReserve (g) maximum reserve [GWh]
pMinReserve (g) minimum reserve [GWh]
pIniReserve (g) initial reserve [GWh]
pEffic (g) pumping efficiency [p.u.]
pInflows (g,n) inflows [GWh]
pENSCost energy not served cost [M€ per GWh]
pCO2Cost CO2 emission cost [ € per tCO2]

Good Optimization Modeling Practices with GAMS. July 2016 54

Stochastic Daily Unit Commitment of Thermal and
Hydro Units (SDUC) (iii)
variables
vTotalVCost total system variable cost [M€]

binary variables
vCommitt ( n,g) commitment of the unit [0-1]
vStartup ( n,g) startup of the unit [0-1]
vShutdown ( n,g) shutdown of the unit [0-1]

positive variables
vProduct (sc,n,g) output of the unit [GW]
vProduct1 (sc,n,g) output of the unit > min load [GW]
vConsump (sc,n,g) consumption of the unit [GW]
vENS (sc,n ) energy not served [GW]
vIG (sc,n ) intermittent generation [GW]
vReserve (sc,n,g) reserve at the end of period [GWh]
vSpillage (sc,n,g) spillage [GWh]

equations
eTotalVCost total system variable cost [M€]
eBalance (sc,n ) load generation balance [GW]
eOpReserve( n ) operating reserve [GW]
eReserveUp(sc,n ) operating reserve upwards [GW]
eReserveDw(sc,n ) operating reserve downwards [GW]
eMaxOutput(sc,n,g) max output of a committed unit [GW]
eMinOutput(sc,n,g) min output of a committed unit [GW]
eTotOutput(sc,n,g) tot output of a committed unit [GW]
eRampUp (sc,n,g) bound on ramp up [GW]
eRampDw (sc,n,g) bound on ramp down [GW]
eUCStrShut( n,g) relation among commitment startup and shutdown
eWtReserve(sc,n,g) water reserve [GWh]
eMinTUp (sc,n,g) minimum up time ( committed)
eMinTDown (sc,n,g) minimum down time (not committed) ;

Good Optimization Modeling Practices with GAMS. July 2016 55

 Stochastic Daily Unit Commitment of Thermal and
Hydro Units (SDUC) (iv)
* mathematical formulation

eTotalVCost .. vTotalVCost =e= sum[(sc,n ), pENSCost *vENS (sc,n )*pScenProb(sc)] +

sum[(sc,n,t), pSlopeVarCost(t)*vProduct (sc,n,t)*pScenProb(sc)] +
sum[(sc,n,t), pEmissionCost(t)*vProduct (sc,n,t)*pScenProb(sc)] +
sum[( n,t), pInterVarCost(t)*vCommitt ( n,t)] +
sum[( n,t), pStartupCost (t)*vStartup ( n,t)] +
sum[( n,t), pShutdownCost(t)*vShutdown( n,t)] ;

eBalance (sc,n ) $ pScenProb(sc) .. sum[t, vProduct(sc,n,t)] + sum[h, vProduct(sc,n,h)] - sum[h, vConsump(sc,n,h)/pEffic(h)] + vIG(sc,n) + vENS(sc,n) =e= pDemand(n) ;

eOpReserve( n ) .. sum[t, pMaxProd(t) * vCommitt(n,t)] + sum[h, pMaxProd(h)] =g= pOperReserve (n) + pDemand(n) ;
eReserveUp(sc,n ) $ pScenProb(sc) .. sum[t, pMaxProd(t) * vCommitt(n,t) - vProduct(sc,n,t)] =g= pOperReserveUp(n) ;
eReserveDw(sc,n ) $ pScenProb(sc) .. sum[t, pMinProd(t) * vCommitt(n,t) - vProduct(sc,n,t)] =l= - pOperReserveDw(n) ;

eMaxOutput(sc,n,t) $[pScenProb(sc) and pMaxProd(t)] .. vProduct(sc,n,t) / pMaxProd(t) =l= vCommitt(n,t) ;

eMinOutput(sc,n,t) $[pScenProb(sc) and pMinProd(t)] .. vProduct(sc,n,t) / pMinProd(t) =g= vCommitt(n,t) ;

eTotOutput(sc,n,t) $ pScenProb(sc) .. vProduct(sc,n,t) =e= pMinProd(t)*vCommitt(n,t) + vProduct1(sc,n,t) ;

eRampUp (sc,n,t) $ pScenProb(sc) .. vProduct1(sc,n,t) - vProduct1(sc,n-1,t) - pIniOut(t) $n1(n) =l= pRampUp(t) ;

eRampDw (sc,n,t) $ pScenProb(sc) .. vProduct1(sc,n,t) - vProduct1(sc,n-1,t) - pIniOut(t) $n1(n) =g= - pRampDw(t) ;

eUCStrShut( n,t) .. vCommitt(n,t) - vCommitt(n-1,t) - pIniUC(t) $n1(n) =e= vStartup(n,t) - vShutdown(n,t) ;

eMinTUp (sc,n,t) $ pScenProb(sc) .. sum[nn $(ord(nn) >= ord(n)+1-pMinTU(t) and ord(nn) <= ord(n)), vStartup (nn,t)] =l= vCommitt(n,t) ;
eMinTDown (sc,n,t) $ pScenProb(sc) .. sum[nn $(ord(nn) >= ord(n)+1-pMinTD(t) and ord(nn) <= ord(n)), vShutdown(nn,t)] =l= 1 - vCommitt(n,t) ;

eWtReserve(sc,n,h) $ pScenProb(sc) .. vReserve(sc,n-1,h) + pIniReserve(h) $n1(n) + pInflows(h,n) - vSpillage(sc,n,h) - vProduct(sc,n,h) + vConsump(sc,n,h) =e=
vReserve(sc,n,h) ;

model SDUC / all / ;

SDUC.OptFile = 1 ; SDUC.solprint = 2 ; SDUC.holdfixed = 1 ;

Good Optimization Modeling Practices with GAMS. July 2016 56

Stochastic Daily Unit Commitment of Thermal and
Hydro Units (SDUC) (v)
* read input data from Excel and include into the model

file TMP / tmp.txt /

$OnEcho > tmp.txt
r1= indices
o1=tmp_indices.txt
r2= param
o2=tmp_param.txt
r3= demand
o3=tmp_demand.txt
r4= oprres
o4=tmp_oprres.txt
r5= oprresup
o5=tmp_oprresup.txt
r6= oprresdw
o6=tmp_oprresdw.txt
r7= IGgen
o7=tmp_IGgen.txt
r8= thermalgen
o8=tmp_thermalgen.txt
r9= hydrogen
o9=tmp_hydrogen.txt
r10= inflows
o10=tmp_inflows.txt
$OffEcho
$ call xls2gms m i="%gams.user1%.xlsm“ @"tmp.txt"

sets
$include tmp_indices.txt
;
$include tmp_param.txt
parameter pDemand(n) hourly load [MW] /
$include tmp_demand.txt
/
parameter pOperReserve(n) hourly operating reserve [MW] /
$include tmp_oprres.txt
/
parameter pOperReserveUp(n) hourly operating reserve [MW] /
$include tmp_oprresup.txt
/
parameter pOperReserveDw(n) hourly operating reserve [MW] /
$include tmp_oprresdw.txt
/
table pIntermGen(n,sc) stochastic IG generation [MW]
$include tmp_IGgen.txt
table pThermalGen(g,*)
$include tmp_thermalgen.txt
table pHydroGen (g,*)
$include tmp_hydrogen.txt
table pInflows (g,n)
$include tmp_inflows.txt

execute 'del tmp.txt tmp_indices.txt tmp_param.txt tmp_demand.txt tmp_oprres.txt tmp_oprresup.txt tmp_oprresdw.txt tmp_IGgen.txt tmp_thermalgen.txt tmp_hydrogen.txt tmp_inflows.txt' ;

Good Optimization Modeling Practices with GAMS. July 2016 57

Stochastic Daily Unit Commitment of Thermal and
Hydro Units (SDUC) (vi)
* determine the first hour of the day

n1(n) $[ord(n) = 1] = yes ;

* assignment of thermal units, storage hydro and pumped storage hydro plants

t (g) $[pThermalGen(g,'MaxProd' ) and pThermalGen(g,'FuelCost')] = yes ;

h (g) $[pHydroGen (g,'MaxProd' ) ] = yes ;

* scaling of parameters to GW and M€

pDemand (n ) = pDemand (n ) * 1e-3 ;

pOperReserve (n ) = pOperReserve (n ) * 1e-3 ;
pOperReserveUp(n ) = pOperReserveUp(n ) * 1e-3 ;
pOperReserveDw(n ) = pOperReserveDw(n ) * 1e-3 ;
pIntermGen (n,sc) $pScenProb(sc) = pIntermGen (n,sc) * 1e-3 ;

pENSCost = pENSCost * 1e-3 ;

pMaxProd (t) = pThermalGen(t,'MaxProd' ) * 1e-3 ;
pMinProd (t) = pThermalGen(t,'MinProd' ) * 1e-3 ;
pIniOut (t) = pThermalGen(t,'IniProd' ) * 1e-3 ;
pRampUp (t) = pThermalGen(t,'RampUp' ) * 1e-3 ;
pRampDw (t) = pThermalGen(t,'RampDw' ) * 1e-3 ;
pMinTU (t) = pThermalGen(t,'MinTU' ) ;
pMinTD (t) = pThermalGen(t,'MinTD' ) ;
pSlopeVarCost(t) = pThermalGen(t,'OMVarCost' ) * 1e-3 +
pThermalGen(t,'SlopeVarCost') * 1e-3 * pThermalGen(t,'FuelCost') ;
pEmissionCost(t) = pThermalGen(t,'EmissionRate') * 1e-3 * pCO2Cost ;
pInterVarCost(t) = pThermalGen(t,'InterVarCost') * 1e-6 * pThermalGen(t,'FuelCost') ;
pStartupCost (t) = pThermalGen(t,'StartupCost' ) * 1e-6 * pThermalGen(t,'FuelCost') ;
pShutdownCost(t) = pThermalGen(t,'ShutdownCost') * 1e-6 * pThermalGen(t,'FuelCost') ;

pMaxProd (h) = pHydroGen (h,'MaxProd' ) * 1e-3 ;

pMinProd (h) = pHydroGen (h,'MinProd' ) * 1e-3 ;
pMaxCons (h) = pHydroGen (h,'MaxCons' ) * 1e-3 ;
pEffic (h) = pHydroGen (h,'Efficiency' ) ;
pMaxReserve (h) = pHydroGen (h,'MaxReserve' ) * 1e-3 ;
pMinReserve (h) = pHydroGen (h,'MinReserve' ) * 1e-3 ;
pIniReserve (h) = pHydroGen (h,'IniReserve' ) * 1e-3 ;

* if the initial output of the unit is above its minimum load then the unit is committed, otherwise it is not committed
pIniUC (g) = 1 $[pIniOut(g) >= pMinProd(g)] ;

* if the efficiency of a hydro plant is 0, it is changed to 1

pEffic (h) $[pEffic (h) = 0] = 1 ;

* if the minimum up or down times are 0, they are changed to 1

pMinTU (t) $[pMinTU (t) = 0] = 1 ;
pMinTD (t) $[pMinTD (t) = 0] = 1 ; Good Optimization Modeling Practices with GAMS. July 2016 58
Stochastic Daily Unit Commitment of
Thermal and Hydro Units (SDUC) (vii)
* bounds on variables

vProduct.up (sc,n,g) $pScenProb(sc) = pMaxProd (g ) ;

vConsump.up (sc,n,g) $pScenProb(sc) = pMaxCons (g ) ;
vProduct1.up(sc,n,t) $pScenProb(sc) = pMaxProd (t ) - pMinProd(t) ;
vIG.up (sc,n ) $pScenProb(sc) = pIntermGen (n,sc) ;
vENS.up (sc,n ) $pScenProb(sc) = pDemand (n ) ;
vReserve.up (sc,n,g) $pScenProb(sc) = pMaxReserve(g ) ;
vReserve.lo (sc,n,g) $pScenProb(sc) = pMinReserve(g ) ;

vCommitt.up ( n,g) = 1 ;
vStartup.up ( n,g) = 1 ;
vShutdown.up( n,g) = 1 ;

* solve stochastic daily unit commitment model

solve SDUC using MIP minimizing vTotalVCost ;

* scaling of the results

pCommitt( t,n) = vCommitt.l( n,t) + eps ;

pProduct(sc,g,n) $pScenProb(sc) = vProduct.l(sc,n,g)*1e3 + eps ;
pSRMC (sc, n) $pScenProb(sc) = eBalance.m(sc,n )*1e3 + eps ;

* data output to xls file

put TMP putclose 'par=pCommitt rdim=1 rng=UC!a1' / 'par=pProduct rdim=2 rng=Output!a1' / 'par=pSRMC rdim=1 rng=SRMC!a1'
execute_unload 'tmp.gdx' pProduct pCommitt pSRMC
execute 'gdxxrw.exe tmp.gdx SQ=n EpsOut=0 O="tmp.xlsx" @tmp.txt'
execute 'del tmp.gdx '
execute 'del tmp.txt'

$OnListing

Good Optimization Modeling Practices with GAMS. July 2016 59

Scaling
• Solvers are powerful but not magic
• Input data and output results must be in commonly used units
• But internally variables, equations, parameters must be around 1 (i.e., from 0.01 to 100)
• Scaling can be done:
• Manually (e.g., from MW to GW, from euros to Meuros)
• Automatically by the language (ModelName.ScaleOpt=1) or
• By the solver (ScaInd 1 in CPLEX)

• Specially useful in large-scale LP problems or NLP problems and/or when willing to use the dual
variables
• The condition number is a measure for the sensitivity of the solution of a system of linear equations
to errors in the data.
• It is the ratio between largest and smallest eigenvalues
• Condition numbers below 106 are good enough. Numerical problems arise for condition numbers
above 1010
• Quality 1 in CPLEX
• Kappa 1 in GUROBI
Good Optimization Modeling Practices with GAMS. July 2016 60
Rule of thumb for selecting an LP algorithm
• Simplex (or dual simplex) method can the best choice for
moderate size (up to 100000 x 100000)
• Interior point method is usually the most efficient for
huge problems
• It can be threaded quite efficiently (compared to simplex)

• Difference in solution time can reach

one order of magnitude

Good Optimization Modeling Practices with GAMS. July 2016 62

Good Optimization Modeling Practices with GAMS. July 2016 63
Options
Options
LimRow Number of rows to show
LimCol Number of columns to show
SolPrint Solution output
SolveOpt Replace
Decimals Number of decimals in displaying values
IterLim Maximum number of solver iterations
ResLim Maximum solution time
Profile Time profiling
ProfileTol Profile threshold
Seed Initialize seed for random numbers

Good Optimization Modeling Practices with GAMS. July 2016 64

$ Directives
$ Directives
$OnEmpty Allow introduction of empty sets
$OnMulti Allow redeclaration of sets
$OffListing Suppress listing of the code

Good Optimization Modeling Practices with GAMS. July 2016 65

Variable attributes
Attributes
VarName.lo lower bound
VarName.up upper bound
VarName.fx fixes the variable to a constant
VarName.Range range of the variable
VarName.l initial value before and optimal value after
VarName.m marginal value (reduced cost)
VarName.Scale numerical scale factor
VarName.Prior branching priority in a MIP model (∞ →not discrete)
VarName.SlackUp slack from upper bound
VarName.SlackLo slack from lower bound
VarName.Infeas infeasibility out of bounds

Good Optimization Modeling Practices with GAMS. July 2016 66

Equation attributes
Attributes
EquationName.lo lower bound
EquationName.up upper bound
EquationName.l initial value before and optimal value after
EquationName.m marginal value (dual variable or shadow price)
EquationName.Scale numerical scaling factor

Good Optimization Modeling Practices with GAMS. July 2016 67

Model
model ModelName / EquationNames /
model ModelName / all /
model ModelName / all – EquationName1 + EquationName2 /

Good Optimization Modeling Practices with GAMS. July 2016 68

 Model attributes
Attributes
ModelName.ResLim Resource limit
ModelName.SolveOpt Replace/merge/clear in consecutive solves
ModelName.SolSlack Show slack variables
ModelName.SolvePrint 0, 1, 2 (to remove the detailed solution from the .lst file)
ModelName.TryLinear Try linear model first
ModelName.ModelStat Model status
Attributes
ModelName.SolveStat Solve status
ModelName.IterUsd Number of iterations
ModelName.NumEqu Number of equations
ModelName.NumVar Number of variables ModelName.ResUsd Resource used

ModelName.NumDVar Number of discrete variables ModelName.BRatio Basis ration controls the use of previous basis
ModelName.NumNz Number of non zeros ModelName.HoldFixed Fix and eliminate variables
ModelName.IterLim Iteration limit
ModelName.NodLim Node limit
ModelName.OptCA Absolute optimality tolerance
ModelName.OptCR Relative optimality tolerance
ModelName.OptFile Use of an option file
Good Optimization Modeling Practices with GAMS. July 2016 69
ModelName.PriorOpt Use of priority
GAMS Call Options
GAMS Options
Suppress Suppress echo of the code listing
PW Page width
PS Page size
Charset Allows international characters
For example, InterfaceName, SolverSelection,
U1..U10 User parameter SkipExcelInput, SkipExcelOutput,

u1="E50_Nodal_450_2_op_wodcpst_ptdf" u2=0 u3=0 u4=1 --NumberCores=4

Good Optimization Modeling Practices with GAMS. July 2016 70

Solvers

Good Optimization Modeling Practices with GAMS. July 2016 71

Solvers

Good Optimization Modeling Practices with GAMS. July 2016 72

Boosting performance
• Threads

• GUSS (Gather-Update-Solve-Scatter)
• Use of sensitivity analysis for solving many similar problems

• Grid and Multi-Threading Solve Facility

• Use of multiple cores of a computer

• You can launch several GAMS process simultaneously

being careful with conflicting filenames

Good Optimization Modeling Practices with GAMS. July 2016 73

 Scenario analysis of the transportation problem solved
with GUSS
set scen_dem stochastic demand scenario dictionary /
sc . scenario . ''
sets
* update the LHS with values of the RHS
I origins
J destinations pB . param . pBS
SC scenarios
* store in the RHS with values of the LHS
parameters vX . level . pX
pA (i ) origin capacity vCost . level . pCost
pB (j ) destination demand eDemand . marginal . pPrice /
pC (i,j ) per unit transportation cost
pBS (sc, j) stochastic destination demand *************************************************************
pX (sc,i,j) stochastic units transported
sets
pCost (sc ) stochastic transportation cost
pPrice(sc, j) stochastic spot price I origins / VIGO, ALGECIRAS /
J destinations / MADRID, BARCELONA, VALENCIA /
variables SC scenarios / sc000*sc999 /
vX(i,j) units transported
vCost transportation cost parameters
pA(i) origin capacity
positive variable vX / VIGO 350
ALGECIRAS 700 /
equations
eCost transportation cost
eCapacity(i) maximum capacity of each origin pBS(sc,j) stochastic destination demand
eDemand (j) demand supply at destination ; / sc000 . MADRID 400
sc000 . BARCELONA 450
eCost .. sum[(i,j), pC(i,j) * vX(i,j)] =e= vCost ; sc000 . VALENCIA 150 / ;
eCapacity(i) .. sum[j, vX(i,j)] =l= pA(i) ;
eDemand (j) .. sum[i, vX(i,j)] =g= pB(j) ; * lazy input, feeding data for all the scenarios with random demand
pBS(sc,j) = pBS('sc000',j) * [1+uniform(-0.05,0.05)] ;
model mTransport / all /
table pC(i,j) per unit transportation cost
MADRID BARCELONA VALENCIA
VIGO 0.06 0.12 0.09
ALGECIRAS 0.05 0.15 0.11 ;
*************************************************************

* initialization of the destination demand

pB(j) = pBS('sc000',j) ;

solve mTransport using LP minimizing vCost scenario scen_dem

Good Optimization Modeling Practices with GAMS. July 2016 74
Scenario analysis of the transportation problem solved
with Grid computing and GUSS (i)
sets
gs(sc) scenarios per GUSS run
sh solution headers / System.GUSSModelAttributes /
sets scen_dem stochastic demand scenario dictionary /
I origins sc . scenario . ''
J destinations scen_optn . opt . st_report_o
SC scenarios
* update the LHS with values of the RHS
parameters pB . param . pBS
pA (i ) origin capacity
pB (j ) destination demand * store in the RHS with values of the LHS
pC (i,j ) per unit transportation cost
vX . level . pX
pBS (sc, j) stochastic destination demand
pX (sc,i,j) stochastic units transported vCost . level . pCost
pCost (sc ) stochastic transportation cost eDemand . marginal . pPrice /
pPrice(sc, j) stochastic spot price
*************************************************************
variables sets
vX(i,j) units transported I origins / VIGO, ALGECIRAS /
vCost transportation cost
J destinations / MADRID, BARCELONA, VALENCIA /
positive variable vX SC scenarios / sc000*sc999 /

equations parameters
eCost transportation cost pA(i) origin capacity
eCapacity(i) maximum capacity of each origin / VIGO 350
eDemand (j) demand supply at destination ; ALGECIRAS 700 /
eCost .. sum[(i,j), pC(i,j) * vX(i,j)] =e= vCost ; pBS(sc,j) stochastic destination demand
eCapacity(i) .. sum[j, vX(i,j)] =l= pA(i) ;
/ sc000 . MADRID 400
eDemand (j) .. sum[i, vX(i,j)] =g= pB(j) ;
sc000 . BARCELONA 450
model mTransport / all / sc000 . VALENCIA 150 / ;

* lazy input, feeding data for all the scenarios with random demand
pBS(sc,j) = pBS('sc000',j) * [1+uniform(-0.05,0.05)] ;

table pC(i,j) per unit transportation cost

MADRID BARCELONA VALENCIA
VIGO 0.06 0.12 0.09
ALGECIRAS 0.05 0.15 0.11 ;

Good Optimization Modeling Practices with GAMS. July 2016 75

Scenario analysis of the transportation problem solved
with Grid computing and GUSS (ii)
sets
* using four cores and assignment of scenarios to cores
core grid jobs to run / core001*core004 /
coresc(core,sc) cores to scenario / core001.(sc000*sc249)
core002.(sc250*sc499)
core003.(sc500*sc749)
core004.(sc750*sc999) /
parameter
scen_optn scenario options / OptFile 2, LogOption 1, SkipBaseCase 1,
UpdateType 1, RestartType 1, NoMatchLimit 999 /
st_report_o(sc,sh) status report
pGridHandle(core) grid handles ;
* initialization of the destination demand
pB(j) = pBS('sc000',j) ;
mTransport.SolveLink = %SolveLink.AsyncGrid% ;
loop (core,
gs(sc) = coresc(core,sc)
if (sum[gs(sc), 1] > 0,
solve mTransport using LP minimizing vCost scenario scen_dem ;
pGridHandle(core) = mTransport.Handle ;
) ;
) ;
repeat
loop (core $HandleCollect(pGridHandle(core)),
display $HandleDelete (pGridHandle(core)) 'Trouble deleting handles' ;
pGridHandle(core) = 0 ;
) ;
until card(pGridHandle) = 0 or TimeElapsed > 1000 ;
mTransport.SolveLink = %SolveLink.LoadLibrary% ;
display st_report_o

Good Optimization Modeling Practices with GAMS. July 2016 76

How to write multiple language versions
* Language for Excel headings: Spanish 0 English 1 French 2

$set language 1

$ifthene %language% == 0

$ set iTitle MiModelo Versión 6.19 --- 11 Julio 2016

$ include ModelEs.gms

$elseif %language% == 1

$ set iTitle MyModel Release 6.19 --- July 11, 2016

$ include ModelEn.gms

$elseif %language% == 2

$ set iTitle MonModel Version 6.19 --- 11 Juillet 2016

$ include ModelFr.gms

$endif

* File ModelEs.gms
$setglobal Lema Más sabe el diablo por viejo que por diablo

* File ModelEn.gms
$setglobal Lema The devil knows many things because he is old

* File ModelFr.gms
$setglobal Lema Le diable sait beaucoup parce qu'il est vieux
Good Optimization Modeling Practices with GAMS. July 2016 77
Introducing flexibility
$SetGlobal x 100
parameter pDimension / %x% /
set u / unit1*unit%x% /
display pDimension, u

Good Optimization Modeling Practices with GAMS. July 2016 78

Releasing memory Clear Data
• Define a dummy solve
• Clear parameters
• Run dummy model
option profile=10

set i / 1 * 10000000 /
parameter pp(i) ;

pp(i) = 33 ;

* dummy optimization problem used for releasing memory

variable vDummy
equation eDummy ; eDummy .. vDummy =e= 0 ;
model mDummy / eDummy / ;

* only parameters that are no longer used can be cleared

option Clear=pp

* solve a dummy optimization problem to release memory usage

solve mDummy using LP minimizing vDummy
Good Optimization Modeling Practices with GAMS. July 2016 79
 GAMS to LaTeX
sets
I origins / VIGO, ALGECIRAS /
J destinations / MADRID, BARCELONA, VALENCIA / Generate the doc file
parameters
gams transport.gms DocFile=transport
pA(i) origin capacity
/ VIGO 350
Write the tex file
ALGECIRAS 700 / model2tex transport
pB(j) destination demand
/ MADRID 400
BARCELONA 450 min
VALENCIA 150 /

≤ ∀#
table pC(i,j) per unit transportation cost
MADRID BARCELONA VALENCIA
VIGO 0.06 0.12 0.09
ALGECIRAS 0.05 0.15 0.11

variables
≥ $ ∀%
vX(i,j) units transported
vCost transportation cost ≥0
positive variable vX

equations
eCost transportation cost
eCapacity(i) maximum capacity of each origin
eDemand (j) demand supply at destination ;

eCost .. sum[(i,j), pC(i,j) * vX(i,j)] =e= vCost ;

eCapacity(i) .. sum[ j , vX(i,j)] =l= pA(i) ;
eDemand (j) .. sum[ i , vX(i,j)] =g= pB(j) ;

model mTransport / all /

solve mTransport using LP minimizing vCost

Good Optimization Modeling Practices with GAMS. July 2016 80

Each tool has a specific purpose
• GDX (GAMS Data eXchange) utilities to interface with other applications
• Interfaces with other programs 1000 685

• Access (gdx2access, mdb2gms) 900

680

Excel (gdx2xls, xls2gms)

800

• 700
675
Perc 95 %

• SQL (gdx2sqlite, sql2gms)

Volumen últil o Vertido [hm3]

Perc 5 %
Mínimo
600
670 Media

• Matlab (gdxmrw)

Cota [m]
RMn
500 RMx
Vertido

• R (gdxrrw) 400
665
TunelInf
TunelSup
TnlMn

• APIs
300
660 TnlMx
Cota
200

• .Net 100
655

• Java 0 650

s0840
s0843
s0846
s0849
s0852
s0903
s0906
s0909
s0912
s0915
s0918
s0921
s0924
s0927
s0930
s0933
s0936
s0939
s0942
s0945
s0948
s0951
s1002
s1005
s1008
s1011
s1014
s1017
s1020
s1023
s1026
s1029
s1032
s1035
s1038
• Python
• Input/output data and simple graphs (Excel, Access)
• Advanced graphs (GNUPlot, Matlab)

Good Optimization Modeling Practices with GAMS. July 2016 81

 1. Programming Style
2. GAMS Code
3. Look and Feel
4. Mathematical Formulation
5. Advanced Algorithms
3

Look and Feel

Look and feel
In software design, look and feel is used in respect of a
graphical user interface and comprises aspects of its
design, including elements such as colors, shapes, layout,
and typefaces (the “look”), as well as the behavior of
dynamic elements such as buttons, boxes, and menus (the
“feel”) (http://en.wikipedia.org/wiki/Look_and_feel)

Good Optimization Modeling Practices with GAMS. July 2016 83

Interface (http://iit.upcomillas.es/aramos/StarGenLite_SDUC.zip)
• Microsoft Excel (with macros1) is the most popular application in business
environment for data manipulation and is used as interface to
• Introduce the input data,
• Run the model and
• Present the output results numerically and graphically
• Menu sheet contains buttons to
• Run the model
• Load the results
• Green sheets are used to introduce input data in named ranges
• Green cells can be modified by the user
• The other ones it is suggested not to modify them.

• Salmon sheets are used to present the output data

• Blue sheets are used to plot the results

1. Enable the use of macros when opening it in Excel

Good Optimization Modeling Practices with GAMS. July 2016 84
Input / Output

• Input
• Use simple input formats. Make input easy to prepare and adapted to the user.
• Check input for validity and plausibility. This can be done in input/output interface
(e.g., using conditional formatting in Excel) or in GAMS
• Identify bad input, recover if possible, if not abort the execution giving
information to the user
• Document your data layouts
• Definitions with physical dimensions
• References to the data origin
• Output
• It is difficult to determine the infeasibility of an optimization program. It is better
to check bad data.
• Produce self-explanatory output.

Good Optimization Modeling Practices with GAMS. July 2016 85

Input Sheets
• In each sheet one o several named ranges are included to
input the data.
• You can introduce or delete data but the named ranges
have to be kept. To avoid errors it is important to
insert/delete rows/columns out of the first or last
row/column, i.e., in any intermediate row/column.
• The range names can not be changed by the user, are
predefined by the GAMS model.
• The named range input in included directly as GAMS
input files so it must be grammatically correct.

Good Optimization Modeling Practices with GAMS. July 2016 86

Output Sheets
• Their names are predefined by the model and can not be
changed by the user.

Good Optimization Modeling Practices with GAMS. July 2016 87

Running the model
• From the Excel interface
• Write the GAMS directory
• Press the button Run Run

• Remember to look for errors in StarGenLite_XXX.lst

• Press the button Load results Load Results

• From the GAMSIde

• Write u1="StarGenLite_XXX" in the upper right window
• Run GAMS (F9)
• Remember to look for errors in StarGenLite_XXX.lst
• Results are written to Excel file tmp.xls
• Press the button Load results form the Excel interface

Good Optimization Modeling Practices with GAMS. July 2016 88

 1. Programming Style
2. GAMS Code
3. Look and Feel
4. Mathematical Formulation
5. Advanced Algorithms
4

Mathematical Formulation
Reformulation in MIP problems
• Most MIP problems can be formulated in different ways
• In MIP problems, a good formulation is crucial to solve
the model
• How good is a MIP formulation?
• Integrality gap difference between the objective function of the MIP and
LP relaxation solutions
• Given two equivalent MIP formulations, one is stronger
(tighter/better) than the other, if the feasible region of
the linear relaxation is strictly contained in the feasible
region of the other one. Integrality gap is lower.
Good Optimization Modeling Practices with GAMS. July 2016 90
Warehouse location problem (no limits) (i)
• Choose where to locate warehouses among a set of
locations and assign clients to the warehouses minimizing
the total cost. No limits means that there is no limit in the
number of clients assigned to a warehouse.
• Data j locations
i clients
c j localization cost in j
hij cost of satisfying the demand of client i from j

• Variables y =
1 warehouse located in j
j
0 otherwise
xij fraction of demand of client i met from j

Good Optimization Modeling Practices with GAMS. July 2016 91

Warehouse location problem (no limits) (ii)
Formulation #1 Formulation #2
min ∑c j y j + ∑hij xij min ∑c j y j + ∑hij xij
j ij
j ij

∑x ij = 1 ∀i ∑x
j
ij = 1 ∀i
j

xij ≤ y j ∀ij ∑x
i
ij ≤ My j ∀j

y j ∈{0,1} , xij ∈[ 0,1] y j ∈{0,1} , xij ∈[ 0,1]

Number of constraints: I+IJ Number of constraints: I+J

• Both formulations are MIP equivalent. However,

formulation #1 is stronger
• Intuitively the fewer constraints the better. That’s true in
LP. However, in many MIP problems the more constraints
the better. Good Optimization Modeling Practices with GAMS. July 2016 92
Production problem with fixed and inventory costs (i)
• Data t time period
ct fixed cost, pt variable cost, ht inventory cost
dt demand
1 to produce
• Variables yt = 
0 not produce
xt amount produced
st inventory at the end of the period
• Formulation #1
min ∑( ct yt + pt xt + ht st )
t

st −1 + xt = dt + st ∀t
Number of constraints: 2T
xt ≤ Myt ∀t
Number of variables: 3T
s0 = sT = 0
xt , st ≥ 0, yt ∈{0,1}

Good Optimization Modeling Practices with GAMS. July 2016 93

Production problem with fixed and inventory costs (ii)
• Variables
1 to produce
yt = 
0 not produce
qit quantity produced in period i to meet the demand in period t ≥ i

• Formulation #2
T t T
min ∑∑( pi + hi + hi+1 + ⋯+ ht −1 ) qit + ∑ct yt
t =1 i =1 t =1
t Number of constraints: T+T2/2
∑q
i =1
it = dt ∀t
Number of variables: T+T2/2
qit ≤ dt yi ∀it
qit ≥ 0, yt ∈{0,1}
• Formulation #2 is better. However, it has a greater
number of constraints and variables.

Good Optimization Modeling Practices with GAMS. July 2016 94

Some tips for MIP
• It can be interesting increase the number of variables if they
can be used in the branching strategy of B&B. For example,
“artificial” division of a zone in regions N, S, E and W to
branch first in these zonal regions.
• Alternatively, introduce lazy constraints (only in
GAMS/GUROBI)
• Avoid the use of big M parameters or put tight (lowest upper
bound) values for the big M
• GAMS/CPLEX/GUROBI supports the use of an indicator
constraint ≤ ()
*#+,) + - *#+,) + - Write in the file cplex.opt
≤ () ≤0 indic constraint$y 0
≥0 ≥0
) ∈ 0,1 ) ∈ 0,1

Good Optimization Modeling Practices with GAMS. July 2016 95

Tight and compact unit commitment
• G. Gentile, G. Morales-España and A. Ramos A Tight MIP Formulation of the Unit Commitment
Problem with Start-up and Shut-down Constraints EURO Journal on Computational Optimization
5 (1), 177–201 March 2017 10.1007/s13675-016-0066-y
• G. Morales-España, C.M. Correa-Posada, A. Ramos Tight and Compact MIP Formulation of
Configuration-Based Combined-Cycle Units IEEE Transactions on Power Systems 31 (2), 1350-
1359, March 2016 10.1109/TPWRS.2015.2425833
• G. Morales-España, J.M. Latorre, and A. Ramos Tight and Compact MILP Formulation for the
Thermal Unit Commitment Problem IEEE Transactions on Power Systems 28 (4): 4897–4908, Nov
2013 10.1109/TPWRS.2012.2222938
• G. Morales-España, J.M. Latorre, and A. Ramos Tight and Compact MILP Formulation of Start-Up
and Shut-Down Ramping in Unit Commitment IEEE Transactions on Power Systems 28 (2): 1288-
1296, May 2013 10.1109/TPWRS.2012.2222938

Good Optimization Modeling Practices with GAMS. July 2016 96

Reformulation of an NLP problem
n
n n n n min ∑ x iwi
min ∑ ∑q xx ij i j min ∑ x i ∑q x ij j
i =1
n
i =1 j =i +1 i =1 j =i +1
n n wi = ∑q x ij j

∑x j =1 ∑x
j =1
j =1
n
j =i +1

j =1
n n ∑x j =1
∑r x
j =1
j j = r0 ∑r x
j =1
j j = r0 j =1
n

∑r x
j =1
j j = r0
• Formulation #2 is better than the #1. The evaluation of the objective function in #1
requires 2n2/2 multiplications. In #2 only n+n2/2
• Formulation #3 has essentially the same number of multiplications but they appear in
linear constraints. Number of constraints is bigger but all of them are linear. Linear
algebra is much more efficient. Formulation #3 is the most efficient
Good Optimization Modeling Practices with GAMS. July 2016 97
Reformulation of an NLP problem
0
*#+
1
+) 0= +)
*#+
∑
1=

1≥2

• Formulation #1 has a lot of nonlinear variables and it not protected against division by
zero
• Formulation #2 has only 2 nonlinear variables, the remaining ones appear in linear
equations, and the denominator is lower bounded to avoid division by zero. Model is
easier to solve and more robust

Good Optimization Modeling Practices with GAMS. July 2016 98

 Product of two variables
=) −) Quadratic separable form

) =( + )/2
) =( − )/2
7 ≤ ≤0
7 ≤ ≤0
(7 + 7 ) ≤ ) ≤ (0 + 0 )
(7 − 0 ) ≤ ) ≤ (0 − 7 )
AIMMS Modeling Guide - Integer Programming Tricks
Good Optimization Modeling Practices with GAMS. July 2016 99
Solving real problems
• MIP
• Solve with a sensible relative optimality tolerance
• Provide an initial solution based on specific knowledge of the model
or use the solution from a previous solve

• NLP
• Introduce sensible bounds on variables AND
• Provide a good enough starting point AND
• Scale the problem

• Symmetry breaking constraint ≥ 8

Good Optimization Modeling Practices with GAMS. July 2016 100

CPLEX Performance Tuning for MIP
(http://www-01.ibm.com/support/docview.wss?uid=swg21400023)
Names no
Pure branch and bound
NodeFileInd 3
Cuts -1
NodeSel 0
HeurFreq -1
VarSel 3
StartAlg 4 Presolve
MemoryEmphasis 1 PreInd, PrePass
WorkMem 1000
MIPEmphasis 2 Solution method of LP problem
MIPSearch 2 First iteration (interior point or simplex method)
SolveFinal 0 Successive iterations (primal or dual simplex)
Solution Polishing Priority for variable selection
Solution pool Select variables that impact the most in the o.f. (e.g.,
FlowCovers investment vs. operation variables)
FeasOptMode 2 Initial cutoff or incumbent
FeasOpt 1 Initial valid bound of the o.f. estimated by the user
tuning cplex.opt
RINSHeur 100
FpHeur 2
Good Optimization Modeling Practices with GAMS. July 2016 101
GUROBI (http://www.gurobi.com/documentation/6.5/refman/mip_models.html)
• Most important parameters Heuristics
• Threads, MIPFocus Heuristics, SubMIPNodes,
• Solution Improvement MinRelNodes, PumpPasses,
• ImproveStartTime, ZeroObjNodes
ImproveStartGap RINS 100
• Termination Cutting Planes
• TimeLimit Cuts, FlowCoverCuts, MIRCuts
• MIPGap, MIPGapAbs Presolve
• NodeLimit, IterationLimit, Presolve, PrePasses, AggFill,
SolutionLimit PreSparsify
• Cutoff
• Speeding Up The Root Relaxation
• Method

Good Optimization Modeling Practices with GAMS. July 2016 102

 1. Programming Style
2. GAMS Code
3. Look and Feel
4. Mathematical Formulation
5. Advanced Algorithms
5

Advanced Algorithms
 Transportation problem solved as MCP
sets
I origins / VIGO, ALGECIRAS / min sets
J destinations / MADRID, BARCELONA, VALENCIA / I origins / VIGO, ALGECIRAS /
J destinations / MADRID, BARCELONA, VALENCIA /
parameters
pA(i) origin capacity ≤ ∀# parameters
pA(i) origin capacity
/ VIGO 350
ALGECIRAS 700 / / VIGO 350
ALGECIRAS 700 /
pB(j) destination demand ≥ $ ∀% pB(j) destination demand
/ MADRID 400 / MADRID 400
BARCELONA 450
≥0
BARCELONA 450
VALENCIA 150 / VALENCIA 150 /

table pC(i,j) per unit transportation cost table pC(i,j) per unit transportation cost
MADRID BARCELONA VALENCIA MADRID BARCELONA VALENCIA
VIGO 0.06 0.12 0.09 VIGO 0.06 0.12 0.09
ALGECIRAS 0.05 0.15 0.11 ℒ = +9 − +: $ − ALGECIRAS 0.05 0.15 0.11

variables variables
vX(i,j) units transported vX(i,j) units transported
vCost transportation cost vA(i ) Lagrange multiplier of capacity constraint
vB( j) Lagrange multiplier of demand constraint
=> +9 ≥ : : ∀#%
→
positive variable vX
positive variables vX, vA, vB
equations =?@A equations
eCost transportation cost
eCapacity(i) maximum capacity of each origin − ≥− :9 ∀# eProfit(i,j) marginal cost >= marginal profit
eCapacity(i) maximum capacity of each origin
eDemand (j) demand supply at destination ; eDemand (j) demand supply at destination ;

≥$ :: ∀%
eCost .. sum[(i,j), pC(i,j) * vX(i,j)] =e= vCost ; eProfit(i,j) .. vA(i) + pC(i,j) =g= vB(j) ;
eCapacity(i) .. sum[ j , vX(i,j)] =l= pA(i) ; eCapacity(i) .. -sum[j, vX(i,j)] =g= -pA(i) ;
eDemand (j) .. sum[ i , vX(i,j)] =g= pB(j) ; eDemand (j) .. sum[i, vX(i,j)] =g= pB(j) ;

model mTransport / all / ,9 ,: ≥ 0 model mTransport / eProfit.vX eCapacity.vA eDemand.vB /

solve mTransport using LP minimizing vCost solve mTransport using MCP

Good Optimization Modeling Practices with GAMS. July 2016 104

 Fixed-Charge Transportation Problem (FCTP)
• Bd Relaxed Master
Flows
(second stage) Complete problem min G + C )
B ,F
Investment decisions
(first stage) HIG + DI ( ) ≥ DI ( ) 7 = 1, … , K
min C ) +
,B ) ∈ 0,1
Capacity of each
origin ≤ ∀#
min
Demand of each
destination ≤$ ∀%
• Bd Subproblem ≤ ∀#
Flow can pass only
for installed
≤( ) ∀#%
≥ 0, ) ∈ 0,1
≤$ ∀%
connections

≤ ( )E ∀#% : DE
≥0
Good Optimization Modeling Practices with GAMS. July 2016 105
FCTP solved by Benders decomposition (i)
$title Fixed-charge transportation problem (FCTP) solved by Benders' decomposition
* relative optimality tolerance in solving MIP problems
option OptcR = 0
sets
L iterations / l1 * l20 /
LL(l) iterations subset
I origins / i1 * i4 /
J destinations / j1 * j3 /
* Begin problem data
parameters
A(i) product offer
/ i1 20, i2 30, i3 40, i4 20 /
B(j) product demand
/ j1 20, j2 50, j3 30 /
table C(i,j) per unit variable transportation cost
j1 j2 j3
i1 1 2 3
i2 3 2 1
i3 2 3 4
i4 4 3 2
table F(i,j) fixed transportation cost
j1 j2 j3
i1 10 20 30
i2 20 30 40
i3 30 40 50
i4 40 50 60
* End problem data
abort $(sum[i, A(i)] < sum[j, B(j)]) 'Infeasible problem'
parameters
BdTol relative Benders' tolerance / 1e-6 /
Z_Lower lower bound / -inf /
Z_Upper upper bound / inf /
Y_L (l,i,j) first stage variables values in iteration l
PI_L (l,i,j) dual variables of second stage constraints in iteration l
Delta(l) cut type (feasibility 0 optimality 1) in iteration l
Z2_L (l) subproblem objective function value in iteration l
positive variable
X(i,j) arc flow
binary variable
Y(i,j) arc investment decision
variables
Z1 first stage objective function
Z2 second stage objective function
Theta recourse function
Good Optimization Modeling Practices with GAMS. July 2016 106
FCTP solved by Benders decomposition (ii)
equations
EQ_Z1 first stage objective function
EQ_Z2 second stage objective function
EQ_OBJ complete problem objective function
Offer (i ) offer at origin
Demand ( j) demand at destination
FlowLimit(i,j) arc flow limit
Bd_Cuts (l) Benders' cuts ;

EQ_Z1 .. Z1 =e= sum[(i,j), F(i,j)*Y(i,j)] + Theta ;

EQ_Z2 .. Z2 =e= sum[(i,j), C(i,j)*X(i,j)] ;

EQ_OBJ .. Z1 =e= sum[(i,j), F(i,j)*Y(i,j)] + sum[(i,j), C(i,j)*X(i,j)] ;

Offer (i ) .. sum[j, X(i,j)] =l= A(i) ;

Demand ( j) .. sum[i, X(i,j)] =g= B(j) ;

FlowLimit(i,j) .. X(i,j) =l= min[A(i),B(j)] * Y(i,j) ;

Bd_Cuts(ll) .. Delta(ll) * Theta =g= Z2_L(ll) -

sum[(i,j), PI_L(ll,i,j) * min[A(i),B(j)] * (Y_L(ll,i,j) - Y(i,j))] ;

model Master_Bd / EQ_Z1 Bd_Cuts /

model Subproblem_Bd / EQ_Z2 Offer Demand FlowLimit /
model Complete / EQ_OBJ Offer Demand FlowLimit / ;

X.up(i,j) = min[A(i),B(j)]

* to allow CPLEX correctly detect rays in an infeasible problem

* only simplex method can be used and no preprocessing neither scaling options
* optimality and feasibility tolerances are very small to avoid primal degeneration

file COPT / cplex.opt /

put COPT putclose 'ScaInd -1' / 'LPMethod 1' / 'PreInd 0' / 'EpOpt 1e-9' / 'EpRHS 1e-9' / ;

Subproblem_Bd.OptFile = 1 ;

Good Optimization Modeling Practices with GAMS. July 2016 107

FCTP solved by Benders decomposition (iii)
* parameters initialization

LL (l) = no ;
Theta.fx = 0 ;
Delta (l) = 0 ;
Z2_L (l) = 0 ;
PI_L(l,i,j) = 0 ;
Y_L (l,i,j) = 0 ;

* Benders' algorithm iterations

loop (l $(abs(1-Z_Lower/Z_Upper) > BdTol),

* solving master problem

solve Master_Bd using MIP minimizing Z1 ;

* storing the master solution

Y_L(l,i,j) = Y.l(i,j) ;

* fixing first-stage variables and solving subproblem

Y.fx (i,j) = Y.l(i,j) ;

* solving subproblem
solve Subproblem_Bd using RMIP minimizing Z2 ;

* storing parameters to build a new Benders' cut

if (Subproblem_Bd.ModelStat = 4,
Delta(l) = 0 ;
Z2_L (l) = Subproblem_Bd.SumInfes ;
else
* updating lower and upper bound
Z_Lower = Z1.l ;
Z_Upper = min(Z_Upper, Z1.l - Theta.l + Z2.l) ;

Delta(l) = 1 ;
Z2_L (l) = Subproblem_Bd.ObjVal ;
) ;

PI_L(l,i,j) = FlowLimit.m(i,j) ;

Y.lo( i,j) = 0 ;
Y.up( i,j) = 1 ;

Theta.lo = -inf ;
Theta.up = inf ;

* increase the set of Benders' cuts

LL(l) = yes ;
) ;

solve Complete using MIP minimizing Z1 Good Optimization Modeling Practices with GAMS. July 2016 108
 Stochastic FCTP
Flows
(second stage) Complete problem
Investment decisions
(first stage)
L
min
N
C ) +
,B
L
Capacity of each
L
origin ≤ ∀#M
Demand of each
L
≤$ ∀%M
destination

Flow can pass only

L
for installed ≤( ) ∀#%M
L
connections
≥ 0, ) ∈ 0,1

Good Optimization Modeling Practices with GAMS. July 2016 109

Det & Stoch FCTP $title Stochastic fixed-charge transportation problem (SFCTP)

option OptcR = 0

sets
$title Deterministic fixed-charge transportation problem (DFCTP)
I origins / i1 * i4 /
* relative optimality tolerance in solving MIP problems J destinations / j1 * j3 /
option OptcR = 0 S scenarios / s1 * s3 /

sets parameters
I origins / i1 * i4 / A(i) product offer / i1 20, i2 30, i3 40, i4 20 /
J destinations / j1 * j3 / P(s) scenario probability / s1 0.5, s2 0.3, s3 0.2 /
table B(s,j) product demand
parameters j1 j2 j3
A(i) product offer
/ i1 20, i2 30, i3 40, i4 20 / s1 21 51 31
B(j) product demand s2 32 22 52
/ j1 20, j2 50, j3 30 / s3 53 33 23
table C(i,j) per unit variable transportation cost
table C(i,j) per unit variable transportation cost j1 j2 j3
j1 j2 j3 i1 1 2 3
i1 1 2 3 i2 3 2 1
i2 3 2 1 i3 2 3 4
i3 2 3 4 i4 4 3 2
i4 4 3 2
table F(i,j) fixed transportation cost
table F(i,j) fixed transportation cost j1 j2 j3
j1 j2 j3 i1 10 20 30
i1 10 20 30 i2 20 30 40
i2 20 30 40 i3 30 40 50
i3 30 40 50 i4 40 50 60
i4 40 50 60
loop (s, abort $(sum[i, A(i)] < sum[j, B(s,j)]) 'Infeasible problem' )
abort $(sum[i, A(i)] < sum[j, B(j)]) 'Infeasible problem'

positive variable positive variable

X(i,j) arc flow X(s,i,j) arc flow
binary variable
binary variable Y(i,j ) arc investment decision
Y(i,j) arc investment decision variables
Z1 objective function
variables
Z1 objective function equations
EQ_OBJ complete problem objective function
equations
EQ_OBJ complete problem objective function Offer (s,i ) offer at origin
Offer (i ) offer at origin Demand (s, j) demand at destination
Demand ( j) demand at destination FlowLimit(s,i,j) arc flow limit ;
FlowLimit(i,j) arc flow limit ;
EQ_OBJ .. Z1 =e= sum[(i,j), F(i,j)*Y(i,j)] + sum[(s,i,j), P(s)*C(i,j)*X(s,i,j)] ;
EQ_OBJ .. Z1 =e= sum[(i,j), F(i,j)*Y(i,j)] + sum[(i,j), C(i,j)*X(i,j)] ; Offer (s,i ) .. sum[j, X(s,i,j)] =l= A( i) ;
Demand (s, j) .. sum[i, X(s,i,j)] =g= B(s,j) ;
Offer (i ) .. sum[j, X(i,j)] =l= A(i) ; FlowLimit(s,i,j) .. X(s,i,j) =l= min[A(i),B(s,j)] * Y(i,j) ;
Demand ( j) .. sum[i, X(i,j)] =g= B(j) ;
model Complete / EQ_OBJ Offer Demand FlowLimit / ;
FlowLimit(i,j) .. X(i,j) =l= min[A(i),B(j)] * Y(i,j) ;
X.up(s,i,j) = min[A(i),B(s,j)]
model Complete / EQ_OBJ Offer Demand FlowLimit / ;
solve Complete using MIP minimizing Z1
X.up(i,j) = min[A(i),B(j)]

solve Complete using MIP minimizing Z1

display X.l, Y.l
Good Optimization Modeling Practices with GAMS. July 2016 110
 Stochastic FCTP with EMP
$title Deterministic fixed-charge transportation problem (FCTP)

* relative optimality tolerance in solving MIP problems

option OptcR = 0
set S scenarios / s1 * s3 /
sets parameter P(s) scenario probability / s1 0.5, s2 0.3, s3 0.2 /
I origins / i1 * i4 / YS(s,i,j) arc investment decision
J destinations / j1 * j3 / XS(s,i,j) arc flow
parameters
A(i) product offer table BS(s,j) product demand
/ i1 20, i2 30, i3 40, i4 20 / j1 j2 j3
B(j) product demand s1 21 51 31
/ j1 11, j2 44, j3 66 / s2 32 22 52
s3 53 33 23
table C(i,j) per unit variable transportation cost
j1 j2 j3
i1 1 2 3 file emp / '%emp.info%' / ; emp.pc=2
i2 3 2 1 put emp
i3 2 3 4 put '* problem %gams.i%' / "jrandvar "
i4 4 3 2
loop (j,
table F(i,j) fixed transportation cost put B.tn(j) " "
j1 j2 j3 )
i1 10 20 30 loop (s,
i2 20 30 40 put P(s)
i3 30 40 50
loop (j,
i4 40 50 60
put BS(s,j)
positive variable )
X(i,j) arc flow )
put / "stage 2 B X Offer Demand FlowLimit"
binary variable
putclose emp
Y(i,j) arc investment decision

variables set dict / s . scenario . ''

Z1 objective function B . randvar . BS
X . level . XS
equations
Y . level . YS /
EQ_OBJ complete problem objective function
Offer (i ) offer at origin
Demand ( j) demand at destination loop (s, abort $(sum[i, A(i)] < sum[j, BS(s,j)]) 'Infeasible problem' ) ;
FlowLimit(i,j) arc flow limit ;
solve Complete minimizing Z1 using emp scenario dict
EQ_OBJ .. Z1 =e= sum[(i,j), F(i,j)*Y(i,j)] + sum[(i,j), C(i,j)*X(i,j)] ;

Offer (i ) .. sum[j, X(i,j)] =l= A(i) ; display XS, YS

Demand ( j) .. sum[i, X(i,j)] =g= B(j) ;

FlowLimit(i,j) .. X(i,j) =l= 100 * Y(i,j) ;

model Complete / all / ;

X.up(i,j) = 100
Good Optimization Modeling Practices with GAMS. July 2016 111
 Cournot model: formulation
e Company
• Objective function: profits e* Other companies
Be = p ⋅ qe - C e Be Profits
p Price
• Inverse demand function qe Output
 
p = f ∑ qe  Ce Variable costs
 e  MCe Marginal cost
• Cournot conjecture: vertical supply MRe Marginal revenue
p-MCe Price mark-up
System of ∂p ∂p ∂(qe + qe * ) ∂q
equations = = p ′ ↔ e* = 0
∂qe ∂q ∂qe ∂qe Own output decision will not
have an effect on the
decisions of the competitors
• Optimality conditions
∂ Be p − MC e (qe )
= 0 → MRe = p + qe p ′ = MC e (qe ) → qe =
∂qe −p ′

Good Optimization Modeling Practices with GAMS. July 2016 113

Cournot model: example (I)
• Perfect competition
• Company 1: MC1 = 2 €/MW q1 MAX = 5 MW
• Company 2: MC2 = 3 €/MW q2 MAX = 5 MW
• Inverse demand function: p = 10 - (q1 + q2)
€/MW p
Inverse demand function
10

Equilibrium under
perfect competition

MC2 = 3; q2 MAX = 5
3
Competitive
MC1 = 2 ; q1 MAX = 5 supply function
2
D = q1 + q2
5 7 10 MW
Good Optimization Modeling Practices with GAMS. July 2016 114
Cournot model: example (II)
• Duopoly ∂ B1
= 0 → p + q1 ⋅ (−1) = 2
• Company 1: MC1 = 2 €/MW Cournot
∂q 1
• Company 2: MC2 = 3 €/MW ∂ B2
= 0 → p + q2 ⋅ (−1) = 3
∂q 2
• IDF: p = 10 - (q1 + q2); p’=-1
€/MW p Solving

10

q1 = 3, q2 = 2
p’=-1
p = 10 - (q1 + q2) = 5
5

3
q2 = 2
2 q1 = 3
D = q1 + q2
5 7 10 MW
Good Optimization Modeling Practices with GAMS. July 2016 115
Cournot Market Equilibrium
$Title Cournot Market Equilibrium Model

$ontext

Developed by

Andrés Ramos
Instituto de Investigación Tecnológica
Escuela Técnica Superior de Ingeniería - ICAI
UNIVERSIDAD PONTIFICIA COMILLAS
Alberto Aguilera 23
28015 Madrid, Spain
Andres.Ramos@comillas.edu

August 5, 2016

$offtext

sets
p periods / p1 /
t thermal units / t1 * t2 /

parameters
pMarginalCost(t) marginal costs [€ per MW] / t1 2, t2 3 /
pMaxThrOutput(t) maximum output [MW] / t1 5, t2 5 /

pPriceInter (p) price intercept for each period [€ per MW ] / p1 10 /

pPriceSlope (p) price slope for each period [€ per MW^2] / p1 -1 /

positive variables
vSMP ( p) system marginal price [€ per MW]
vThrOutput(t,p) thermal unit output [MW]

binary variables
vCommitment(t,p) thermal unit commitment [p.u.]

Good Optimization Modeling Practices with GAMS. July 2016 116

Cournot Market Equilibrium
equations
eDemandFunction( p) price as a function of the load [€ per MW]
eDerivRevenue (t,p) derivative of the net revenue [€ per MW] ;

eDemandFunction( p) .. vSMP(p) =g= pPriceInter(p) + pPriceSlope(p)*sum[t, vThrOutput(t,p)] ;

eDerivRevenue (t,p) .. pMarginalCost(t) =g= vSMP(p) + pPriceSlope(p)* vThrOutput(t,p) ;

model mCournot / eDemandFunction.vSMP eDerivRevenue.vThrOutput /

variable
vTotalCost total investment and operation cost

equations
eTotalCost total investment and operation cost [€]
eLoadBalance(p) balance between generation and load [€ per MW] ;

eTotalCost .. vTotalCost =e= - sum[(t,p), pMarginalCost(t)*vThrOutput(t,p)] ;

eLoadBalance(p) .. pMarginalCost('t2') =e= pPriceInter(p) + pPriceSlope(p)*sum[t, vThrOutput(t,p)] ;

model mPerfectCompetition / eTotalCost eLoadBalance /

variable
vSocialWelfare social welfare [€]

equations
eSocialWelfare social welfare [€] ;

eSocialWelfare .. vSocialWelfare =e= sum[p, (pPriceInter(p) + pPriceSlope(p)*sum[t, vThrOutput(t,p)]/2)*sum[t, vThrOutput(t,p)]]

- sum[(t,p), pMarginalCost(t)*vThrOutput(t,p)] ;

model mSocialWelfare / eSocialWelfare /

vThrOutput.up(t,p) = pMaxThrOutput(t) ;

* Cournot equilibrium model

solve mCournot using MCP ;

* perfect competition model (cost minimization)

solve mPerfectCompetition using MIP maximizing vTotalCost

* perfect competition model (social welfare maximization)

solve mSocialWelfare using QCP maximizing vSocialWelfare

Good Optimization Modeling Practices with GAMS. July 2016 117

UC model competition
• Formulate mathematically the minimum up and down
time of a thermal unit
• Introduce the constraints in the UC model
• Compare the efficiency of the formulation for a large-
scale model
• Tuning the solver parameters
• Stronger formulation of the constraints

Good Optimization Modeling Practices with GAMS. July 2016 118

Antonio Machado. Cantares

“Todo pasa y todo queda, “Everything passes and everything

pero lo nuestro es pasar, stays, but our fate is to pass, to pass
pasar haciendo caminos, making paths, paths on the sea.”
caminos sobre el mar.”

“All things pass and stay forever, yet we

pass eternally, drawing footpaths in our
passing, footpaths on the restless sea.”
Good Optimization Modeling Practices with GAMS. July 2016 120
Formulating, writing and solving
optimization models

Thank you for your attention

Prof. Andres Ramos

http://www.iit.comillas.edu/aramos/
Andres.Ramos@comillas.edu
Good Optimization Modeling Practices with GAMS. July 2016 121