## Warehouse Location

This is an example of solving the warehouse problem with Numberjack and SCIP

In the Warehouse Location problem, a company considers opening warehouses

at some candidate locations in order to supply its existing stores. Each possible warehouse

has the same maintenance cost, and a capacity designating the maximum number of stores

that it can supply. Each store must be supplied by exactly one open warehouse.

The supply cost to a store depends on the warehouse. The objective is to determine which

warehouses to open, and which of these warehouses should supply the various stores, such

that the sum of the maintenance and supply costs is minimized. See prob034 in the CSPLib.

def model_warehouse_planning(data):

WareHouseOpen = VarArray(data.NumberOfWarehouses)

ShopSupplied = Matrix(data.NumberOfWarehouses,

data.NumberOfShops)

# Cost of running warehouses

warehouseCost = Sum(WareHouseOpen, data.WareHouseCosts)

# Cost of shops using warehouses

transpCost = Sum([ Sum(varRow, costRow)

for (varRow, costRow) in zip(ShopSupplied, data.SupplyCost)])

obj = warehouseCost + transpCost

model = Model(

# Objective function

Minimise(obj),

# Channel from store opening to store supply matrix

[[var <= store for var in col]

for (col, store) in zip(ShopSupplied.col, WareHouseOpen)],

# Make sure every shop if supplied by one store

[Sum(row) == 1 for row in ShopSupplied.row],

# Make sure that each store does not exceed it's supply capacity

[Sum(col) <= cap

for (col, cap) in zip(ShopSupplied.col, data.Capacity)]

)

return (obj, WareHouseOpen, ShopSupplied, model)

def solve_warehouse_planning(data, param):

(obj, WareHouseOpen, ShopSupplied, model) = model_warehouse_planning(data)

solver = model.load(param['solver'])

solver.setVerbosity(1)

solver.solve()

print obj.get_value()

print "",WareHouseOpen

print ShopSupplied

class WareHouseData:

def __init__(self):

self.NumberOfWarehouses = 5

self.NumberOfShops = 10

self.FixedCost = 30

self.WareHouseCosts = [30, 30, 30, 30, 30]

self.Capacity = [1,4,2,1,3]

self.SupplyCost = supplyCost = [

[ 20, 24, 11, 25, 30 ],

[ 28, 27, 82, 83, 74 ],

[ 74, 97, 71, 96, 70 ],

[ 2, 55, 73, 69, 61 ],

[ 46, 96, 59, 83, 4 ],

[ 42, 22, 29, 67, 59 ],

[ 1, 5, 73, 59, 56 ],

[ 10, 73, 13, 43, 96 ],

[ 93, 35, 63, 85, 46 ],

[ 47, 65, 55, 71, 95 ]

]

solve_warehouse_planning(WareHouseData(), input({'solver':'SCIP'}))

## Quasigroup existence

This is an example showing how to solve the Golomb ruler problem using numberjack.

See prob003 in the CSPLib.. An order m quasigroup is a Latin square of size m. That is, a m by m multiplication table in which each element occurs once in every row and column. For example,

1 2 3 4

4 1 2 3

3 4 1 2

2 3 4 1

is an order 4 quasigroup. A quasigroup can be specified by a set and a binary multiplication opertor, * defined over this set.

Quasigroup existence problems determine the existence or non-existence of quasigroups of a given size with additional properties. Certain existence problems are of sufficient interest that a naming scheme has been invented for them.

QG3.m problems are order m quasigroups for which (a*b)*(b*a) = a.

QG4.m problems are order m quasigroups for which (b*a)*(a*b) = a.

QG5.m problems are order m quasigroups for which ((b*a)*b)*b = a.

QG6.m problems are order m quasigroups for which (a*b)*b = a*(a*b).

QG7.m problems are order m quasigroups for which (b*a)*b = a*(b*a).

For each of these problems, we may additionally demand that the quasigroup is idempotent. That is, a*a=a for every element a.

if t == 3: #(a*b)*(b*a) = a

return (matrix[matrix[a,b],matrix[b,a]] == a)

elif t == 4: #(b*a)*(a*b) = a

return (matrix[matrix[b,a],matrix[a,b]] == a)

elif t == 5: #((b*a)*b)*b = a

return (matrix[matrix[matrix[b,a],b],b] == a)

elif t == 6: #(a*b)*b = a*(a*b)

return (matrix[matrix[a,b],b] == matrix[a,matrix[a,b]])

elif t == 7: #(b*a)*b = a*(b*a)

return (matrix[matrix[b,a],b] == matrix[a,matrix[b,a]])

else:

return []

def quasigroup(T,N):

model = Model()

matrix = Matrix(N,N,N)

print matrix, '\n'

model.add( [AllDiff(row) for row in matrix.row] )

model.add( [AllDiff(col) for col in matrix.col] )

model.add( [matrix[a,a] == a for a in range(N)] ) # idempotency

model.add( [QG(T,a,b) for a in range(N) for b in range(N)] )

print matrix, '\n'

solver = Solver( model )

if solver.solve():

print matrix

else:

print 'no such quasigroup'

## Langford's number

This is an example of solving the Langford's number problem with Numberjack.

Consider two sets of the numbers from 1 to 4. The problem is to arrange the eight numbers in the two sets into a single sequence in which the two 1's appear one number apart, the two 2's appear two numbers apart, the two 3's appear three numbers apart, and the two 4's appear four numbers apart.

The problem generalizes to the L(k,n) problem, which is to arrange k sets of numbers 1 to n, so that each appearance of the number m is m numbers on from the last. For example, the L(3,9) problem is to arrange 3 sets of the numbers 1 to 9 so that the first two 1's and the second two 1's appear one number apart, the first two 2's and the second two 2's appear two numbers apart, etc. See prob024 in the CSPLib.

def model_langford_sequence(M,N):

X = [Variable(1, N*M-(i+2)*(M-1)) for i in range(N)]

model = Model(

AllDiff( [X[i]+((i+2)*j) for j in range(M) for i in range(N)] ),

X[0] > X[1] # Break symmetry

)

return (X,model)

def compute_langford_number(param):

(X,model) = model_langford_sequence(param['M'],param['N'])

solver = model.load(param['solver'])

solver.startNewSearch();

langford_number = 0

while solver.getNextSolution() == SAT:

printLangford(param['M'],X)

langford_number += 1

print 'L('+str(param['M'])+','+str(param['N'])+') =', langford_number

print 'Nodes:', solver.getNodes(), ' Time:', solver.getTime()

def printLangford(M,X):

N = len(X)

sequence = [0]*(M*N)

for i in range(N):

for j in range(M):

sequence[X[i].get_value()+j*(i+2)-1] = (i+1)

print sequence

compute_langford_number(input({'solver':'Mistral', 'N':10, 'M':3}))

## n-Queens

This is the numberjack solution to the famous n-Queens problem.

The n queens puzzle is the problem of putting n chess queens on an nÃ—n chessboard such that none of them is able to capture any other using the standard chess queen's moves. The queens must be placed in such a way that no two queens would be able to attack each other. Thus, a solution requires that no two queens share the same row, column, or diagonal.

def model_queens(N):

queens = [Variable(N) for i in range(N)]

model = Model(

AllDiff( queens ),

AllDiff( [queens[i] + i for i in range(N)] ),

AllDiff( [queens[i] - i for i in range(N)] )

)

return (queens,model)

def solve_queens(param):

(queens,model) = model_queens(param['N'])

solver = model.load(param['solver'])

solver.solve()

print_chessboard(queens)

print 'Nodes:', solver.getNodes(), ' Time:', solver.getTime()

def print_chessboard(queens):

separator = '+---'*len(queens)+'+'

for queen in queens:

print separator

print '| '*queen.get_value()+'| Q |'+' |'*(len(queens)-1-queen.get_value())

print separator

solve_queens(input({'solver':'Mistral', 'N':10}))