(!*******************************************************
Mosel Example Problems
======================
file benders_main.mos
`````````````````````
Benders decomposition for solving a simple MIP.
- Main model -
*** ATTENTION: This model will return an error if ***
*** no more than one Xpress licence is available. ***
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, Jan. 2006, rev. Sep. 2022
*******************************************************!)
model "Benders (main model)"
uses "mmxprs", "mmjobs", "mmsystem"
parameters
NCTVAR = 3 ! 1
NINTVAR = 3 ! 2
NC = 4 ! 3
BIGM = 1000
ALG = 1 ! 1: Use Benders decomposition (dual)
! 2: Use Benders decomposition (primal)
DATAFILE = "bprob33.dat" ! "bprob12.dat"
end-parameters
forward procedure start_solution
forward procedure solve_primal_int(ct: integer)
forward procedure solve_cont
forward function eval_solution: boolean
forward procedure print_solution
declarations
STEP_0=2 ! Event codes sent to submodels
STEP_1=3
STEP_2=4
EVENT_SOLVED=6 ! Event codes sent by submodels
EVENT_INFEAS=7
EVENT_READY=8
CtVars = 1..NCTVAR ! Continuous variables
IntVars = 1..NINTVAR ! Discrete variables
Ctrs = 1..NC ! Set of constraints (orig. problem)
A: array(Ctrs,CtVars) of integer ! Coeff.s of continuous variables
B: array(Ctrs,IntVars) of integer ! Coeff.s of discrete variables
b: array(Ctrs) of integer ! RHS values
C: array(CtVars) of integer ! Obj. coeff.s of continuous variables
D: array(IntVars) of integer ! Obj. coeff.s of discrete variables
sol_u: array(Ctrs) of real ! Solution of dual problem
sol_x: array(CtVars) of real ! Solution of primal prob. (cont.)
sol_y: array(IntVars) of real ! Solution of primal prob. (integers)
sol_obj: real ! Objective function value (primal)
RM: range ! Model indices
stepmod: dynamic array(RM) of Model ! Submodels
end-declarations
initializations from DATAFILE
A B b C D
end-initializations
! **** Submodels ****
initializations to "bin:shmem:probdata" ! Save data for submodels
A B b C D
end-initializations
! Compile + load all submodels
if compile("benders_int.mos")<>0: exit(1)
create(stepmod(1)); load(stepmod(1), "benders_int.bim")
if compile("benders_dual.mos")<>0: exit(2)
if ALG=1 then
create(stepmod(2)); load(stepmod(2), "benders_dual.bim")
else
create(stepmod(0)); load(stepmod(0), "benders_dual.bim")
if compile("benders_cont.mos")<>0: exit(3)
create(stepmod(2)); load(stepmod(2), "benders_cont.bim")
run(stepmod(0), "NCTVAR=" + NCTVAR + ",NINTVAR=" + NINTVAR + ",NC=" + NC)
end-if
! Start the execution of the submodels
run(stepmod(1), "NINTVAR=" + NINTVAR + ",NC=" + NC)
run(stepmod(2), "NCTVAR=" + NCTVAR + ",NINTVAR=" + NINTVAR + ",NC=" + NC)
forall(m in RM) do
wait ! Wait for "Ready" messages
ev:= getnextevent
if getclass(ev) <> EVENT_READY then
writeln("Error occurred in a subproblem")
exit(4)
end-if
end-do
! **** Solution algorithm ****
start_solution ! 0. Initial solution for cont. var.s
ct:= 1
repeat
writeln("\n**** Iteration: ", ct)
solve_primal_int(ct) ! 1. Solve problem with fixed cont.
solve_cont ! 2. Solve problem with fixed int.
ct+=1
until eval_solution ! Test for optimality
print_solution ! 3. Retrieve and display the solution
! **** Cleaning up temporary files ****
fdelete("benders_int.bim")
fdelete("benders_dual.bim")
if ALG<>1: fdelete("benders_cont.bim")
fdelete("shmem:probdata")
fdelete("shmem:sol")
!-----------------------------------------------------------
! Produce an initial solution for the dual problem using a dummy objective
procedure start_solution
if ALG=1 then ! Start the problem solving
send(stepmod(2), STEP_0, 0)
else
send(stepmod(0), STEP_0, 0)
end-if
wait ! Wait for the solution
ev:=getnextevent
if getclass(ev)=EVENT_INFEAS then
writeln("Problem is infeasible")
exit(6)
end-if
end-procedure
!-----------------------------------------------------------
! Solve a modified version of the primal problem, replacing continuous
! variables by the solution of the dual
procedure solve_primal_int(ct: integer)
send(stepmod(1), STEP_1, ct) ! Start the problem solving
wait ! Wait for the solution
ev:=getnextevent
sol_obj:= getvalue(ev) ! Store objective function value
initializations from "bin:shmem:sol" ! Retrieve the solution
sol_y
end-initializations
end-procedure
!-----------------------------------------------------------
! Solve the Step 2 problem (dual or primal depending on value of ALG)
! for given solution values of y
procedure solve_cont
send(stepmod(2), STEP_2, 0) ! Start the problem solving
wait ! Wait for the solution
dropnextevent
initializations from "bin:shmem:sol" ! Retrieve the solution
sol_u
end-initializations
end-procedure
!-----------------------------------------------------------
function eval_solution: boolean
write("Test optimality: ", sol_obj - sum(i in IntVars) D(i)*sol_y(i),
" >= ", sum(j in Ctrs) sol_u(j)* (b(j) -
sum(i in IntVars) B(j,i)*sol_y(i)) )
returned:= ( sol_obj - sum(i in IntVars) D(i)*sol_y(i) >=
sum(j in Ctrs) sol_u(j)* (b(j) - sum(i in IntVars) B(j,i)*sol_y(i)) )
writeln(if(returned, " : true", " : false"))
end-function
!-----------------------------------------------------------
procedure print_solution
! Retrieve results
initializations from "bin:shmem:sol"
sol_x
end-initializations
forall(m in RM) stop(stepmod(m)) ! Stop all submodels
write("\n**** Solution (Benders): ", sol_obj, "\n x: ")
forall(i in CtVars) write(sol_x(i), " ")
write(" y: ")
forall(i in IntVars) write(sol_y(i), " ")
writeln
end-procedure
end-model
|
(!*******************************************************
Mosel Example Problems
======================
file benders_dual.mos
`````````````````````
Benders decomposition for solving a simple MIP.
- Dual problem in the continuous variables
for start solution and Step 2 of the algorithm -
*** Not intended to be run standalone - run from benders_main.mos ***
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, Jan. 2006, rev. Jan. 2013
*******************************************************!)
model "Benders (dual problem)"
uses "mmxprs", "mmjobs"
parameters
NCTVAR = 3
NINTVAR = 3
NC = 4
BIGM = 1000
end-parameters
forward procedure save_solution
declarations
STEP_0=2 ! Event codes sent to submodels
STEP_2=4
EVENT_SOLVED=6 ! Event codes sent by submodels
EVENT_INFEAS=7
EVENT_READY=8
CtVars = 1..NCTVAR ! Continuous variables
IntVars = 1..NINTVAR ! Discrete variables
Ctrs = 1..NC ! Set of constraints (orig. problem)
A: array(Ctrs,CtVars) of integer ! Coeff.s of continuous variables
B: array(Ctrs,IntVars) of integer ! Coeff.s of discrete variables
b: array(Ctrs) of integer ! RHS values
C: array(CtVars) of integer ! Obj. coeff.s of continuous variables
sol_u: array(Ctrs) of real ! Solution of dual problem
sol_x: array(CtVars) of real ! Solution of primal prob. (cont.)
sol_y: array(IntVars) of real ! Solution of primal prob. (integers)
u: array(Ctrs) of mpvar ! Dual variables
end-declarations
initializations from "bin:shmem:probdata"
A B b C
end-initializations
forall(i in CtVars) CtrD(i):= sum(j in Ctrs) u(j)*A(j,i) <= C(i)
send(EVENT_READY,0) ! Model is ready (= running)
! (Re)solve this model until it is stopped by event "STEP_3"
repeat
wait
ev:= getnextevent
Alg:= getclass(ev)
if Alg=STEP_0 then ! Produce an initial solution for the
! dual problem using a dummy objective
maximize(sum(j in Ctrs) u(j))
if(getprobstat = XPRS_INF) then
writeln("Problem is infeasible")
send(EVENT_INFEAS,0) ! Problem is infeasible
else
write("**** Start solution: ")
save_solution
end-if
else ! STEP 2: Solve the dual problem for
! given solution values of y
initializations from "bin:shmem:sol"
sol_y
end-initializations
Obj:= sum(j in Ctrs) u(j)* (b(j) - sum(i in IntVars) B(j,i)*sol_y(i))
maximize(XPRS_PRI, Obj)
if(getprobstat=XPRS_UNB) then
BigM:= sum(j in Ctrs) u(j) <= BIGM
maximize(XPRS_PRI, Obj)
end-if
write("Step 2: ")
save_solution ! Write solution to memory
BigM:= 0 ! Reset the 'BigM' constraint
end-if
until false
!-----------------------------------------------------------
! Process solution data
procedure save_solution
! Store values of u and x
forall(j in Ctrs) sol_u(j):= getsol(u(j))
forall(i in CtVars) sol_x(i):= getdual(CtrD(i))
initializations to "bin:shmem:sol"
sol_u sol_x
end-initializations
send(EVENT_SOLVED, getobjval)
write(getobjval, "\n u: ")
forall(j in Ctrs) write(sol_u(j), " ")
write("\n x: ")
forall(i in CtVars) write(getdual(CtrD(i)), " ")
writeln
fflush
end-procedure
end-model
|
(!*******************************************************
Mosel Example Problems
======================
file benders_cont.mos
`````````````````````
Benders decomposition for solving a simple MIP.
- Solve the primal problem for given solution values
of the integer variable (Step 2 of the algorithm) -
*** Not intended to be run standalone - run from benders_main.mos ***
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, Jan. 2006, rev. Jan. 2013
*******************************************************!)
model "Benders (continuous problem)"
uses "mmxprs", "mmjobs"
parameters
NCTVAR = 3
NINTVAR = 3
NC = 4
BIGM = 1000
end-parameters
declarations
STEP_0=2 ! Event codes sent to submodels
STEP_2=4
STEP_3=5
EVENT_SOLVED=6 ! Event codes sent by submodels
EVENT_READY=8
CtVars = 1..NCTVAR ! Continuous variables
IntVars = 1..NINTVAR ! Discrete variables
Ctrs = 1..NC ! Set of constraints (orig. problem)
A: array(Ctrs,CtVars) of integer ! Coeff.s of continuous variables
B: array(Ctrs,IntVars) of integer ! Coeff.s of discrete variables
b: array(Ctrs) of integer ! RHS values
C: array(CtVars) of integer ! Obj. coeff.s of continuous variables
Ctr: array(Ctrs) of linctr ! Constraints of orig. problem
sol_u: array(Ctrs) of real ! Solution of dual problem
sol_x: array(CtVars) of real ! Solution of primal prob. (cont.)
sol_y: array(IntVars) of real ! Solution of primal prob. (integers)
x: array(CtVars) of mpvar ! Continuous variables
end-declarations
initializations from "bin:shmem:probdata"
A B b C
end-initializations
Obj:= sum(i in CtVars) C(i)*x(i)
send(EVENT_READY,0) ! Model is ready (= running)
! (Re)solve this model until it is stopped by event "STEP_3"
repeat
wait
dropnextevent
initializations from "bin:shmem:sol"
sol_y
end-initializations
forall(j in Ctrs)
Ctr(j):= sum(i in CtVars) A(j,i)*x(i) +
sum(i in IntVars) B(j,i)*sol_y(i) >= b(j)
minimize(Obj) ! Solve the problem
! Store values of u and x
forall(j in Ctrs) sol_u(j):= getdual(Ctr(j))
forall(i in CtVars) sol_x(i):= getsol(x(i))
initializations to "bin:shmem:sol"
sol_u sol_x
end-initializations
send(EVENT_SOLVED, getobjval)
write("Step 2: ", getobjval, "\n u: ")
forall(j in Ctrs) write(sol_u(j), " ")
write("\n x: ")
forall(i in CtVars) write(getsol(x(i)), " ")
writeln
fflush
until false
end-model
|
(!*******************************************************
Mosel Example Problems
======================
file benders_int.mos
````````````````````
Benders decomposition for solving a simple MIP.
- Solve the primal problem for given solution values
of the continuous variables (Step 1 of the algorithm) -
*** Not intended to be run standalone - run from benders_main.mos ***
(c) 2008 Fair Isaac Corporation
author: S. Heipcke, Jan. 2006, rev. Jan. 2013
*******************************************************!)
model "Benders (integer problem)"
uses "mmxprs", "mmjobs"
parameters
NINTVAR = 3
NC = 4
BIGM = 1000
end-parameters
declarations
STEP_0=2 ! Event codes sent to submodels
STEP_1=3
EVENT_SOLVED=6 ! Event codes sent by submodels
EVENT_READY=8
IntVars = 1..NINTVAR ! Discrete variables
Ctrs = 1..NC ! Set of constraints (orig. problem)
B: array(Ctrs,IntVars) of integer ! Coeff.s of discrete variables
b: array(Ctrs) of integer ! RHS values
D: array(IntVars) of integer ! Obj. coeff.s of discrete variables
MC: array(range) of linctr ! Constraints generated by alg.
sol_u: array(Ctrs) of real ! Solution of dual problem
sol_y: array(IntVars) of real ! Solution of primal prob.
y: array(IntVars) of mpvar ! Discrete variables
z: mpvar ! Objective function variable
end-declarations
initializations from "bin:shmem:probdata"
B b D
end-initializations
z is_free ! Objective is a free variable
forall(i in IntVars) y(i) is_integer ! Integrality condition
forall(i in IntVars) y(i) <= BIGM ! Set (large) upper bound
send(EVENT_READY,0) ! Model is ready (= running)
repeat
wait
ev:= getnextevent
ct:= integer(getvalue(ev))
initializations from "bin:shmem:sol"
sol_u
end-initializations
! Add a new constraint
MC(ct):= z >= sum(i in IntVars) D(i)*y(i) +
sum(j in Ctrs) sol_u(j)*(b(j) - sum(i in IntVars) B(j,i)*y(i))
minimize(z)
! Store solution values of y
forall(i in IntVars) sol_y(i):= getsol(y(i))
initializations to "bin:shmem:sol"
sol_y
end-initializations
send(EVENT_SOLVED, getobjval)
write("Step 1: ", getobjval, "\n y: ")
forall(i in IntVars) write(sol_y(i), " ")
write("\n Slack: ")
forall(j in 1..ct) write(getslack(MC(j)), " ")
writeln
fflush
until false
end-model
|
(!*******************************************************
Mosel Example Problems
======================
file benders_single.mos
```````````````````````
Benders decomposition for solving a simple MIP.
- Single model containing several subproblems -
(c) 2009 Fair Isaac Corporation
author: S. Heipcke, Jan. 2009, rev. Jul. 2020
*******************************************************!)
model "Benders (single model)"
uses "mmxprs", "mmjobs"
parameters
NCTVAR = 3 ! 1
NINTVAR = 3 ! 2
NC = 4 ! 3
BIGM = 1000
ALG = 1 ! 1: Use Benders decomposition (dual)
! 2: Use Benders decomposition (primal)
DATAFILE = "bprob33.dat" ! "bprob12.dat"
end-parameters
forward procedure start_solution
forward procedure solve_primal_int(ct: integer)
forward procedure solve_cont
forward function eval_solution: boolean
forward procedure print_solution
forward procedure define_intprob(prob:mpproblem)
forward function solve_intprob(prob:mpproblem, ct:integer): real
forward procedure define_dualprob(prob:mpproblem)
forward function solve_dualprob(prob:mpproblem, Alg:integer): real
forward procedure define_contprob(prob:mpproblem)
forward function solve_contprob(prob:mpproblem): real
declarations
STEP_0=2 ! Codes sent to subproblems
STEP_1=3
STEP_2=4
STAT_SOLVED=6 ! Status codes returned by subproblems
STAT_INFEAS=7
STAT_READY=8
CtVars = 1..NCTVAR ! Continuous variables
IntVars = 1..NINTVAR ! Discrete variables
Ctrs = 1..NC ! Set of constraints (orig. problem)
A: array(Ctrs,CtVars) of integer ! Coeff.s of continuous variables
B: array(Ctrs,IntVars) of integer ! Coeff.s of discrete variables
b: array(Ctrs) of integer ! RHS values
C: array(CtVars) of integer ! Obj. coeff.s of continuous variables
D: array(IntVars) of integer ! Obj. coeff.s of discrete variables
Ctr: array(Ctrs) of linctr ! Constraints of orig. problem
CtrD: array(CtVars) of linctr ! Constraints of dual problem
MC: array(range) of linctr ! Constraints generated by alg.
Obj: linctr ! Continuous objective
y: array(IntVars) of mpvar ! Discrete variables
z: mpvar ! Objective function variable
u: array(Ctrs) of mpvar ! Dual variables
x: array(CtVars) of mpvar ! Continuous variables
sol_u: array(Ctrs) of real ! Solution of dual problem
sol_x: array(CtVars) of real ! Solution of primal prob. (cont.)
sol_y: array(IntVars) of real ! Solution of primal prob. (integers)
sol_obj: real ! Objective function value (primal)
RM: range ! Problem indices
stepprob: dynamic array(RM) of mpproblem ! Subproblems
status: dynamic array(mpproblem) of integer ! Subproblem status
end-declarations
initializations from DATAFILE
A B b C D
end-initializations
! **** Subproblems ****
! Create and define submodels
create(stepprob(1)); define_intprob(stepprob(1))
if ALG=1 then
create(stepprob(2)); define_dualprob(stepprob(2))
else
create(stepprob(0)); define_dualprob(stepprob(0))
create(stepprob(2)); define_contprob(stepprob(2))
end-if
! **** Solution algorithm ****
start_solution ! 0. Initial solution for cont. var.s
ct:= 1
repeat
writeln("\n**** Iteration: ", ct)
solve_primal_int(ct) ! 1. Solve problem with fixed cont.
solve_cont ! 2. Solve problem with fixed int.
ct+=1
until eval_solution ! Test for optimality
print_solution ! 3. Retrieve and display the solution
!-----------------------------------------------------------
! Produce an initial solution for the dual problem using a dummy objective
procedure start_solution
num:= if(ALG=1, 2, 0)
res:=solve_dualprob(stepprob(num), STEP_0) ! Start the problem solving
if status(stepprob(num))=STAT_INFEAS then
writeln("Problem is infeasible")
exit(6)
end-if
end-procedure
!-----------------------------------------------------------
! Solve a modified version of the primal problem, replacing continuous
! variables by the solution of the dual
procedure solve_primal_int(ct: integer)
sol_obj:= solve_intprob(stepprob(1), ct)
end-procedure
!-----------------------------------------------------------
! Solve the Step 2 problem (dual or primal depending on value of ALG)
! for given solution values of y
procedure solve_cont
if ALG=1 then ! Start the problem solving
res:= solve_dualprob(stepprob(2), STEP_2)
else
res:= solve_contprob(stepprob(2))
end-if
end-procedure
!-----------------------------------------------------------
function eval_solution: boolean
write("Test optimality: ", sol_obj - sum(i in IntVars) D(i)*sol_y(i),
" >= ", sum(j in Ctrs) sol_u(j)* (b(j) -
sum(i in IntVars) B(j,i)*sol_y(i)) )
returned:= ( sol_obj - sum(i in IntVars) D(i)*sol_y(i) >=
sum(j in Ctrs) sol_u(j)* (b(j) - sum(i in IntVars) B(j,i)*sol_y(i)) )
writeln(if(returned, " : true", " : false"))
end-function
!-----------------------------------------------------------
procedure print_solution
! Reset all submodels
forall(m in RM) reset(stepprob(m))
write("\n**** Solution (Benders): ", sol_obj, "\n x: ")
forall(i in CtVars) write(sol_x(i), " ")
write(" y: ")
forall(i in IntVars) write(sol_y(i), " ")
writeln
end-procedure
!-----------------------------------------------------------
! Define the integer problem
procedure define_intprob(prob:mpproblem)
with prob do
z is_free ! Objective is a free variable
forall(i in IntVars) y(i) is_integer ! Integrality condition
forall(i in IntVars) y(i) <= BIGM ! Set (large) upper bound
end-do
status(prob):= STAT_READY
end-procedure
! Solve the integer problem
function solve_intprob(prob:mpproblem, ct:integer): real
with prob do
status(prob):= STAT_READY
! Add a new constraint
MC(ct):= z >= sum(i in IntVars) D(i)*y(i) +
sum(j in Ctrs) sol_u(j)*(b(j) - sum(i in IntVars) B(j,i)*y(i))
minimize(z)
! Store solution values of y
forall(i in IntVars) sol_y(i):= getsol(y(i))
returned:=getobjval
status(prob):= STAT_SOLVED
write("Step 1: ", getobjval, "\n y: ")
forall(i in IntVars) write(sol_y(i), " ")
write("\n Slack: ")
forall(j in 1..ct) write(getslack(MC(j)), " ")
writeln
fflush
end-do
end-function
!-----------------------------------------------------------
! Define the dual problem
procedure define_dualprob(prob:mpproblem)
with prob do
forall(i in CtVars) CtrD(i):= sum(j in Ctrs) u(j)*A(j,i) <= C(i)
end-do
status(prob):= STAT_READY
end-procedure
! Process dual solution data
procedure save_dualsolution(prob:mpproblem)
! There is no need for 'with mpproblem do' here since this subroutine
! is called from within a subproblem
! Store values of u and x
forall(j in Ctrs) sol_u(j):= getsol(u(j))
forall(i in CtVars) sol_x(i):= getdual(CtrD(i))
status(prob):= STAT_SOLVED
write(getobjval, "\n u: ")
forall(j in Ctrs) write(sol_u(j), " ")
write("\n x: ")
forall(i in CtVars) write(getdual(CtrD(i)), " ")
writeln
fflush
end-procedure
! (Re)solve the dual problem
function solve_dualprob(prob:mpproblem, Alg:integer): real
with prob do
status(prob):= STAT_READY
if Alg=STEP_0 then ! Produce an initial solution for the
! dual problem using a dummy objective
maximize(sum(j in Ctrs) u(j))
if(getprobstat = XPRS_INF) then
writeln("Problem is infeasible")
status(prob):= STAT_INFEAS ! Problem is infeasible
else
write("**** Start solution: ")
save_dualsolution(prob)
returned:= getobjval
end-if
else ! STEP 2: Solve the dual problem for
! given solution values of y
Obj:= sum(j in Ctrs) u(j)* (b(j) - sum(i in IntVars) B(j,i)*sol_y(i))
maximize(XPRS_PRI, Obj)
if(getprobstat=XPRS_UNB) then
BigM:= sum(j in Ctrs) u(j) <= BIGM
maximize(XPRS_PRI, Obj)
end-if
write("Step 2: ")
save_dualsolution(prob) ! Save solution values
returned:= getobjval
BigM:= 0 ! Reset the 'BigM' constraint
end-if
end-do
end-function
!-----------------------------------------------------------
! Define the continuous primal problem
procedure define_contprob(prob:mpproblem)
with prob do
Obj:= sum(i in CtVars) C(i)*x(i)
end-do
status(prob):= STAT_READY
end-procedure
! Solve the continuous problem
function solve_contprob(prob:mpproblem): real
with prob do
status(prob):= STAT_READY
! (Re)define and solve this problem
forall(j in Ctrs)
Ctr(j):= sum(i in CtVars) A(j,i)*x(i) +
sum(i in IntVars) B(j,i)*sol_y(i) >= b(j)
minimize(Obj) ! Solve the problem
! Store values of u and x
forall(j in Ctrs) sol_u(j):= getdual(Ctr(j))
forall(i in CtVars) sol_x(i):= getsol(x(i))
returned:=getobjval
status(prob):= STAT_SOLVED
write("Step 2: ", getobjval, "\n u: ")
forall(j in Ctrs) write(sol_u(j), " ")
write("\n x: ")
forall(i in CtVars) write(getsol(x(i)), " ")
writeln
fflush
end-do
end-function
!-----------------------------------------------------------
end-model
|
