Skip to content

Instantly share code, notes, and snippets.

@270ajay
Created December 23, 2020 16:15
Show Gist options
  • Save 270ajay/96c343165693df579729b1c628e453a6 to your computer and use it in GitHub Desktop.
Save 270ajay/96c343165693df579729b1c628e453a6 to your computer and use it in GitHub Desktop.
Advanced constraint programming modeling
from ortools.sat.python import cp_model
''' course: https://www.coursera.org/learn/advanced-modeling#about
These models have different classes of symmetries.
Week6 is about symmetries.
Currently ortools don't have methods for symmetry breaking '''
def crossbowTrapProblem(): # nqueens problem
model = cp_model.CpModel()
sizeOfField = 7
trapVarDict = {}
for i in range(sizeOfField):
for j in range(sizeOfField):
trapVarDict[(i, j)] = model.NewBoolVar("trapIsOn"+str(i)+str(j))
model.Add(sum([trapVarDict[(i,j)] for i in range(sizeOfField) for j in range(sizeOfField)]) == sizeOfField)
for i in range(sizeOfField):
varList = []
for j in range(sizeOfField):
varList.append(trapVarDict[(i,j)])
model.Add(sum(varList) <= 1)
for j in range(sizeOfField):
varList = []
for i in range(sizeOfField):
varList.append(trapVarDict[(i,j)])
model.Add(sum(varList) <= 1)
for k in range(1-sizeOfField, sizeOfField): # forward diagonals
varList = []
for i in range(sizeOfField):
for j in range(sizeOfField):
if k == i-j:
varList.append(trapVarDict[(i,j)])
model.Add(sum(varList) <= 1)
for k in range((2*sizeOfField)-1): # backward diagonals
varList = []
for i in range(sizeOfField):
for j in range(sizeOfField):
if k == i+j:
varList.append(trapVarDict[(i,j)])
model.Add(sum(varList) <= 1)
solver = cp_model.CpSolver()
status = solver.Solve(model)
if status == cp_model.INFEASIBLE:
print("Infeasible model")
if status == cp_model.OPTIMAL:
for i,j in trapVarDict:
if solver.Value(trapVarDict[i,j]) > 0.0:
print("trap is placed at "+str(i)+", "+str(j))
print("----------------------------")
#----------------------------------------------------------------------------------------------------------------------#
def theLampModel():
model = cp_model.CpModel()
numOfRows = 7
numOfColumns = 56
numOfLampsOnEachRow = 24
numOfLampsOnEachColumn = 3
numOfCommonLitLampsInAnyTwoRows = 8
lampsLitVarDict = {}
for row in range(numOfRows):
for col in range(numOfColumns):
lampsLitVarDict[(row, col)] = model.NewBoolVar("lampLitAt"+str(row)+str(col))
for row in range(numOfRows):
varList = []
for col in range(numOfColumns):
varList.append(lampsLitVarDict[(row, col)])
model.Add(sum(varList) == numOfLampsOnEachRow)
for col in range(numOfColumns):
varList = []
for row in range(numOfRows):
varList.append(lampsLitVarDict[(row, col)])
model.Add(sum(varList) == numOfLampsOnEachColumn)
for row1 in range(numOfRows-1):
for row2 in range(row1+1, numOfRows):
varList = []
for col in range(numOfColumns):
isCommonLitLamp = model.NewBoolVar("isCommonLitLampAt"+str(col)+"for"+str(row1)+str(row2))
model.AddBoolAnd([lampsLitVarDict[(row1, col)], lampsLitVarDict[(row2, col)]]).OnlyEnforceIf(isCommonLitLamp)
varList.append(isCommonLitLamp)
model.Add(sum(varList) == numOfCommonLitLampsInAnyTwoRows)
solver = cp_model.CpSolver()
status = solver.Solve(model)
if status == cp_model.INFEASIBLE:
print("Infeasible model")
if status == cp_model.OPTIMAL:
for row in range(numOfRows):
for col in range(numOfColumns):
if solver.Value(lampsLitVarDict[(row, col)]) > 0.0:
print("Lamp lit at row_"+str(row)+", col_"+str(col))
print("----------------------------")
#----------------------------------------------------------------------------------------------------------------------#
def windPrayingModel():
model = cp_model.CpModel()
numOfVertices = 6
numOfPolygon = 3
tributeList = ["LETTUCE", "TURNIP", "CARROT", "BANANA", "PEACH", "GRAPE", "SNAKE", "CHICKEN", "FISHHEAD", "RICEWINE",
"GAOLIANG", "GRAPEWINE", "GOLD", "SILVER", "BRONZE", "GREEN", "RED", "BLUE"]
rankList = [1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4, 4, 5, 5, 5, 6, 6, 6]
rankTributesDict = {1: [0,1,2], 2: [3,4,5], 3:[6,7,8], 4:[9,10,11], 5:[12,13,14], 6:[15,16,17]}
numOfTributes = len(tributeList)
numOfRanks = max(rankList)
tributeVarList = []
polygonVertexIndexDict, counter = {}, 0
for polygon in range(numOfPolygon):
for vertex in range(numOfVertices):
tributeVarList.append(model.NewIntVar(0, numOfTributes-1, "tributeFor"+str(polygon)+str(vertex)))
polygonVertexIndexDict[(polygon, vertex)] = counter
counter += 1
numCoinsVarList = []
for polygon in range(numOfPolygon-1):
for vertex in range(numOfVertices):
numCoinsVarList.append(model.NewIntVar(0, numOfRanks, "numCoinsfor"+str(polygon)+str(vertex)))
model.AddAllDifferent(tributeVarList)
rankVarList = []
for polygon in range(numOfPolygon):
for vertex in range(numOfVertices):
rankVarList.append(model.NewIntVar(1, numOfRanks, "rankFor"+str(polygon)+str(vertex)))
index = polygonVertexIndexDict[(polygon, vertex)]
model.AddElement(tributeVarList[index], rankList, rankVarList[index])
positionVarList = []
for tribute in range(numOfTributes):
positionVarList.append(model.NewIntVar(0, (numOfPolygon*numOfVertices)-1, "positionForTribute"+str(tribute)))
model.AddInverse(tributeVarList, positionVarList)
for polygon in range(numOfPolygon-1):
for vertex in range(numOfVertices):
diff = model.NewIntVar(-numOfRanks, numOfRanks, "diffInRanksFor"+str(polygon)+str(polygon+1)+str(vertex))
index = polygonVertexIndexDict[(polygon, vertex)]
index2 = polygonVertexIndexDict[(polygon+1, vertex)]
model.Add(diff == rankVarList[index] - rankVarList[index2])
model.AddAbsEquality(numCoinsVarList[index], diff)
for vertex in range(numOfVertices-1): # removing variable symmetry
index = polygonVertexIndexDict[(0, vertex)]
index2 = polygonVertexIndexDict[(0, vertex+1)]
model.Add(tributeVarList[index] < tributeVarList[index2])
index = polygonVertexIndexDict[(0, 0)]
index2 = polygonVertexIndexDict[(numOfPolygon-1, 0)]
model.Add(tributeVarList[index] < tributeVarList[index2]) # removing variable symmetry
for rank, tributeList in rankTributesDict.items(): # removing variable symmetry/value symmetry
for tribute in range(len(tributeList)-1):
model.Add(positionVarList[tribute] < positionVarList[tribute+1])
model.Maximize(sum(numCoinsVarList))
solver = cp_model.CpSolver()
status = solver.Solve(model)
if status == cp_model.INFEASIBLE:
print("Infeasible model")
if status == cp_model.OPTIMAL:
print('Max number of coins: %i' % solver.ObjectiveValue())
print("----------------------------")
#----------------------------------------------------------------------------------------------------------------------#
def generalAllocationProblem ():
model = cp_model.CpModel()
numOfGenerals = 14
numOfRoads = 14
coopList = [[0, 15, 3, 13, 4, 1, 5, 4, 4, 1, 1, 11, 10, 1],
[15, 0, 3, 13, 4, 11, 10, 9, 11, 7, 1, 5, 4, 1],
[3, 3, 0, 13, 4, 11, 10, 9, 10, 7, 1, 5, 4, 1],
[13, 13, 13, 0, 4, 11, 10, 9, 9, 7, 1, 5, 4, 1],
[4, 4, 4, 4, 0, 1, 5, 4, 10, 1, 1, 11, 10, 1],
[1, 11, 11, 11, 1, 0, 10, 9, 8, 7, 1, 5, 4, 1],
[5, 10, 10, 10, 5, 10, 0, 9, 8, 7, 1, 5, 4, 1],
[4, 9, 9, 9, 4, 9, 9, 0, 8, 7, 1, 5, 4, 1],
[4, 11, 10, 9, 10, 7, 8, 8, 0, 7, 1, 5, 4, 1],
[1, 7, 7, 7, 1, 7, 7, 7, 7, 0, 1, 5, 4, 1],
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 5, 4, 1],
[11, 5, 5, 5, 11, 5, 5, 5, 5, 5, 5, 0, 4, 4],
[10, 4, 4, 4, 10, 4, 4, 4, 4, 4, 4, 4, 0, 3],
[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 4, 3, 0]]
coopListFlattened = sum(coopList, [])
maxCoop = max(coopListFlattened)
roadVarList = []
for road in range(numOfRoads):
roadVarList.append(model.NewIntVar(0, numOfGenerals-1, "generalForRoad"+str(road)))
model.AddAllDifferent(roadVarList)
coopVarList = []
for road in range(numOfRoads-1):
indexVar = model.NewIntVar(0, (numOfGenerals*numOfGenerals)-1, "indexFor"+str(road)+str(road+1))
coopVarList.append(model.NewIntVar(0, maxCoop, "coopFor"+str(road)+str(road+1)))
model.Add(indexVar == (roadVarList[road] * numOfGenerals) + roadVarList[road+1])
model.AddElement(indexVar, coopListFlattened, coopVarList[road])
model.Add(roadVarList[0] < roadVarList[numOfRoads-1]) # variable symmetry breaking
for road1 in range(1, numOfRoads-2): # dominance breaking constraint
for road2 in range(road1+1, numOfRoads-1):
coop1Var = model.NewIntVar(0, maxCoop, "coop1For"+str(road1-1)+str(road2))
coop2Var = model.NewIntVar(0, maxCoop, "coop2For"+str(road1)+str(road2+1))
index1Var = model.NewIntVar(0, (numOfGenerals*numOfGenerals)-1, "index1For"+str(road1-1)+str(road2))
index2Var = model.NewIntVar(0, (numOfGenerals*numOfGenerals)-1, "index2For" +str(road1) + str(road2+1))
model.Add(index1Var == (roadVarList[road1-1] * numOfGenerals) + roadVarList[road2])
model.Add(index2Var == (roadVarList[road1] * numOfGenerals) + roadVarList[road2+1])
model.AddElement(index1Var, coopListFlattened, coop1Var)
model.AddElement(index2Var, coopListFlattened, coop2Var)
model.Add(coopVarList[road1-1] + coopVarList[road2] >= coop1Var + coop2Var)
# solList = [10, 12, 4, 8, 9, 7, 6, 5, 2, 3, 1, 0, 11, 13]
# for i,var in enumerate(roadVarList):
# model.Add(var == solList[i])
model.Maximize(sum(coopVarList))
solver = cp_model.CpSolver()
status = solver.Solve(model)
if status == cp_model.INFEASIBLE:
print("Infeasible model")
if status == cp_model.OPTIMAL:
print('Max coop: %i' % solver.ObjectiveValue())
print("----------------------------")
#----------------------------------------------------------------------------------------------------------------------#
if __name__ == '__main__':
crossbowTrapProblem()
theLampModel()
windPrayingModel()
generalAllocationProblem()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment