Initializing help system before first use

Benders decomposition: sequential solving of several different submodels


Type: Programming
Rating: 5 (difficult)
Description: Benders decomposition is a method for solving large MIP problems. The model implementation shows the following features:
  • iterative sequence of concurrent solving of a set of subproblems,
  • data exchange between several models via shared memory, and
  • coordination of several models via events.
An implementation using a single model is also presented (benders_single.mos).
File(s): benders_main.mos, benders_dual.mos, benders_cont.mos, benders_int.mos, benders_single.mos
Data file(s): bprob12.dat, bprob33.dat


benders_main.mos
(!*******************************************************
   Mosel Example Problems 
   ======================

   file benders_main.mos
   `````````````````````
   Benders decomposition for solving a simple MIP.
   - Master 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. Dec. 2017
*******************************************************!)
        
model "Benders (master 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
  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.

  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 then exit(1); end-if
 create(stepmod(1)); load(stepmod(1), "benders_int.bim")
 if compile("benders_dual.mos")<>0 then exit(2); end-if

 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 then exit(3); end-if
  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 then fdelete("benders_cont.bim"); end-if
 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

benders_dual.mos
(!*******************************************************
   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 -

   (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

benders_cont.mos
(!*******************************************************
   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) -

   (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

benders_int.mos
(!*******************************************************
   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) -

   (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

benders_single.mos
(!*******************************************************
   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. Dec. 2017
*******************************************************!)
        
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  
      ! TODO: Reset 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

!----------------------------------------------------------- 
! 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

© 2001-2019 Fair Isaac Corporation. All rights reserved. This documentation is the property of Fair Isaac Corporation (“FICO”). Receipt or possession of this documentation does not convey rights to disclose, reproduce, make derivative works, use, or allow others to use it except solely for internal evaluation purposes to determine whether to purchase a license to the software described in this documentation, or as otherwise set forth in a written software license agreement between you and FICO (or a FICO affiliate). Use of this documentation and the software described in it must conform strictly to the foregoing permitted uses, and no other use is permitted.