(!**************************************************************************
   Mosel Example Problems
   ====================== 

   file GAPbp3.mos
   ```````````````   
   TYPE:         Branch-and-Price for Generalized Assignment Problem (GAP)
   DIFFICULTY:   5
   FEATURES:     mmjobs, Mosel 2.0 records.
   DESCRIPTION:  A set of production batches is assigned to machines with the
                 objective to maximize profit. The model implements a branch-
                 and-price algorithm that solves a disaggregated formulation of
                 the GAP where columns represent feasible assignments of batches
                 to machines. Column generation is applied at every node of the
                 branch-and-bound tree. Branching is based on fixing compact 
                 formulation variables that represent assignments of production
                 batches to machines. The branch-and-bound logic is 
                 implemented in Mosel, and the optimizer is used only to solve
                 the LP relaxation at each node. 
                
   FURTHER INFO: Similar Problems:
                 Mosel Parallel Whitepaper
                   3.1 Example problem: cutting stock
                   5.1 Example problem: multi-item, multi-period production 
                       planning
                 References for Branch-and-Price:
                 Savelsbergh, M.W.P.: "A Branch-and-Price Algorithm for the 
                   Generalized Assignment Problem" Oper. Res., v45, pp. 831-
                   841, 1997. 
                 Savelsbergh, M.W.P.: "Branch-and-price: Integer programming 
                   with column generation." In: Encyclopedia of Optimization,
                   C. Floudas, pp. Pardalos, eds, 2001

   *** ATTENTION: This model will return an error if no   ***
   *** sufficient number of Xpress licences is available. ***

   The input data required by this model is generated by executing
   model "genGapDataD.mos" prior to running this model.
   
   (c) 2008 Fair Isaac Corporation
       author: Hernan Wurgaft, 2007, rev. S. Heipcke Sep. 2022
***************************************************************************!)

model GAPmas
uses "mmxprs", "mmjobs", "mmsystem" 

parameters
  DATAFILE = "Data/Dtestx5_30_1.dat"  
  TOL = 0.00001
  DO_BRANCHING=true
  DO_HEURISTIC=false
  CUTOFF=0
  BRANCH_STRATEGY=0   ! 0 is for variable branching (this is a placeholder for
                      !  other branching strategies)
   
  BIG=1000000
  SUBMOD="Data/GAPsubDP.mos"   
  MESSLEV=4           ! 0: silent, 1-5: progressivly more detailed output
end-parameters

setparam("XPRS_PRESOLVE",0)

declarations                 
  EVENT_GEN_INI_COLS=2        ! Event codes sent to submodels
  EVENT_GEN_COLS=3
  EVENT_STOP_SUBMOD=5
  EVENT_NEW_NODE=9
  
  EVENT_SOLVED=6              ! Event codes sent by submodels
  EVENT_FAILED=7
  EVENT_READY=8
  EVENT_INFEAS=11
end-declarations   


! Declare and read input data
forward procedure read_data_size
forward procedure read_prob_data

declarations 
  NM, NP : integer
end-declarations
 
read_data_size
    
declarations   
  MACH = 1..NM                          ! Set of machines
  PRODS = 1..NP                         ! Set of production batches
  CAP: array(MACH) of integer           ! Machine capacities
  DUR: array(MACH,PRODS) of integer     ! Production durations
  PROFIT: array(MACH,PRODS) of integer  ! Production profit
end-declarations

read_prob_data
  
! Make data available for submodels to read
initializations to "bin:shmem:data"
   DUR
   CAP
   PROFIT
end-initializations
 
   
 !**** Branch-and-Price declarations****
declarations
  excessA: mpvar                        ! Violation of Assignment
  weight: array(MACH,range) of mpvar    ! Weight for each column
  nPROP: array(MACH) of integer         ! Number of columns from each machine
  Assign: array(PRODS) of linctr        ! Assignment constraints
  Convex: array(MACH) of linctr         ! Convexity constraints
  Mobjective: linctr                    ! Objective function of main problem
  USE: array(MACH,range,PRODS) of integer ! Stores columns as they are generated
  Price_convex: array(MACH) of real     ! Dual prices for convexity constraints
  Price_assign: array(PRODS) of real    ! Dual prices for assignment constraints
 
  branchvartype=
   record        ! Data structure representing compact formulation variables:
    m: integer     ! machine
    p: integer     ! product
  end-record    
  
  branchvar: branchvartype
  
  Node=
   record          ! Data structure representing a node in the tree:
                   ! Pointers refer to indexes in dynamic array nodes
     value: real           ! optimal value of relaxation associated with node
     left,right,parent: integer   ! pointers to parent and children in node tree
     feasible_sol: boolean
     integer_sol: boolean              
     solbasis: basis       ! basis associated with optimal value of relaxation
     prev_val, next_val: integer  ! pointers to nodes in linked queue of
                                  ! active nodes (queue is maintained sorted by
                                  ! value)
     branchvar: branchvartype     ! compact formulation variable to
                                  ! branch at (variable branching)
   end-record
  
  R: range  
  nodes:dynamic array(R) of Node  ! Array containing all nodes
  node: integer    ! Current node being solved by column generation
  head: integer    ! Head of queue of active nodes (0 if queue is empty)
  tail: integer    ! Tail of queue of active nodes (0 if queue is empty)
  Select_node:integer             ! Active node selected for branching
  
  Fixed=
   record          ! Data structure representing fixed variables:
     var:    branchvartype        ! compact formulation variable that is fixed
     to_zero: boolean             ! true if fixed to zero, false if fixed to one
   end-record
  fixed_vars: array(1..NP) of Fixed     ! Array with all fixed variables
  
  fixedcolumns: set of integer
    
  Nfixed: integer  ! Number of compact formulation variables fixed
  NWfix: integer   ! Count of constraints to fix variables
  consfix: array(range) of linctr       ! Constraints to fix variables
   
  sol_BEST: array(PRODS) of integer     ! Incumbent solution in terms of compact
                                        ! formulation variables
  sol_GAP: array(MACH,PRODS) of real    ! Solution in terms of compact
                                        ! formulation variables
  
  ubound: real                          ! Upper bound from Lagrangian relaxation
                                        ! calculated at every iteration
  Bestubound: real                      ! Best upper bound from Lagrangian
                                        ! relaxation
  Best_Value: real                      ! Incumbent value
  masbasis: basis                       ! Used to save basis at every column
                                        ! generation iteration

  Totiter: integer                      ! Count of column generation iterations
  countCols: integer                    ! Count of columns in disaggregated
                                        ! formulation
  integer_solution: boolean             ! true if latest node solution is
                                        ! integer (set by translate_solution)
end-declarations

! Column generation procedures
forward procedure solve_node 
forward procedure process_sub_result 
forward procedure solve_main
forward procedure process_main_result  

! Branch-and-Price supporting procedures
forward procedure translate_solution
forward procedure branching_variable 
forward procedure variable_branching_reset
forward procedure fixvar(branchvar:branchvartype,to_zero:boolean)   
! (functions used to manage set fixedcolumns):
forward function encode(m:integer,k:integer):integer
forward function getm(e:integer):integer
forward function getk(e:integer):integer

! Placeholder for heuristic
forward procedure heuristic    
        
! Reporting procedures
forward procedure save_best_solution
forward procedure write_sol_GAP         ! Writes node solution in terms of
                                        ! compact representation variables
forward procedure write_sol_BEST        ! Writes incumbent solution in terms of
                                        ! compact representation variables
forward procedure write_new_col(Machine:integer)  
forward procedure write_msg(msglevel:integer, msg:string)

!*************************************************
!*************************************************
!***************  Get the Sub-models running
!*************************************************
!*************************************************
 
declarations
  submod: array(MACH) of Model             ! One submodel per machine
  modid: array(set of integer) of integer  ! Model indices
end-declarations 
 
res:= compile("g",SUBMOD,"tmp:gapsub.bim") ! Compile the submodel file

forall(m in MACH) do                       ! Load & run submodels
  load(submod(m), "tmp:gapsub.bim")
  modid(getid(submod(m))):= m
  run(submod(m), "Machine=" + m + ",NM=" + NM + ",NP=" + NP )
                !!All submodels are started in sequence
                !!They will wait until they receive signal to continue
    
  wait                                     ! Wait for submodels to be ready
  ev:=getnextevent
  if ev.class=EVENT_END then
   writeln("*** Impossible to start all necessary models - aborting ***")
   exit(1)
  end-if
end-do
 
!*************************************************
!*************************************************
!***************  Solve Root Node
!*************************************************
!*************************************************
declarations
  Machine: integer
  sol_use: array(PRODS) of integer
  sol_Profit: real
end-declarations

StartTime:= gettime
 
! Generate initial columns
!*************************
forall(m in MACH) send(submod(m),EVENT_GEN_INI_COLS,0) 
                ! All submodels are sent signal to run in parallel
  
forall(m in MACH) do 
                ! This piece uses the memory pipe to wait and then processes
                ! the solution of each sub-model to generate initial columns
  initializations from "mempipe:noindex,sol"
     Machine
     sol_use
     sol_Profit
  end-initializations
  
 ! Add the new column to the main problem
  nPROP(Machine)+=1
  countCols+=1
  create(weight(Machine,nPROP(Machine)))
  weight(Machine,nPROP(Machine)) is_binary
  forall(p in PRODS) USE(Machine,nPROP(Machine),p):=sol_use(p)
        
  write_new_col(Machine)
end-do
  
 ! Catches event that submodels completed generation of initial columns
forall(m in MACH) do 
  wait;   dropnextevent;  
end-do
  
! Define and solve main problem with initial columns
!***************************************************
forall(m in MACH)
  Convex(m):= sum (k in 1..nPROP(m)) weight(m,k) <= 1  ! Convexity constraints
 
forall(p in PRODS)
  Assign(p):= sum(m in MACH, k in 1..nPROP(m)) USE(m,k,p) * 
                            weight(m,k)+ excessA = 1   ! Assignment constraints
                                                            
Mobjective:= sum(m in MACH, k in 1..nPROP(m))
                (sum (p in PRODS) USE(m,k,p)* PROFIT(m,p))* 
                 weight(m,k)- 100 * excessA            ! Objective function
  
write_msg(2, "Solving initial main problem: " + (gettime-StartTime) + "sec")
maximize(XPRS_LIN, Mobjective)  
savebasis(masbasis)  
  
write_msg(1, "Iter: 0  sol: " + getobjval)    
 
process_main_result
Best_Value:=CUTOFF
  
solve_node
 
!*************************************************
!*************************************************
!*************************************************
!*************************************************
!*************************************************
!*************************************************
!***************  Begin Branch-and-Price
!*************************************************
!*************************************************
!*************************************************
!*************************************************
!*************************************************
!*************************************************
 
if DO_BRANCHING then 
 !*************************************************
 !*************************************************
 !***  Initialize root node in search tree active nodes
 !*************************************************
 !*************************************************

  write_msg(2, "Start branching: " + (gettime-StartTime) + "sec")
  node                   := 0
  NWfix                  := 0
  Nfixed                 := 0
  translate_solution
  
  if integer_solution then
    save_best_solution
  else
    if BRANCH_STRATEGY=0: branching_variable  
    if DO_HEURISTIC: heuristic
  end-if
 
  if getsol(excessA)>0 then
    nodes(node).feasible_sol:=false
    write_msg(1, "Infeasible")
  elif integer_solution then
    write_msg(1, "Integer Solution")
    nodes(node).integer_sol:=true
  else
    nodes(node).feasible_sol:=true  
    nodes(node).integer_sol:=false
    nodes(node).value:=getobjval
    nodes(node).parent:=-1
    nodes(node).prev_val:=-1
    nodes(node).next_val:=-1
    head:=0
    tail:=0
    Select_node:=0
    nodes(node).left:= -1
    nodes(node).right:=-1
  end-if
    
  if nodes(0).feasible_sol and not nodes(0).integer_sol then
 
!*****************************************************
!***************  Begin looping over search tree nodes
!*****************************************************
 
    search_done:=false
    repeat        ! until search_done
      node+=1
      nodes(node).parent:=Select_node
      nodes(node).left:= -1
      nodes(node).right:=-1
 
 !*****************************************************
 !** Reset the main problem for the current node.
 !*****************************************************
(! Main problem columns that are fixed to zero for current node
   are first put in set fixedcolumns. Constraints consfix are
   used to fix the columns to zero. Also, information about
   fixed variables is passed to knapsack submodels.
   Procedures used:
    variable_branching_reset
    fixvar(branchvar:branchvartype,to_zero:boolean)
   Functions used (to manage set fixedcolumns):
    function encode(m:integer,k:integer):integer
    function getm(e:integer):integer
    function getk(e:integer):integer 
 !)
 
      forall(c in 1..NWfix) consfix(c):=0
      NWfix:=0
      Nfixed:=0
      fixedcolumns:={}        
 
      if BRANCH_STRATEGY=0: variable_branching_reset
      initializations to "bin:shmem:fixed"
        Nfixed
        fixed_vars
      end-initializations
 
      forall(col in fixedcolumns) do
        NWfix+=1
        consfix(NWfix):=weight(getm(col),getk(col))=0
      end-do
 
      forall(m in MACH) send(submod(m), EVENT_NEW_NODE, 0)
  
 !*****************************************************
 !** Optimize main problem
 !*****************************************************
      loadprob(Mobjective)                    ! Reload the problem
      loadbasis(nodes(Select_node).solbasis)  ! Load saved basis
      maximize(XPRS_LIN,Mobjective)  
      savebasis(masbasis)
      process_main_result
      solve_node 
      translate_solution
 
      if not integer_solution then
   !Determine how to branch off the current node        active nodes
        if BRANCH_STRATEGY=0:  branching_variable  
        if DO_HEURISTIC: heuristic
      end-if
 
 !*****************************************************
 !** Update search tree based on solution of current node
 !** (Update status of all nodes using new information)
 !*****************************************************
      if getsol(excessA)>0 then
        nodes(node).feasible_sol:=false
        write_msg(3, "Infeasible")
      else
        nodes(node).feasible_sol:=true       
        nodes(node).value:=Bestubound
        savebasis(nodes(node).solbasis)
        if integer_solution  then
          write_msg(2, "Integer Solution Found")
          nodes(node).value:=getobjval
          nodes(node).integer_sol:=true
          if nodes(node).value>Best_Value then
            save_best_solution
            Best_Value:=nodes(node).value
            write_msg(2, "Improved integer value: " + Best_Value)          
            if head=-1 or nodes(head).value-Best_Value<1.0 then
              search_done:=true
            end-if
          end-if
        else   !Current node becomes active in the tree
          if nodes(node).value>=Best_Value+1 then
            position:=head
       !Tail
            repeat
              if position<>-1 and nodes(node).value < nodes(position).value  then
                position:=nodes(position).next_val
              else
          !p1 will precede node and position will follow
          !p1 --- node --- position
                if position<>-1 then
                  p1:=nodes(position).prev_val
                  nodes(position).prev_val:=node
                else
                  p1:=tail
                  tail:=node
                end-if                  
                if p1<>-1 then
                  nodes(p1).next_val:=node
                else
                  head:=node    
                end-if
                                                
                nodes(node).prev_val:=p1
                nodes(node).next_val:=position
                        
               break
              end-if
            until false
          else
            write_msg(3,"Cutoff node: " + node)
          end-if  
        end-if       
      end-if
      if MESSLEV>3 then
        writeln("************NODE REPORT*****************")
        sn:=node
        writeln("node ",sn," Value ",nodes(sn).value) 
        sn:=Select_node
        writeln("parent ", sn, " Value ", nodes(sn).value, " P ",
	        nodes(sn).parent," L ",
         nodes(sn).left," R ",nodes(sn).right) 
        writeln("branch m ", nodes(sn).branchvar.m, " branch p ",
	        nodes(sn).branchvar.p) 
        writeln("time elapsed ", gettime, " total iter: ", Totiter,
	        " total cols: ", countCols)
        writeln("******************************************")
      end-if
  
 !*****************************************************
 !** Select next node using Best Bound Rule
 !*****************************************************
      if head=-1 or nodes(head).value-Best_Value<1.0 then
        search_done:=true
      else
        Select_node:=head
      end-if
 
    until search_done
  end-if    !if root node is not infeasible or integer
end-if      !if DO_BRANCHING
!*************************************************
!*************************************************
  
  forall(m in MACH) send(submod(m), EVENT_STOP_SUBMOD,0)  ! Stop all models
   
  !Catch event indicating that all submodels stopped.
  forall(m in MACH) do 
    wait;   dropnextevent;  
  end-do


  write_msg(1, "Total time: " + (gettime-StartTime) + "sec")
  write_msg(1, "Optimal objective value: " + Best_Value)
  write_sol_BEST

!****************************************************
!**************PROCEDURES************************
!****************************************************

procedure solve_node
 (!  Solve node relaxation applying column generation
  ***************************************************
  At each iteration, a knapsack submodel is solved for each
  machine. Knapsack solutions produce columns for the main
  problem and they are also used to calculate a Lagrangean
  upper bound. The loop is broken when no knapsack solution
  produces a column that improves the objective function 
  (solved to optimality), when the main problem objective value 
  is close enough to the upper bound (solved to almost optimality),
  or when the upper bound is smaller than the cutoff value
  (Best bound implies cutoff).
  
  Uses procedures:
   process_sub_result
   solve_main  
   process_main_result
  !)
 Bestubound:=BIG
 Niter:=0
 
 repeat
  Niter+=1
  forall(m in MACH)                        ! Start solving all submodels
    send(submod(m), EVENT_GEN_COLS, 0)
  improve:=0
  
  forall(m in MACH) do 
    wait
    ev:=getnextevent
    if getclass(ev)<>EVENT_INFEAS then
      ubound+=getvalue(ev)
    else
      write_msg(3, "Infeasible knapsack")
    end-if       !!(useful message)
    if getclass(ev)=EVENT_SOLVED then
      improve+=1
      process_sub_result      
    end-if
  end-do
  
    
  if ubound<Bestubound: Bestubound:=ubound
  
  if MESSLEV>2: writeln("Iter: ", Niter, " sol: ", getobjval, " ubound: ",
                        ubound, " Bestubound ", Bestubound)
  
  if improve = 0  then
    write_msg(3, "solved to optimality")
    break
  end-if
  
  if Bestubound<Best_Value+1 then
    write_msg(3, "Bestbound implies cutoff")
    break
  end-if 
  
  !If bound and value are close enough
  if  integer(Bestubound)=integer(getobjval) then
    write_msg(3, "solved to almost optimality")
    break
  end-if
  solve_main  
  process_main_result
 until false
 Totiter+=Niter
end-procedure  

!****************************************************
procedure process_sub_result
 ! Add the new column to the main problem
 
  initializations from "mempipe:noindex,sol"
    Machine
    sol_use
    sol_Profit
  end-initializations
    
  nPROP(Machine)+=1
  countCols+=1
  create(weight(Machine,nPROP(Machine)))
  weight(Machine,nPROP(Machine)) is_binary
  
  forall(p in PRODS) USE(Machine,nPROP(Machine),p):=sol_use(p)
    
  Convex(Machine)+=weight(Machine,nPROP(Machine))
  forall(p in PRODS)
        Assign(p)+= USE(Machine,nPROP(Machine),p) * 
                                  weight(Machine,nPROP(Machine))  
  Mobjective+=(sum (p in PRODS)USE(Machine,nPROP(Machine),p)* 
                         PROFIT(Machine,p))* weight(Machine,nPROP(Machine))
   
  write_new_col(Machine)
end-procedure

!****************************************************
procedure solve_main
 ! Optimize main problem
  loadprob(Mobjective)                    ! Reload the problem
  loadbasis(masbasis)                     ! Load the saved basis
  maximize(XPRS_LIN,Mobjective) 
end-procedure

!****************************************************
procedure process_main_result
 (! Update dual pricing data for subproblems
    Add dual prices to upper bound calculation 
 !)
  ubound:=0
  savebasis(masbasis) 
  forall(p in PRODS) do
    Price_assign(p):=getdual(Assign(p))
    ubound+=Price_assign(p)
  end-do
  forall(m in MACH) do
    Price_convex(m):=getdual(Convex(m))
    ubound+=Price_convex(m)
  end-do
   
  initializations to "bin:shmem:Price"
    Price_assign
    Price_convex
  end-initializations
end-procedure

!****************************************************
procedure translate_solution
 (! Translate solution from disaggregated representation variables (weight)  
    to compact representation variables (sol_GAP). 
    Check whether current solution is integer (integer_solution)
 !)  
  integer_solution:=true
  if MESSLEV>1: writeln("Node ", node, " Objective: ",getobjval)
  forall(m in MACH) do                         ! Loop over machines
    forall(p in PRODS) sol_GAP(m,p):=0.
    forall(k in 1..nPROP(m)) do                ! Loop over columns for machine
      if (getsol(weight(m,k))> TOL) then       ! Finds active columns
        if abs(getsol(weight(m,k))-1.)>TOL then
          integer_solution:=false
        end-if
    
        forall(p in PRODS) do
          if (USE(m,k,p)=1) then
            sol_GAP(m,p)+=getsol(weight(m,k))  ! Translates variables
          end-if
         end-do
      end-if
    end-do
  end-do
end-procedure
 
!****************************************************
procedure branching_variable
 (! Determine the compact formulation variable to branch on for 
    the current node   
 !)
  branch_use    :=1.
  branch_price  :=0.
  forall(m in MACH) do                         ! Loop over machines
    forall(p in PRODS) do    
      if (abs(sol_GAP(m,p)-0.5)<abs(branch_use-0.5)) then
        nodes(node).branchvar.m:=m
        nodes(node).branchvar.p:=p
        branch_use:=sol_GAP(m,p)
        branch_price:=PROFIT(m,p)
      elif (abs(sol_GAP(m,p)-0.5)=abs(branch_use-0.5)) and
           (PROFIT(m,p)>branch_price) then
        nodes(node).branchvar.m:=m
        nodes(node).branchvar.p:=p
        branch_use:=sol_GAP(m,p)
        branch_price:=PROFIT(m,p)
      end-if
    end-do
  end-do
end-procedure

!****************************************************
procedure variable_branching_reset
 (! Determine all compact formulation variables that need to be fixed
    at the current node. Right children are fixed to zero and left 
    children are fixed to one.
     uses procedure:
     fixvar
 !)
  branchvar:= nodes(Select_node).branchvar
          
 !Fix variables
  if nodes(Select_node).left=-1 then           ! Left node
    fixvar(branchvar,false)          
    nodes(Select_node).left:=node
  else                                         ! Right node
    fixvar(branchvar,true)
    nodes(Select_node).right:=node
    !!and take selected node out of queue
    head:=nodes(Select_node).next_val
    nodes(head).prev_val:=-1
  end-if
  
!Fix variable for all the parent nodes in the tree.
  visit:=Select_node
  while (nodes(visit).parent<>-1) do
    parent:=nodes(visit).parent
    branchvar:=nodes(parent).branchvar
    if visit=nodes(parent).left then
      fixvar(branchvar,false)
    else
      fixvar(branchvar,true) 
    end-if
    visit:=parent
  end-do
end-procedure

!****************************************************
procedure fixvar(branchvar:branchvartype, to_zero: boolean)
 (! Build array fixed_vars and set fixedcolumns with compact and disaggregated
   formulation variables to be fixed. Array fixed_vars is to be later passed to 
   knapsack subproblems. The set fixedcolumns is to be used to create consfix
   constraints. Implementation avoids duplicating consfix constraints by using
   function encode to place unique pairs m-k in set fixedcolumns.
 !)
  
  Nfixed+=1
  fixed_vars(Nfixed).var:=branchvar
  
  case to_zero of
    true: do                                    ! Fix to zero
      forall(k in 1..nPROP(branchvar.m)) do
        if USE(branchvar.m,k,branchvar.p)=1 then
          fixedcolumns+={encode(branchvar.m,k)}
        end-if 
      end-do 
      fixed_vars(Nfixed).to_zero:=true
    end-do
 
    false: do                                   ! Fix to one
      fixed_vars(Nfixed).to_zero:=false
      forall(k in 1..nPROP(branchvar.m)) do
        if USE(branchvar.m,k,branchvar.p)<>1 then
          fixedcolumns+={encode(branchvar.m,k)}
        end-if 
      end-do 
   
   !Fix the rest to zero
      forall(m in MACH|m<>branchvar.m) do
        forall(k in 1..nPROP(m)) do
          if USE(m,k,branchvar.p)=1 then
            fixedcolumns+={encode(m,k)}
          end-if 
        end-do
      end-do
    end-do
 
  end-case  
end-procedure

!****************************************************
procedure heuristic
  writeln("Placeholder for heuristic")
end-procedure

!****************************************************
procedure write_sol_GAP
  forall(m in MACH) do
    write("Machine ", m, ": ")
    forall(p in PRODS) do  
      if sol_GAP(m,p)>0 then
        write(p, "(",sol_GAP(m,p),")")
      end-if
    end-do
    writeln
  end-do 
end-procedure  

!****************************************************
procedure write_sol_BEST
  forall(m in MACH) do
    write("Machine ", m, ": ")
    forall(p in PRODS | sol_BEST(p)=m) write(p, " ")
    writeln
  end-do 
end-procedure  

!****************************************************
procedure save_best_solution
  forall(m in MACH,p in PRODS | sol_GAP(m,p)>0) sol_BEST(p):=m
end-procedure 

!****************************************************
procedure write_new_col(Machine:integer)
  if MESSLEV>4 then
    write("Machine ", Machine, ": ")
    forall(p in PRODS | sol_use(p)=1) write(p, " ")
    writeln(" (total dur profit: ", sum(p in PRODS) DUR(Machine,p)*sol_use(p),
                " ",sol_Profit," ", CAP(Machine), ")")
  end-if
end-procedure

!****************************************************
procedure write_msg(msglevel:integer, msg:string)
  if msglevel<= MESSLEV: writeln(msg)
end-procedure 

!****************************************************
function encode(m:integer,k:integer):integer
  returned:=m+k*1024
end-function

function getm(e:integer):integer
  returned:=bittest(e,1023)
end-function

function getk(e:integer):integer
  returned:=e div 1024
end-function
 
!****************************************************
procedure read_data_size
  initializations from DATAFILE
    NM NP
  end-initializations
end-procedure
 
procedure read_prob_data
  initializations from DATAFILE
    DUR CAP PROFIT
  end-initializations
end-procedure

end-model