Warehouse Location

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

Description: 

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.

Code: 

from Numberjack import *

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.

Description: 

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.

Code: 

def QG(t,a,b):
    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.

Description: 

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.

Code: 

from Numberjack import *

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.

Description: 

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.

Code: 

from Numberjack import *

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}))

Syndicate content