Mixed Integer Linear Programming application problem

Question

I have a Mixed-Integer Linear Programming question.

There's this game called Islanders. Pretty graphics, it's about maximizing your score building a small city in a confined island.

enter image description here

Here's a quick rundown of the rules:

1. There are various building types
2. Every building type has a different size
3. Every building type has a different range (explained later)
4. Every building may give a non-zero base score when placed
5. Every building may have give a positive, zero or negative score based on what buildings are within its range when placed
6. Buildings already placed are not affected by newer buildings (their score contribution is not updated)

I thought it would be a fun exercise to formulate these rules into an MILP and see what comes out of it.

I have written the attached Python 3 program, using the PuLP library (link). I have used the default solver (CBC), provided with the installation of PuLP.

There are currently 4 building types encoded, with their base score declared at line 61, the 4x4 cross-score matrix at line 64 and their dimension and range in lines 72 and 75.

In line 82 there is currently stated that there is 1 building of the first type and 3 buildings of the 2nd type to be placed optimally (N=5)

There are several decision variables:

1. px and py are each building coordinates
2. t is the build order (integer)
3. placed_before is a NxN array saying if building i is to be placed before building j (binary)
4. east_of, north_of have similar logic (binary)
5. covers states whether building j is in range of building i (binary)
6. x_overlaps and y_overlaps check whether two buildings are built on top of each other, which is restricted (binary)
7. gives_score is a NxN array saying if building i contributes to the score of building j (is whthin range is is placed before it, binary)

The constraints are as follows:

1. Lines 140-147 state the building order constraints, which form variable t
2. 149-161 decide if buildings are east-west and north-east, used for later
3. 163-183 establish whether buildings overlaps each other, in the east-west (x) direction and the north-south (y) direction
4. 185-219 establish if building i covers building j, filling in the covers variable
5. 222-237 enforces the separation of two buildings, if they overlap in both x and y directions
6. 239-250 finds out if a building i contributes to the score of building j

Currently, the cost function under optimization is gives_score, to maximize the number of buildings which give score to another.

The issue:

• Currently with 1 of the first building and 4 of the second, the optimization does not finish in reasonable time. 380 constraints are generated.
• If instead you go for 3 of the second building type, the problem becomes tractable again, finishing in approx. 10s. 228 constraints are generated.
• Also, if you disable the constraints of lines 222-237, the optimization manages to run, but the buildings overlap.

Any suggestions on the problem formulation or explanations on why this is such a hard problem?

Thank you in advance.

Bonus: Here is an example optimal result for 1 of type-1 building and 3 of type-2: Example optimal result

• Stuart Mitchell: Hi does it take a long time to generate the constraints or in the solve step?
• George ZP: It doesn't take a long time to generate constraints. For the solvable case of 4 buildings, 228 constraints are generated. With the non-terminating case, there are 380 constraints. I use the default solver provided by PuLP, the CBC solver: github.com/coin-or/Cbc
• Stuart Mitchell: Interesting you could try installing an updated version of cbc, or trying to solve it with gurobi.

Answer 1

While I can't answer if your model is okay or if it has some subtle flaw, I think I can confirm that general approach works and the issue is with CBC solver.

I reimplemented the game on my own: https://gist.github.com/afish/56d710cf440bf652fa2ea51f356f7378

CPLEX finds solution for city center and 4 houses in 3 seconds. CBC can't do it (I killed it after couple minutes). I imagine the same situation might be with your code.

I don't know what solvers you can use but seems like CBC is not a good choice here. Just for the record, I also tried MIPCL, SCIP, and GLPK - none of them could solve the problem in a reasonable (~seconds) time.