There is an area as shown below.
--A sensor will be installed in the center of this triangle to collect data. -** The purpose is to maximize the total amount of acquired data **. ――Only one sensor can be placed in one triangle. ――If you put one sensor, you can get the amount of data ** 3 **. ――The place where the sensor should be installed and the place where the sensor should not be installed are decided in advance. ――As shown in the table below, if there is another sensor next to it, it will cause interference and the amount of acquired data will be reduced by ** 1 ** respectively. (Total ** 2 ** decrease) ――Therefore, if there are two or more sensors around, it is better not to install them.
Number of sensors around | Amount of acquired data | Amount of interference data | The amount of data that increases as a whole |
---|---|---|---|
0 | 3 | 0 | 3 |
1 | 3 | 2 | 1 |
2 | 3 | 4 | -1 |
3 | 3 | 6 | -3 |
Let's model it as a Mixed Integer optimization Problem (MIP). Please refer to Use Combinatorial Optimization for formulation and ideas.
Maximize td> | $ Amount of acquired data \\ = 3 \ times Whether to install-2 \ times Whether to interfere $ td> tr> |
Variables td> | $ Whether to install (sensor) \ in \ {0, 1 \} \ forall Installation location $ td> tr> |
$ Whether to interfere \ ge 0 \ forall Adjacent $ td> tr> | |
Constraints td> | $ Whether to install _1 + Whether to install _2 --Whether to interfere \ le 1 $ td> tr> |
[Minimum cut problem] in graph theory (https://ja.wikipedia.org/wiki/%E3%82%AB%E3%83%83%E3%83%88_(%E3%82%B0%E3%83] % A9% E3% 83% 95% E7% 90% 86% E8% AB% 96) # .E6.9C.80.E5.B0.8F.E3.82.AB.E3.83.83.E3.83.88) There is something.
Consider the following s-t graph. Finding the s-t minimum cut of this directed graph will give you the best solution to the original sensor installation problem. The original problem is the maximization problem, but by considering "3-acquired data amount" as the cost, it can be regarded as the minimization problem. For details, see the reference site below.
-The above figure shows three examples of (0,0), (0,1), and (1,0). --From the start node "s", draw an edge to each part (triangle). Also, draw an edge from each location to the terminal node "t". --In the upper triangle ((0,0) and (1,0) in the above figure), the side from s means "installed" and the side to t means "not installed". --In the lower triangle ((0,1) in the above figure), the side from s means "not installed" and the side to t means "installed". --The weight of the installed side is 0, and the weight of the non-installed side is 3. (The weight is "3-Amount of acquired data".) --For all adjacent points, draw an edge with a weight of 2 from the lower triangle to the upper triangle. -When fixing with "Installation", delete (cut) the "Installation" side and make the weight of the "Non-installation" side sufficiently large. -When fixing with "Non-installed", delete (cut) the "Non-installed" side and make the weight of the "Installed" side sufficiently large.
--In order to cut s-t, either the installed or non-installed side must be cut at each location. It is considered that the cut one is selected. --In order for the s-t cut to work, if both adjacent locations are installed, it is necessary to cut the side of interference as well, and the amount of interference can be expressed. (It works well by dividing it into upper and lower triangles.)
The minimum cut problem is an easy-to-solve problem that can be solved in polynomial time. In the program described below, it will be solved using Python networkx. By the way, the minimum cut problem dual problem is the maximum flow problem.
Get ready.
python3
import numpy as np, networkx as nx
from pulp import *
def addvar(lowBound=0, var_count=[0], *args, **kwargs):
"""Variable creation utility"""
var_count[0] += 1
return LpVariable('v%d' % var_count[0], lowBound=lowBound, *args, **kwargs)
def calc(a, r):
"""Calculate the amount of acquired data with r as the installation location"""
b = a.copy()
b[b > 1] = 0
for x, y in r:
b[y, x] = 1
s = b.sum() * 3
for y in range(0, b.shape[0], 2):
for x in range(b.shape[1]):
s -= 2 * b[y, x] * b[y+1,x]
if x:
s -= 2 * b[y, x] * b[y+1,x-1]
if y:
s -= 2 * b[y, x] * b[y-1,x]
return s
solve_by_mip is a MIP solution. Returns the installation location.
python3
def solve_by_mip(a):
"""Solve the problem with MIP and return the installation location"""
nm, nn = a.shape
b = a.astype(object)
vv1 = [addvar(cat=LpBinary) for _ in range((b > 1).sum())]
b[b > 1] = vv1
vv2 = []
m = LpProblem(sense=LpMaximize)
for y in range(0, nm, 2):
for x in range(nn):
chk(m, vv2, b[y,x] + b[y+1,x])
if x: chk(m, vv2, b[y,x] + b[y+1,x-1])
if y: chk(m, vv2, b[y,x] + b[y-1,x])
m += 3 * lpSum(vv1) - 2 * lpSum(vv2)
m.solve()
return [(x, y) for x in range(nn) for y in range(nm)
if isinstance(b[y,x], LpVariable) and value(b[y, x]) > 0.5]
def chk(m, vv2, e):
"""If e contains a variable, add a constraint that reduces the objective function by 2 if both are 1."""
if isinstance(e, LpAffineExpression):
v = addvar()
vv2.append(v)
m += e - v <= 1
solve_by_graph is a minimum cut solution. Also return the installation location.
python3
def solve_by_graph(a):
"""Solve the problem with the minimum cut problem and return the installation location"""
nm, nn = a.shape
g = nx.DiGraph()
for y in range(0, nm, 2):
for x in range(nn):
if a[y, x] == 0: # off
g.add_edge('s', (x,y), capacity=7)
elif a[y, x] == 1: # on
g.add_edge((x,y), 't', capacity=7)
else:
g.add_edge('s', (x,y), capacity=0)
g.add_edge((x,y), 't', capacity=3)
if a[y+1, x] == 0: # off
g.add_edge((x,y+1), 't', capacity=7)
elif a[y+1, x] == 1: # on
g.add_edge('s', (x,y+1), capacity=7)
else:
g.add_edge('s', (x,y+1), capacity=3)
g.add_edge((x,y+1), 't', capacity=0)
g.add_edge((x,y+1), (x,y), capacity=2)
if x:
g.add_edge((x-1,y+1), (x,y), capacity=2)
if y:
g.add_edge((x,y-1), (x,y), capacity=2)
r = []
for s in nx.minimum_cut(g, 's', 't')[1]:
b = 's' in s
for t in s:
if isinstance(t, str): continue
x, y = t
if a[y, x] > 1 and b == (y%2 != 0):
r.append((x, y))
return sorted(r)
Let's compare the results by randomly setting fixed points in a size of 40 x 80.
python3
nn, nm = 40, 80 #Horizontal, vertical
np.random.seed(1)
a = np.random.randint(0, 6, (nm, nn)) # 0; fix off, 1: fix on, ow:select
%time rmip = calc(a, solve_by_mip(a))
%time rgrp = calc(a, solve_by_graph(a))
print(rmip == rgrp)
>>>
Wall time: 455 ms
Wall time: 185 ms
True
――You can see that both methods have the same amount of acquired data (rmip == rgrp). --MIP is more than twice as slow. --In general, the dedicated solver is faster than the general-purpose solver. ――When I confirmed it separately, the calculation time of the MIP solver alone was slightly longer than the calculation time of the minimum cut. ――Although the theory is unknown, the amount of acquired data does not change even if MIP is linearly relaxed, and the calculation time is slightly faster.
Reference site
-Problem E. The Year of Code Jam: Referenced issue -Solve the "burn fill problem" using the minimum cut: Referenced solution
that's all
Recommended Posts