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

   file openshop.mos
   `````````````````
   TYPE:         Preemptive open shop scheduling
   DIFFICULTY:   5
   FEATURES:     MIP problem, data preprocessing, algorithm for preemptive 
                 scheduling that involves looping over optimization,
                 `Gantt chart' printing and drawing
   DESCRIPTION:  A satellite routes the traffic from four transmitter 
                 stations to four receiver stations. It has a switch that 
                 allows any permutation (mode) between the four transmitters 
                 and the four receivers. A valid schedule of transmissions 
                 defines a sequence of permutations that routes all the 
                 traffic of the given traffic matrix TRAF. An element of 
                 TRAF may be fragmented, and appear in several modes of 
                 the final decomposition. The objective is to find a 
                 schedule with minimal total duration.
                 The solution algorithm solves an assignment problem for 
                 every permutation.   
   FURTHER INFO: `Applications of optimization with Xpress-MP', 
                 Section 12.5 `Scheduling of telecommunications via satellite'  
   
   (c) 2008 Fair Isaac Corporation
       author: S. Heipcke, 2002, rev. Sep. 2017
*******************************************************!)

model "Satellite scheduling"
 uses "mmxprs", "mmsvg"

 declarations
  TRANSM = 1..4                         ! Set of transmitters
  RECV = 1..4                           ! Set of receivers

  TRAF: array(TRANSM,RECV) of integer   ! Traffic betw. terrestrial stations
  TQBS: array(TRANSM,RECV) of integer   ! Quasi bistochastic traffic matrix
  row: array(TRANSM) of integer         ! Row sums
  col: array(RECV) of integer           ! Column sums
  LB: integer                           ! Maximum of row and column sums
 end-declarations

 initializations from 'openshop.dat'
  TRAF
 end-initializations

! Row and column sums
 forall(t in TRANSM) row(t):= sum(r in RECV) TRAF(t,r)
 forall(r in RECV) col(r):= sum(t in TRANSM) TRAF(t,r)
 LB:=maxlist(max(r in RECV) col(r), max(t in TRANSM) row(t))

! Calculate TQBS
 forall(t in TRANSM,r in RECV) do
  q:= minlist(LB-row(t),LB-col(r))
  TQBS(t,r):= TRAF(t,r)+q
  row(t)+=q
  col(r)+=q
 end-do

 declarations
  MODES: range
  
  flow: array(TRANSM,RECV) of mpvar         ! 1 if transmission from t to r,
                                            ! 0 otherwise
  pmin: mpvar                               ! Minimum exchange
  onerec, minexchg: array(TRANSM) of linctr ! Constraints on transmitters
                                            ! and min exchange
  onetrans: array(RECV) of linctr           ! Constraints on receivers

  solflowt: array(TRANSM,MODES) of integer  ! Solutions of every iteration
  solflowr: array(RECV,MODES) of integer    ! Solutions of every iteration
  solpmin: array(MODES) of integer          ! Objective value per iteration
 end-declarations

 forall(t in TRANSM,r in RECV) flow(t,r) is_binary

 ct:= 0
 while(sum(t in TRANSM,r in RECV) TQBS(t,r) > 0) do
  ct+=1
  
 ! One receiver per transmitter
  forall(t in TRANSM) onerec(t):= sum(r in RECV | TQBS(t,r)>0) flow(t,r) =1
 ! One transmitter per receiver
  forall(r in RECV) onetrans(r):= sum(t in TRANSM | TQBS(t,r)>0) flow(t,r) =1

 ! Minimum exchange
  forall(t in TRANSM) 
   minexchg(t):= sum(r in RECV | TQBS(t,r)>0) TQBS(t,r)*flow(t,r) >= pmin

 ! Solve the problem: maximize the minimum exchange
  maximize(pmin)

 ! Solution printing
  writeln("Round ", ct, " objective: ", getobjval)

 ! Save the solution
  solpmin(ct):= round(getobjval)
  forall(t in TRANSM,r in RECV | TQBS(t,r)>0)
   if(getsol(flow(t,r))>0) then
    solflowt(t,ct):= t
    solflowr(t,ct):= r
   end-if 

 ! Update TQBS
  forall(t in TRANSM)
   TQBS(solflowt(t,ct),solflowr(t,ct)) -= solpmin(ct)

 end-do
 
! Solution printing
 writeln("\nTotal duration: ", sum(m in MODES) solpmin(m))

 write("      ")
 forall(i in 0..ceil(LB/5)) write(strfmt(i*5,5))
 writeln
 forall(t in TRANSM) do
  write("From ", t, " to: ")
  forall(m in MODES)
   forall(i in 1..solpmin(m)) do
    write(if(TRAF(solflowt(t,m),solflowr(t,m))>0, string(solflowr(t,m)),"-"))
    TRAF(solflowt(t,m),solflowr(t,m))-=1
   end-do
  writeln 
 end-do

! Solution drawing
 declarations
  TGraph: array(RECV) of string
 end-declarations

 initializations from 'openshop.dat'
  TRAF
 end-initializations

 FACT:=3
 svgsetgraphviewbox(1, 0, sum(m in MODES) solpmin(m)+1, FACT*TRANSM.size+1)
 svgsetgraphscale(FACT)
 svgsetgraphpointsize(FACT)
 svgsetgraphlabels("Time", "Transmitters")

 forall(t in TRANSM) do
   TGraph(t):="R"+t
   svgaddgroup(TGraph(t), "Receiver "+t) 
 end-do
 forall(t in TRANSM) do
  ct:=0
  forall(m in MODES) do
   forall(i in 1..solpmin(m)) do
    if(TRAF(solflowt(t,m),solflowr(t,m))>0) then
     svgaddpoint(TGraph(round(solflowr(t,m))), ct+i, t*FACT)
    end-if
    TRAF(solflowt(t,m),solflowr(t,m))-=1
   end-do
   ct+=solpmin(m)
  end-do
 end-do

 svgsave("openshop.svg")
 svgrefresh
 svgwaitclose("Close browser window to terminate model execution.", 1)
end-model