paperp.mos |
(!*******************************************************
Mosel User Guide Examples
=========================
file paperp.mos
``````````````
Column generation algorithm for finding (near-)optimal
solutions for a paper cutting example.
Defining several (parallel) models.
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, 2004, rev. May 2019
*******************************************************!)
model "Papermill (multi-prob)"
uses "mmxprs", "mmjobs", "mmsystem"
forward procedure column_gen
forward function knapsack(C:array(range) of real,
A:array(range) of real,
B:real,
D:array(range) of integer,
xbest:array(range) of integer): real
forward procedure show_new_pat(dj:real, vx: array(range) of integer)
(!
The column generation algorithm requires the solution of a knapsack problem.
This problem is defined in a separate model file.
!)
declarations
NWIDTHS = 5 ! Number of different widths
WIDTHS = 1..NWIDTHS ! Range of widths
RP: range ! Range of cutting patterns
MAXWIDTH = 94 ! Maximum roll width
EPS = 1e-6 ! Zero tolerance
WIDTH: array(WIDTHS) of real ! Possible widths
DEMAND: array(WIDTHS) of integer ! Demand per width
PATTERNS: array(WIDTHS, WIDTHS) of integer ! (Basic) cutting patterns
use: dynamic array(RP) of mpvar ! Rolls per pattern
soluse: dynamic array(RP) of real ! Solution values for variables `use'
Dem: array(WIDTHS) of linctr ! Demand constraints
MinRolls: linctr ! Objective function
Knapsack: Model ! Reference to the knapsack model
end-declarations
WIDTH:: [ 17, 21, 22.5, 24, 29.5]
DEMAND:: [150, 96, 48, 108, 227]
! Make basic patterns
forall(j in WIDTHS) PATTERNS(j,j) := floor(MAXWIDTH/WIDTH(j))
(! Alternative starting patterns, particularly for small demand values:
forall(j in WIDTHS)
PATTERNS(j,j) := minlist(floor(MAXWIDTH/WIDTH(j)),DEMAND(j))
forall(j in WIDTHS) PATTERNS(j,j) := 1
!)
forall(j in WIDTHS) do
create(use(j)) ! Create NWIDTHS variables `use'
use(j) is_integer ! Variables are integer and bounded
use(j) <= integer(ceil(DEMAND(j)/PATTERNS(j,j)))
end-do
MinRolls:= sum(j in WIDTHS) use(j) ! Objective: minimize no. of rolls
! Satisfy all demands
forall(i in WIDTHS)
Dem(i):= sum(j in WIDTHS) PATTERNS(i,j) * use(j) >= DEMAND(i)
res:= compile("knapsack.mos") ! Compile the knapsack model
load(Knapsack, "knapsack.bim") ! Load the knapsack model
column_gen ! Column generation at top node
minimize(MinRolls) ! Compute the best integer solution
! for the current problem (including
! the new columns)
writeln("Best integer solution: ", getobjval, " rolls")
write(" Rolls per pattern: ")
forall(i in RP) write(getsol(use(i)),", ")
writeln
fdelete("knapsack.bim") ! Cleaning up
fdelete("shmem:indata"); fdelete("shmem:resdata")
!************************************************************************
! Column generation loop at the top node:
! solve the LP and save the basis
! get the solution values
! generate new column(s) (=cutting pattern)
! load the modified problem and load the saved basis
!************************************************************************
procedure column_gen
declarations
dualdem: array(WIDTHS) of real
xbest: array(WIDTHS) of integer
dw, zbest, objval: real
bas: basis
end-declarations
setparam("XPRS_PRESOLVE", 0) ! Switch presolve off
setparam("zerotol", EPS) ! Set comparison tolerance of Mosel
npatt:=NWIDTHS
npass:=1
while(true) do
minimize(XPRS_LIN, MinRolls) ! Solve the LP
savebasis(bas) ! Save the current basis
objval:= getobjval ! Get the objective value
! Get the solution values
forall(j in 1..npatt) soluse(j):=getsol(use(j))
forall(i in WIDTHS) dualdem(i):=getdual(Dem(i))
! Solve a knapsack problem
zbest:= knapsack(dualdem, WIDTH, MAXWIDTH, DEMAND, xbest) - 1.0
write("Pass ", npass, ": ")
if zbest = 0 then ! Same as zbest <= EPS
writeln("no profitable column found.\n")
break
else
show_new_pat(zbest, xbest) ! Print the new pattern
npatt+=1
create(use(npatt)) ! Create a new var. for this pattern
use(npatt) is_integer
MinRolls+= use(npatt) ! Add new var. to the objective
dw:=0
forall(i in WIDTHS)
if xbest(i) > 0 then
Dem(i)+= xbest(i)*use(npatt) ! Add new var. to demand constr.s
dw:= maxlist(dw, ceil(DEMAND(i)/xbest(i) ))
end-if
use(npatt) <= dw ! Set upper bound on the new var.
! (useful when PRESOLVE is off)
loadprob(MinRolls) ! Reload the problem
loadbasis(bas) ! Load the saved basis
end-if
npass+=1
end-do
writeln("Solution after column generation: ", objval, " rolls, ",
getsize(RP), " patterns")
write(" Rolls per pattern: ")
forall(i in RP) write(soluse(i),", ")
writeln
setparam("XPRS_PRESOLVE", 1) ! Switch presolve on
end-procedure
!************************************************************************
! Solve the integer knapsack problem
! z = max{cx : ax<=b, x<=d, x in Z^N}
! with b=MAXWIDTH, N=NWIDTHS, d=DEMAND
!
! All data is exchanged between the two models via shared memory.
!************************************************************************
function knapsack(C:array(range) of real,
A:array(range) of real,
B:real,
D:array(range) of integer,
xbest:array(range) of integer):real
INDATA := "bin:shmem:indata"
RESDATA := "bin:shmem:resdata"
initializations to INDATA
A B C D
end-initializations
run(Knapsack, "NWIDTHS=" + NWIDTHS + ",INDATA=" + INDATA +
",RESDATA=" + RESDATA) ! Start knapsack (sub)model
wait ! Wait until subproblem finishes
dropnextevent ! Ignore termination message
initializations from RESDATA
xbest returned as "zbest"
end-initializations
end-function
!************************************************************************
procedure show_new_pat(dj:real, vx: array(range) of integer)
declarations
dw: real
end-declarations
writeln("new pattern found with marginal cost ", dj)
write(" Widths distribution: ")
dw:=0
forall(i in WIDTHS) do
write(WIDTH(i), ":", vx(i), " ")
dw += WIDTH(i)*vx(i)
end-do
writeln("Total width: ", dw)
end-procedure
end-model
|
|
knapsack.mos |
(!*******************************************************
Mosel User Guide Examples
=========================
file knapsack.mos
`````````````````
Knapsack (sub)model of the for a paper cutting example,
reading data from shared memory.
*** Can be run standalone - or run from paperp.mos ***
Solve the integer knapsack problem
z = max{Cx : Ax<=B, x<=D, x in Z^N}
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, 2004, rev. Aug. 2011
*******************************************************!)
model "Knapsack"
uses "mmxprs"
parameters
NWIDTHS=5 ! Number of different widths
INDATA = "knapsack.dat" ! Input data file
RESDATA = "knresult.dat" ! Result data file
end-parameters
declarations
WIDTHS = 1..NWIDTHS ! Range of widths
A,C: array(WIDTHS) of real ! Constraint + obj. coefficients
B: real ! RHS value of knapsack constraint
D: array(WIDTHS) of integer ! Variables bounds (demand quantities)
KnapCtr, KnapObj: linctr ! Knapsack constraint+objective
x: array(WIDTHS) of mpvar ! Knapsack variables
xbest: array(WIDTHS) of integer ! Solution values
end-declarations
initializations from INDATA
A B C D
end-initializations
! Define the knapsack problem
KnapCtr:= sum(j in WIDTHS) A(j)*x(j) <= B
KnapObj:= sum(j in WIDTHS) C(j)*x(j)
! Integrality condition and bounds
forall(j in WIDTHS) x(j) is_integer
forall(j in WIDTHS) x(j) <= D(j) ! These bounds can be omitted
maximize(KnapObj)
z:=getobjval
forall(j in WIDTHS) xbest(j):=round(getsol(x(j)))
initializations to RESDATA
xbest z as "zbest"
end-initializations
end-model
|
|
paperpr.mos |
(!*******************************************************
Mosel User Guide Examples
=========================
file paperpr.mos
````````````````
Column generation algorithm for finding (near-)optimal
solutions for a paper cutting example.
Defining several (parallel) models.
Communication of data to/from submodels
using 'shmem' IO driver with 'raw'.
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, 2004, rev. May 2019
*******************************************************!)
model "Papermill (multi-prob)"
uses "mmxprs", "mmjobs", "mmsystem"
forward procedure column_gen
forward function knapsack(C:array(range) of real,
A:array(range) of real,
B:real,
D:array(range) of integer,
xbest:array(range) of integer): real
forward procedure show_new_pat(dj:real, vx: array(range) of integer)
(!
The column generation algorithm requires the solution of a knapsack problem.
This problem is defined in a separate model file.
!)
declarations
NWIDTHS = 5 ! Number of different widths
WIDTHS = 1..NWIDTHS ! Range of widths
RP: range ! Range of cutting patterns
MAXWIDTH = 94 ! Maximum roll width
EPS = 1e-6 ! Zero tolerance
WIDTH: array(WIDTHS) of real ! Possible widths
DEMAND: array(WIDTHS) of integer ! Demand per width
PATTERNS: array(WIDTHS, WIDTHS) of integer ! (Basic) cutting patterns
use: dynamic array(RP) of mpvar ! Rolls per pattern
soluse: dynamic array(RP) of real ! Solution values for variables `use'
Dem: array(WIDTHS) of linctr ! Demand constraints
MinRolls: linctr ! Objective function
Knapsack: Model ! Reference to the knapsack model
end-declarations
WIDTH:: [ 17, 21, 22.5, 24, 29.5]
DEMAND:: [150, 96, 48, 108, 227]
! Make basic patterns
forall(j in WIDTHS) PATTERNS(j,j) := floor(MAXWIDTH/WIDTH(j))
(! Alternative starting patterns, particularly for small demand values:
forall(j in WIDTHS)
PATTERNS(j,j) := minlist(floor(MAXWIDTH/WIDTH(j)),DEMAND(j))
forall(j in WIDTHS) PATTERNS(j,j) := 1
!)
forall(j in WIDTHS) do
create(use(j)) ! Create NWIDTHS variables `use'
use(j) is_integer ! Variables are integer and bounded
use(j) <= integer(ceil(DEMAND(j)/PATTERNS(j,j)))
end-do
MinRolls:= sum(j in WIDTHS) use(j) ! Objective: minimize no. of rolls
! Satisfy all demands
forall(i in WIDTHS)
Dem(i):= sum(j in WIDTHS) PATTERNS(i,j) * use(j) >= DEMAND(i)
setparam("XPRS_verbose",true)
res:= compile("knapsackr.mos") ! Compile the knapsack model
load(Knapsack, "knapsackr.bim") ! Load the knapsack model
column_gen ! Column generation at top node
minimize(MinRolls) ! Compute the best integer solution
! for the current problem (including
! the new columns)
writeln("Best integer solution: ", getobjval, " rolls")
write(" Rolls per pattern: ")
forall(i in RP) write(getsol(use(i)),", ")
writeln
fdelete("knapsackr.bim") ! Cleaning up
!************************************************************************
! Column generation loop at the top node:
! solve the LP and save the basis
! get the solution values
! generate new column(s) (=cutting pattern)
! load the modified problem and load the saved basis
!************************************************************************
procedure column_gen
declarations
dualdem: array(WIDTHS) of real
xbest: array(WIDTHS) of integer
dw, zbest, objval: real
bas: basis
end-declarations
setparam("XPRS_PRESOLVE", 0) ! Switch presolve off
setparam("zerotol", EPS) ! Set comparison tolerance of Mosel
npatt:=NWIDTHS
npass:=1
while(true) do
minimize(XPRS_LIN, MinRolls) ! Solve the LP
savebasis(bas) ! Save the current basis
objval:= getobjval ! Get the objective value
! Get the solution values
forall(j in 1..npatt) soluse(j):=getsol(use(j))
forall(i in WIDTHS) dualdem(i):=getdual(Dem(i))
! Solve a knapsack problem
zbest:= knapsack(dualdem, WIDTH, MAXWIDTH, DEMAND, xbest) - 1.0
write("Pass ", npass, ": ")
if zbest = 0 then ! Same as zbest <= EPS
writeln("no profitable column found.\n")
break
else
show_new_pat(zbest, xbest) ! Print the new pattern
npatt+=1
create(use(npatt)) ! Create a new var. for this pattern
use(npatt) is_integer
MinRolls+= use(npatt) ! Add new var. to the objective
dw:=0
forall(i in WIDTHS)
if xbest(i) > 0 then
Dem(i)+= xbest(i)*use(npatt) ! Add new var. to demand constr.s
dw:= maxlist(dw, ceil(DEMAND(i)/xbest(i) ))
end-if
use(npatt) <= dw ! Set upper bound on the new var.
! (useful when PRESOLVE is off)
loadprob(MinRolls) ! Reload the problem
loadbasis(bas) ! Load the saved basis
end-if
npass+=1
end-do
writeln("Solution after column generation: ", objval, " rolls, ",
getsize(RP), " patterns")
write(" Rolls per pattern: ")
forall(i in RP) write(soluse(i),", ")
writeln
setparam("XPRS_PRESOLVE", 1) ! Switch presolve on
end-procedure
!************************************************************************
! Solve the integer knapsack problem
! z = max{cx : ax<=b, x<=d, x in Z^N}
! with b=MAXWIDTH, N=NWIDTHS, d=DEMAND
!
! All data is exchanged between the two models via shared memory.
!************************************************************************
function knapsack(C:array(range) of real,
A:array(range) of real,
B:real,
D:array(range) of integer,
xbest:array(range) of integer):real
initializations to "raw:noindex"
A as "shmem:A" B as "shmem:B" C as "shmem:C" D as "shmem:D"
end-initializations
run(Knapsack, "NWIDTHS="+NWIDTHS) ! Start solving knapsack subproblem
wait ! Wait until subproblem finishes
dropnextevent ! Ignore termination message
initializations from "raw:"
xbest as "shmem:xbest" returned as "shmem:zbest"
end-initializations
end-function
!************************************************************************
procedure show_new_pat(dj:real, vx: array(range) of integer)
declarations
dw: real
end-declarations
writeln("new pattern found with marginal cost ", dj)
write(" Widths distribution: ")
dw:=0
forall(i in WIDTHS) do
write(WIDTH(i), ":", vx(i), " ")
dw += WIDTH(i)*vx(i)
end-do
writeln("Total width: ", dw)
end-procedure
end-model
|
|
knapsackr.mos |
(!*******************************************************
Mosel User Guide Examples
=========================
file knapsackr.mos
``````````````````
Knapsack (sub)model of the for a paper cutting example,
reading data from shared memory.
Communication of data from/to parent model
using 'shmem' IO driver with 'raw'.
*** Not intended to be run standalone - run from paperpr.mos ***
Solve the integer knapsack problem
z = max{Cx : Ax<=B, x<=D, x in Z^N}
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, 2004, rev. July 2010
*******************************************************!)
model "Knapsack"
uses "mmxprs"
parameters
NWIDTHS=5 ! Number of different widths
end-parameters
declarations
WIDTHS = 1..NWIDTHS ! Range of widths
A,C: array(WIDTHS) of real ! Constraint + obj. coefficients
B: real ! RHS value of knapsack constraint
D: array(WIDTHS) of integer ! Variables bounds (demand quantities)
KnapCtr, KnapObj: linctr ! Knapsack constraint+objective
x: array(WIDTHS) of mpvar ! Knapsack variables
xbest: array(WIDTHS) of integer ! Solution values
end-declarations
initializations from "raw:noindex"
A as "shmem:A" B as "shmem:B" C as "shmem:C" D as "shmem:D"
end-initializations
! Define the knapsack problem
KnapCtr:= sum(j in WIDTHS) A(j)*x(j) <= B
KnapObj:= sum(j in WIDTHS) C(j)*x(j)
! Integrality condition and bounds
forall(j in WIDTHS) x(j) is_integer
forall(j in WIDTHS) x(j) <= D(j) ! These bounds can be omitted
maximize(KnapObj)
z:=getobjval
forall(j in WIDTHS) xbest(j):=round(getsol(x(j)))
initializations to "raw:"
xbest as "shmem:xbest" z as "shmem:zbest"
end-initializations
end-model
|
|
papersn.mos |
(!*******************************************************
Mosel User Guide Examples
=========================
file papersn.mos
````````````````
Column generation algorithm for finding (near-)optimal
solutions for a paper cutting example.
- Working with multiple problems within a single model,
recreating a new problem -
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, Nov. 2008, rev. May 2019
*******************************************************!)
model Papermill
uses "mmxprs"
forward procedure column_gen
forward function knapsack(C:array(range) of real,
A:array(range) of real,
B:real,
D:array(range) of integer,
xbest:array(range) of integer): real
forward procedure show_new_pat(dj:real, vx: array(range) of integer)
(!
The column generation algorithm requires the solution of a knapsack problem.
We define the constraints of the knapsack problem globally because Mosel
does not delete them at the end of a function (due to its principle of incremental model formulation).
!)
declarations
NWIDTHS = 5 ! Number of different widths
WIDTHS = 1..NWIDTHS ! Range of widths
RP: range ! Range of cutting patterns
MAXWIDTH = 94 ! Maximum roll width
EPS = 1e-6 ! Zero tolerance
WIDTH: array(WIDTHS) of real ! Possible widths
DEMAND: array(WIDTHS) of integer ! Demand per width
PATTERNS: array(WIDTHS, WIDTHS) of integer ! (Basic) cutting patterns
use: dynamic array(RP) of mpvar ! Rolls per pattern
soluse: dynamic array(RP) of real ! Solution values for variables `use'
Dem: array(WIDTHS) of linctr ! Demand constraints
MinRolls: linctr ! Objective function
end-declarations
WIDTH:: [ 17, 21, 22.5, 24, 29.5]
DEMAND:: [150, 96, 48, 108, 227]
! Make basic patterns
forall(j in WIDTHS) PATTERNS(j,j) := floor(MAXWIDTH/WIDTH(j))
(! Alternative starting patterns, particularly for small demand values:
forall(j in WIDTHS)
PATTERNS(j,j) := minlist(floor(MAXWIDTH/WIDTH(j)),DEMAND(j))
forall(j in WIDTHS) PATTERNS(j,j) := 1
!)
forall(j in WIDTHS) do
create(use(j)) ! Create NWIDTHS variables `use'
use(j) is_integer ! Variables are integer and bounded
use(j) <= integer(ceil(DEMAND(j)/PATTERNS(j,j)))
end-do
MinRolls:= sum(j in WIDTHS) use(j) ! Objective: minimize no. of rolls
! Satisfy all demands
forall(i in WIDTHS)
Dem(i):= sum(j in WIDTHS) PATTERNS(i,j) * use(j) >= DEMAND(i)
column_gen ! Column generation at top node
minimize(MinRolls) ! Compute the best integer solution
! for the current problem (including
! the new columns)
writeln("Best integer solution: ", getobjval, " rolls")
write(" Rolls per pattern: ")
forall(i in RP) write(getsol(use(i)),", ")
writeln
!************************************************************************
! Column generation loop at the top node:
! solve the LP and save the basis
! get the solution values
! generate new column(s) (=cutting pattern)
! load the modified problem and load the saved basis
!************************************************************************
procedure column_gen
declarations
dualdem: array(WIDTHS) of real
xbest: array(WIDTHS) of integer
dw, zbest, objval: real
bas: basis
end-declarations
setparam("XPRS_PRESOLVE", 0) ! Switch presolve off
setparam("zerotol", EPS) ! Set comparison tolerance of Mosel
npatt:=NWIDTHS
npass:=1
while(true) do
minimize(XPRS_LIN, MinRolls) ! Solve the LP
savebasis(bas) ! Save the current basis
objval:= getobjval ! Get the objective value
! Get the solution values
forall(j in 1..npatt) soluse(j):=getsol(use(j))
forall(i in WIDTHS) dualdem(i):=getdual(Dem(i))
! Solve a knapsack problem
zbest:= knapsack(dualdem, WIDTH, MAXWIDTH, DEMAND, xbest) - 1.0
write("Pass ", npass, ": ")
if zbest = 0 then ! Same as zbest <= EPS
writeln("no profitable column found.\n")
break
else
show_new_pat(zbest, xbest) ! Print the new pattern
npatt+=1
create(use(npatt)) ! Create a new var. for this pattern
use(npatt) is_integer
MinRolls+= use(npatt) ! Add new var. to the objective
dw:=0
forall(i in WIDTHS)
if xbest(i) > 0 then
Dem(i)+= xbest(i)*use(npatt) ! Add new var. to demand constr.s
dw:= maxlist(dw, ceil(DEMAND(i)/xbest(i) ))
end-if
use(npatt) <= dw ! Set upper bound on the new var.
! (useful when PRESOLVE is off)
loadprob(MinRolls) ! Reload the problem
loadbasis(bas) ! Load the saved basis
end-if
npass+=1
end-do
writeln("Solution after column generation: ", objval, " rolls, ",
getsize(RP), " patterns")
write(" Rolls per pattern: ")
forall(i in RP) write(soluse(i),", ")
writeln
setparam("XPRS_PRESOLVE", 1) ! Switch presolve on
end-procedure
!************************************************************************
! Solve the integer knapsack problem
! z = max{cx : ax<=b, x<=d, x in Z^N}
! with b=MAXWIDTH, N=NWIDTHS, d=DEMAND
!************************************************************************
function knapsack(C:array(range) of real,
A:array(range) of real,
B:real,
D:array(range) of integer,
xbest:array(range) of integer):real
declarations
x: array(WIDTHS) of mpvar ! Knapsack variables
end-declarations
with mpproblem do ! Create a local subproblem
! Define the knapsack problem
forall(j in WIDTHS) x(j) is_integer
forall(j in WIDTHS) x(j) <= D(j) ! These bounds can be omitted
sum(j in WIDTHS) A(j)*x(j) <= B
maximize(sum(j in WIDTHS) C(j)*x(j))
returned:=getobjval
forall(j in WIDTHS) xbest(j):=round(getsol(x(j)))
end-do
end-function
!************************************************************************
procedure show_new_pat(dj:real, vx: array(range) of integer)
declarations
dw: real
end-declarations
writeln("new pattern found with marginal cost ", dj)
write(" Widths distribution: ")
dw:=0
forall(i in WIDTHS) do
write(WIDTH(i), ":", vx(i), " ")
dw += WIDTH(i)*vx(i)
end-do
writeln("Total width: ", dw)
end-procedure
end-model
|
|