# Trying to solve a puzzle

Just recently I bumped on this: The PuzzlOR and I want to play a bit with that and solve it using JuMP.

My approach is not yet certain, but something like this:

1. load data into DICT(string, int) * actually just 0 or 1
2. total houses * 0.7 will be maximum for the objective function

Now, can I use indicator constraints to set value for my variables?

If cell A2 will hold a value then that value should be equal to sum from all adjusting cells in DICT (including that cell), and then set adjusting cells variables to be = 0.

Is it possible to set variables based with indicator constraints?

`amount_covered_houses >= 0.7*total_houses` should be a constraint, not an objective, and your objective should be to minimize the amount of antennas

I think I solved it with JuMP, using only Binary variables
I used matrices rather than Dicts to store data and variables

to know if a cell is covered, your constraint can be
`iscovered <= sum(isantenna.(neighbors))`
since these are binary variables and the solver will try to maximize `iscovered`, it’s equivalent to
`iscovered == |(isantenna.(neighbors)...)`

3 Likes

First, thanks for your answer, but what is this syntax? Is isantenna also a variable?

this is kinda pseudo code, with `isantenna` represented as a function
in my solution `isantenna` is a Matrix, so this syntax would cleary not work

I suppose this is the way to declare those variables and constraints:

``````     @variable(prmx,
isCovered[1:_ROWS, 1:_COLUMNS],
Bin)

@variable(prmx,
isAntenna[1:_ROWS, 1:_COLUMNS],
Bin)

@objective(
prmx,
Min,
sum(isCovered[i, j] * isHouse[i, j]
for i in 1:_ROWS, j in 1:_COLUMNS
)
)

@constraint(
prmx,
sum(isCovered[i, j] * isHouse[i, j]
for i in 1:_ROWS, j in 1:_COLUMNS
)
>= _TOTAL * 0.7
)

``````

It looks good

Your last set of constraints is useless I think

And the objective should be to minimize the amount of Antennas, you already have the 70% rule as a constraint

Yes, they say “at least 70%” so it has to bi Min with constraint sum(…) >= _TOTAL * 0.7.

Now I can;t quite define constraints for isAntenna → isCovered. Per my logic it should be if there is Antenna it should cover adjusting cells, but still dont know how to write constraints for that.

I think it worked with

``````    @constraint(
prmx,
[i in 1:_ROWS, j in 1:_COLUMNS],
isCovered[i, j]
<=
isAntenna[i, j]
+ sum(
isAntenna[r, c]
for r in i - 1:i + 1, c in j - 1:j + 1
if r != i && r >= 1 && r <= _ROWS && c != j && c >= 1 && c <= _COLUMNS
)
)
``````

you can simplify a bit:

``````    @constraint(
prmx,
[i in 1:_ROWS, j in 1:_COLUMNS],
isCovered[i, j]
<=
sum(
isAntenna[r, c]
for r in i - 1:i + 1, c in j - 1:j + 1
if 1 <= r <= _ROWS && 1 <= c <= _COLUMNS
)
)
``````

and the objective should simply be the total amount of antennas:

``````@objective(prmx, Min, sum(isAntenna))
``````

then it should be good I think

for @objective I agree and understand, but can’t visualize that isCovered and isAntenna relationship, here it says that isCovered <= sum(isAntenna), hm…, should it be other way around like isCovered >= sum(isAntenna)?

if you have `isCovered >= sum(isAntenna)`, your solver will say,
“ok, isAntenna = 0 everywhere and isCovered = 1 everywhere”

which is of course not want we want (although coverage without antennas would be great )

there is proper maths explanations to transform `and` and `or` into linear constraints but I couldn’t explain the proper maths theory

I found this post that explains how to transform boolean logic in linear equations

the proper complete way to do a `or` is (if b is a list):

``````a = |(b...)
``````

is equivalent to

``````a <= sum(b)
a .>= b
``````

in this case the `a .>= b` isn’t strictly needed (I think) since the solver will try to minimize b and maximize a anyway

``````    @constraint(
prmx,
[i in 2:_ROWS-1, j in 2:_COLUMNS-1],
sum(
isCovered[r, c]
for r in i - 1:i + 1, c in j - 1:j + 1
)  - 9 * isAntenna[i, j]  >= 0
)
``````
``````isCovered[i, j] <=
sum(
isAntenna[r, c]
for r in i-1:i+1, c in j-1:j+1
if 1 <= r <= _ROWS && 1 <= c <= _COLUMNS
)
``````

corresponds to

``````isCovered[i, j] ==
|([
isAntenna[r, c]
for r in i-1:i+1, c in j-1:j+1
if 1 <= r <= _ROWS && 1 <= c <= _COLUMNS
]...)
``````

and for completeness you can add:

``````isCovered[i, j] .>=
[
isAntenna[r, c]
for r in i-1:i+1, c in j-1:j+1
if 1 <= r <= _ROWS && 1 <= c <= _COLUMNS
]
``````

this works fine for me and I find the right answer (6 antennas for 29 houses covered)
I don’t understand what you’re trying to do in your last post

1 Like

Here’s a more compact way to write it:

``````using JuMP, Cbc
function main(; houses::Vector{Tuple{Int,Int}}, ratio::Float64)
M = maximum(h for h in houses)
N = maximum(h for h in houses)
model = Model(Cbc.Optimizer)
@variable(model, cell_tower[1:M, 1:N], Bin)
@variable(model, is_covered[houses], Bin)
@objective(model, Min, sum(cell_tower))
stencil(i, j) = cell_tower[max(i-1,1):min(i+1,M), max(j-1,1):min(j+1,N)]
@constraint(model, [h in houses], is_covered[h] <= sum(stencil(h, h)))
@constraint(model, sum(is_covered) >= ratio * length(houses))
optimize!(model)
return value.(cell_tower)
end
houses = [(1, 9), (3, 2), (4, 4), (5, 6), (8, 9), (9, 1)]
main(houses, ratio = 0.7)
``````
2 Likes

The most important, what I failed to understand, is that isCovered is in the function of houses. I was thinking about the grid that one antenna can cover, but I should think about which houses antenna(s) can cover instead.

Here is the model:

``````using JuMP, Cbc, XLSX

M = size(mat, 1)
N = size(mat, 2)
ratio = 0.7

houses = []
for i = 1 : M
for j = 1 : N
if mat[i,j] == 1
push!(houses, (i,j))
end
end
end

model = Model(Cbc.Optimizer)

@variable(model,
isCovered[houses],
Bin)

@variable(model,
isAntenna[1:M, 1:N],
Bin)

@objective(
model,
Min,
sum(isAntenna)
)

@constraint(
model,
sum(isCovered) >= length(houses) * ratio
)

@constraint(
model,
[h in houses],
isCovered[h] <=
sum(
isAntenna[r, c]
for r in h-1:h+1, c in h-1:h+1
if 1 <= r <= M && 1 <= c <= N
)
)

optimize!(model)

res = primal_status(model)

if res == MOI.FEASIBLE_POINT
val = value.(isAntenna)
open("out/result.txt","w") do io
for i in 1:M
for j in 1:N
print(io, floor(Int, val[i,j]))
print(io, '\t')
end
print(io, "\n")
end
end
end

``````

Got it finally! 6 antennas → 29 houses

Great!

To see a Matrix approach, here was my take on the problem:

``````using JuMP, Cbc
m = Model(Cbc.Optimizer)
set_silent(m)

houses = [
0 0 1 0 1 1 0 1 0 0
1 0 0 0 1 0 0 0 0 0
0 0 0 1 0 1 1 0 1 1
0 0 0 0 1 0 1 0 1 0
0 0 1 0 0 1 1 1 1 0
1 0 0 1 1 0 1 0 0 0
1 1 1 1 1 0 0 0 0 1
0 0 0 0 1 0 1 0 0 0
1 0 1 0 0 0 1 1 0 1
1 0 0 0 1 1 0 0 1 1
]

@variables m begin
isantenna[1:12,1:12], Bin
iscovered[1:10,1:10], Bin
end

@constraints m begin
isantenna[  1,  :] .== 0
isantenna[end,  :] .== 0
isantenna[  :,  1] .== 0
isantenna[  :,end] .== 0

iscovered .<= sum(isantenna[i:end-3+i, j:end-3+j] for i in 1:3, j in 1:3)

sum(iscovered .* houses) >= 0.7*sum(houses)
end

@objective m Min sum(isantenna)

optimize!(m)

result = Int.(value.(isantenna[2:end-1, 2:end-1]))
println(sum(result))
display(result)
``````
1 Like

Great! But I am still a conventional programmer though and will need some time to understand and start using those operators like `.` Anyway I am already on another page: The PuzzlOR did you solve that too?

hmmm, that’s not linear (and possibly not convex)

I might try it when I have some time

1 Like

with 1,000,000 possibilities only, it’s totaly reasonable to brute-force the problem
without optimizing my code, I find the solution in 7 seconds

2 Likes