(!*******************************************************
Mosel Example Problems
======================
file els.mos
````````````
TYPE: Lot sizing problem
DIFFICULTY: 4
FEATURES: MIP problem, implementation of Branch-and-Cut and
Cut-and-Branch algorithms, definition of Optimizer callbacks,
Optimizer and Mosel parameter settings, `case', `procedure',
`function', time measurement
DESCRIPTION: Economic lot sizing (ELS) considers production planning
over a given planning horizon. In every period, there is
a given demand for every product that must be satisfied by
the production in this period and by inventory carried over
from previous periods.
A set-up cost is associated with production in a period,
and the total production capacity per period is limited.
The unit production cost per product and time period is
given. There is no inventory or stock-holding cost.
The model implements a configurable cutting plane algorithm
for this problem.
FURTHER INFO: `Applications of optimization with Xpress-MP teaching
material', Section 2.7 `Branch-and-Cut';
`Mosel User Guide', Section 11.1 `Cut generation'
-- This model cannot be run with a Community Licence
for the provided data instance --
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, 2001, rev. July 2023
*******************************************************!)
model "ELS"
uses "mmxprs","mmsystem", "mmsvg"
parameters
ALG = 6 ! Algorithm choice [default settings: 0]
CUTDEPTH = 10 ! Maximum tree depth for cut generation
EPS = 1e-6 ! Zero tolerance
end-parameters
forward function cb_node:boolean
forward procedure tree_cut_gen
declarations
TIMES = 1..15 ! Range of time
PRODUCTS = 1..4 ! Set of products
DEMAND: array(PRODUCTS,TIMES) of integer ! Demand per period
SETUPCOST: array(TIMES) of integer ! Setup cost per period
PRODCOST: array(PRODUCTS,TIMES) of integer ! Production cost per period
CAP: array(TIMES) of integer ! Production capacity per period
D: array(PRODUCTS,TIMES,TIMES) of integer ! Total demand in periods t1 - t2
produce: array(PRODUCTS,TIMES) of mpvar ! Production in period t
setup: array(PRODUCTS,TIMES) of mpvar ! Setup in period t
solprod: array(PRODUCTS,TIMES) of real ! Sol. values for var.s produce
solsetup: array(PRODUCTS,TIMES) of real ! Sol. values for var.s setup
starttime: real
end-declarations
initializations from "els.dat"
DEMAND SETUPCOST PRODCOST CAP
end-initializations
forall(p in PRODUCTS,s,t in TIMES) D(p,s,t):= sum(k in s..t) DEMAND(p,k)
! Objective: minimize total cost
MinCost:= sum(t in TIMES) (SETUPCOST(t) * sum(p in PRODUCTS) setup(p,t) +
sum(p in PRODUCTS) PRODCOST(p,t) * produce(p,t) )
! Satisfy the total demand
forall(p in PRODUCTS,t in TIMES)
Dem(p,t):= sum(s in 1..t) produce(p,s) >= sum (s in 1..t) DEMAND(p,s)
! If there is production during t then there is a setup in t
forall(p in PRODUCTS, t in TIMES)
ProdSetup(p,t):= produce(p,t) <= D(p,t,getlast(TIMES)) * setup(p,t)
! Capacity limits
forall(t in TIMES) Capacity(t):= sum(p in PRODUCTS) produce(p,t) <= CAP(t)
! Variables setup are 0/1
forall(p in PRODUCTS, t in TIMES) setup(p,t) is_binary
! Uncomment to get detailed MIP output
setparam("XPRS_VERBOSE", true)
! All cost data are integer, we therefore only need to search for integer
! solutions
setparam("XPRS_MIPADDCUTOFF", -0.999)
! Set Mosel comparison tolerance to a sufficiently small value
setparam("ZEROTOL", EPS/100)
writeln("**************ALG=",ALG,"***************")
starttime:=gettime
SEVERALROUNDS:=false; TOPONLY:=false
case ALG of
1: setparam("XPRS_CUTSTRATEGY", 0) ! No cuts
2: setparam("XPRS_PRESOLVE", 0) ! No presolve
3: tree_cut_gen ! User branch-and-cut + automatic cuts
4: do ! User branch-and-cut (several rounds),
tree_cut_gen ! no automatic cuts
setparam("XPRS_CUTSTRATEGY", 0)
SEVERALROUNDS:=true
end-do
5: do ! User cut-and-branch (several rounds)
tree_cut_gen ! + automatic cuts
SEVERALROUNDS:=true
TOPONLY:=true
end-do
6: do ! User branch-and-cut (several rounds)
tree_cut_gen ! + automatic cuts
SEVERALROUNDS:=true
end-do
end-case
minimize(MinCost) ! Solve the problem
writeln("Time: ", gettime-starttime, "sec, Nodes: ", getparam("XPRS_NODES"),
", Solution: ", getobjval)
write("Period setup ")
forall(p in PRODUCTS) write(strfmt(p,-7))
forall(t in TIMES) do
write("\n ", strfmt(t,2), strfmt(getsol(sum(p in PRODUCTS) setup(p,t)),8), " ")
forall(p in PRODUCTS) write(getsol(produce(p,t)), " (",DEMAND(p,t),") ")
end-do
writeln
! Draw the solution
forall(p in PRODUCTS) do
svgaddgroup("P"+p, "Product "+p)
svgsetstyle(SVG_FILL, SVG_CURRENT)
svgsetstyle(SVG_STROKEWIDTH, 0)
end-do
forall(t in TIMES) do
cum:=0.0
forall(p in PRODUCTS | produce(p,t).sol>0) do
svgaddrectangle("P"+p,t, cum, 0.9, setup(p,t).sol)
svgsetstyle(svggetlastobj, SVG_OPACITY, 0.5)
cum+=setup(p,t).sol
svgaddrectangle("P"+p,t, cum, 0.9, produce(p,t).sol)
cum+=produce(p,t).sol
end-do
end-do
! svgsetgraphviewbox(0, 0, TIMES.last+1, max(t in TIMES) CAP(t) +1)
svgsetgraphlabels("Time", "Production amounts and setup")
svgsave("els.svg")
svgrefresh
svgwaitclose("Close browser window to terminate model execution.", 1)
!*************************************************************************
! Cut generation loop:
! get the solution values
! identify and set up violated constraints
! load the modified problem and load the saved basis
!*************************************************************************
function cb_node:boolean
declarations
ncut:integer ! Counter for cuts
cut: array(range) of linctr ! Cuts
cutid: array(range) of integer ! Cut type identification
type: array(range) of integer ! Cut constraint type
ds: real
end-declarations
returned:=false ! OPTNODE: This node is not infeasible
depth:=getparam("XPRS_NODEDEPTH")
cnt:=getparam("XPRS_CALLBACKCOUNT_OPTNODE")
if ((TOPONLY and depth<1) or (not TOPONLY and depth<=CUTDEPTH)) and
(SEVERALROUNDS or cnt<=1) then
ncut:=0
! Get the solution values
forall(t in TIMES, p in PRODUCTS) do
solprod(p,t):=getsol(produce(p,t))
solsetup(p,t):=getsol(setup(p,t))
end-do
! Search for violated constraints
forall(p in PRODUCTS,l in TIMES) do
ds:=0
forall(t in 1..l)
if (solprod(p,t) < D(p,t,l)*solsetup(p,t) + EPS) then ds += solprod(p,t)
else ds += D(p,t,l)*solsetup(p,t)
end-if
! Generate the violated inequality
if ds < D(p,1,l) - EPS then
cut(ncut):= sum(t in 1..l)
if(solprod(p,t)<(D(p,t,l)*solsetup(p,t))+EPS, produce(p,t),
D(p,t,l)*setup(p,t)) - D(p,1,l)
cutid(ncut):= 1
type(ncut):= CT_GEQ
ncut+=1
end-if
end-do
! Add cuts to the problem
if ncut>0 then
addcuts(cutid, type, cut);
if getparam("XPRS_VERBOSE")=true then
writeln("Cuts added : ", ncut, " (depth ", depth, ", node ",
getparam("XPRS_NODES"), ", obj. ", getparam("XPRS_LPOBJVAL"), ")")
end-if
end-if
end-if
end-function
! ****Optimizer settings for using the cut manager****
procedure tree_cut_gen
setparam("XPRS_PRESOLVE", 0) ! Switch presolve off
setparam("XPRS_ROOTPRESOLVE", 0) ! Switch B&B root presolve off
setparam("XPRS_EXTRAROWS", 5000) ! Reserve extra rows in matrix
setcallback(XPRS_CB_OPTNODE, ->cb_node) ! Set the optnode callback function
end-procedure
end-model
|
