model "kamele"
uses "mmxprs"
options noimplicit
parameters
	!filename = "karawane2.dat"
	filename = "v150.dat"
	outfile = "mysol.sol"
end-parameters

forward procedure column_gen
forward function knapsack(C:array(range) of real, xbest: array(range) of integer,pass: integer): real
forward procedure show_new_pat(dj:real, vx: array(range) of integer)



declarations

EPS = 1e-6 ! Zero tolerance

PATTERNS: array(range,range) of integer
ITEMS: set of string
DEMANDS: array(range) of integer
WEIGHTS: array(range) of integer
VALUES: array(range) of integer
NUMITEMS: integer
MAXWEIGHT = 100
MAXVALUE = 200
Weights,t1,t2,t3: set of string

use: array(range) of mpvar ! Rolls per pattern
soluse: array(range) of real ! Solution values for variables ‘use’
Dem: array(range) of linctr ! Demand constraints
MinCamels: linctr ! Objective function
KnapCtr1, KnapCtr2, KnapObj: linctr ! Knapsack constraint+objective
x: array(range) of mpvar ! Knapsack variables
npatt,npass: integer
end-declarations

fopen(outfile, F_OUTPUT)

! Read in the set of ITEMS to compute the number, and
! then initialize the dynamic arrays to make them
! big enough

initializations from filename
	ITEMS as 'Items'
end-initializations

NUMITEMS := getsize(ITEMS)

forall (i in 1..NUMITEMS) do
	VALUES(i) := 0
	WEIGHTS(i) := 0
	DEMANDS(i) := 0
end-do

! Now, read in the VALUES, DEMANDS etc. arrays

initializations from filename
	ITEMS as 'Items'
	WEIGHTS as 'Weights'
	VALUES as 'Values'
	DEMANDS as 'Demands'
end-initializations

! forall(j in 1..NUMITEMS) writeln("item ", j, " ",  ITEMS(j), " value: ", VALUES(j), " weight: ", WEIGHTS(j))

! NUMITEMS := getsize(ITEMS)

!writeln("I have read :")
!writeln(NUMITEMS, " items")
!writeln(getsize(WEIGHTS), " weights")
!writeln(getsize(VALUES), " values")
!writeln(getsize(DEMANDS), " demands")


! Make basic patterns
forall(j in 1..NUMITEMS) do
	PATTERNS(j,j):= minlist(floor(MAXWEIGHT/WEIGHTS(j)), floor(MAXVALUE/VALUES(j)))
	!writeln("Basic pattern ", j, " is: ", PATTERNS(j,j))
end-do

forall(i in 1..NUMITEMS) do
	forall (j in 1..NUMITEMS) do
		write(PATTERNS(i,j), " ")
	end-do
	if (i < NUMITEMS) then
		writeln
	end-if
end-do


forall(j in 1..NUMITEMS) do
	create(use(j))
	use(j) is_integer ! Variables are integer and bounded
	use(j) <= integer(ceil(DEMANDS(j)/PATTERNS(j,j)))
	
	create(x(j)) ! create knapsack variable
end-do

MinCamels := sum(j in 1..NUMITEMS) use(j) ! Objective: minimize no. of camels

forall(i in 1..NUMITEMS) ! Satisfy all demands
	Dem(i):= sum(j in 1..NUMITEMS) PATTERNS(i,j) * use(j) >= DEMANDS(i)
	
column_gen ! Column generation at top node
minimize(MinCamels) ! Compute the best integer solution
! for the current problem (including
! the new columns)
!writeln("Best integer solution: ", getobjval, " rolls")
!write(" Rolls per pattern: ")
writeln("Vector:")
forall(i in 1..npatt) write(getsol(use(i))," ")
writeln

! writeln(getobjval)

procedure column_gen
declarations
dualdem: array(1..NUMITEMS) of real
xbest: array(1..NUMITEMS) of integer
dw, zbest, objval: real
end-declarations
setparam("XPRS_CUTSTRATEGY", 0) ! Disable automatic cuts
setparam("XPRS_PRESOLVE", 0) ! Switch presolve off
npatt:=NUMITEMS
npass:=1
while(true) do
	minimize(XPRS_LIN, MinCamels) ! Solve the (relaxed) LP
	savebasis(1) ! 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 1..NUMITEMS) dualdem(i):=getdual(Dem(i))
	! Solve a knapsack problem
	zbest:= knapsack(dualdem, xbest, npass) - 1.0
	!write("Pass ", npass, ": ")
	if(zbest < EPS) then
		!writeln("no profitable column found.\n")
		break
	else
		writeln
		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
		MinCamels+= use(npatt) ! Add new var. to the objective
		dw:=0
		forall(i in 1..NUMITEMS)
			if(xbest(i) > 0) then
				Dem(i)+= xbest(i)*use(npatt) ! Add new var. to demand constr.s
				dw:= maxlist(dw, ceil(DEMANDS(i)/xbest(i) ))
			end-if
		use(npatt) <= dw ! Set upper bound on the new var.
		loadprob(MinCamels) ! Reload the problem
		loadbasis(1) ! 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
end-procedure


function knapsack(C:array(range) of real, xbest:array(range) of integer, pass: integer):real
! Hide the demand constraints
forall(j in 1..NUMITEMS) sethidden(Dem(j), true)

! Define the knapsack problem
KnapCtr1 := sum(j in 1..NUMITEMS) WEIGHTS(j)*x(j) <= MAXWEIGHT
KnapCtr2 := sum(j in 1..NUMITEMS) VALUES(j)*x(j)  <= MAXVALUE

KnapObj := sum(j in 1..NUMITEMS) C(j)*x(j)

!writeln("in knapsack: maximize:")
!forall (i in 1..NUMITEMS) do
	!write("x", i, " * ",C(i), " + ")
!end-do
!writeln

! Integrality condition
if(pass=1) then
	forall(j in 1..NUMITEMS) x(j) is_integer
end-if
maximize(KnapObj)
returned:=getobjval
! write("maximized, value: ", getobjval, " ", C)
forall(j in 1..NUMITEMS) xbest(j):=round(getsol(x(j)))
! Reset knapsack constraint and objective, unhide demand constraints
KnapCtr1 := 0
KnapCtr2 := 0
KnapObj := 0
forall(j in 1..NUMITEMS) sethidden(Dem(j), false)
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(" distribution: ")
!dw:=0
!forall(i in 1..NUMITEMS) do
	!write(WEIGHTS(i), ":", vx(i), " ")
	!dw += WEIGHTS(i)*vx(i)
!end-do
!writeln("Total weight: ", dw)
forall (i in 1..NUMITEMS)
	write(vx(i), " ")
end-procedure
end-model