In this notebook, we'll tackle the problem of trying to protect a parade route by placing a limited number of guards.

The parade route is discritized into $m$ points. There are $n$ possible guard locations with associated decision variable $x \in
\lbrace 0,1\rbrace^n$, where $x_i = 1$ if and only if a guard is placed at location $i$. Associated with guard location $i$ is a *coverage vector* $a_i \in \mathbf{R}^m$, which describes how well a guard placed at location $i$ would 'cover' each point in the parade route. We will assume that guard coverage is additive, so that the vector describing the total coverage of every edge is given by $Ax$, where $A \in \mathbf{R}^{m \times n}$ has $a_i$ as its $i$th column.

The parade route is only as secure as its least well-covered point. Our goal is to place $k$ guards to maximize the minimum coverage over the points in the route.

We can formulate this as the optimization problem

\[ \begin{array}{ll} \mbox{maximize} & t \ \mbox{subject to} & t \leq Ax\ & x \in \lbrace 0,1 \rbrace^n \ & \mathbf{1}^Tx = k. \end{array} \]

This problem is nonconvex and, in general, NP-hard due to the Boolean decision variable.

We can try to approach the problem with convex optimization by first forming the convex relaxation

\[ \begin{array}{ll} \mbox{maximize} & t \ \mbox{subject to} & t \leq Ax\ & 0 \leq x \leq 1 \ & \mathbf{1}^Tx = k, \end{array} \] by constraining $x \in [0,1]^n$.

In general, the solution to this problem, $x^\star$, will have fractional values. As we want a Boolean allocation, we can use an iterated weighted $\ell_1$ heuristic to try to recover a Boolean solution.

To try and recover a Boolean solution, we will solve a sequence of convex problems where we add a linear term $-w^Tx$ to the objective, picking the weight vector $w \in \mathbf{R}^n_+$ at each iteration to try and induce a sparse solution vector $x^\star$. The details can be found in the Stanford EE364B lecture notes.

The algorithm consists of initializing $w = 0$ and repeating the two steps

- \[ \begin{array}{ll} \mbox{maximize} & t - w^Tx\ \mbox{subject to} & t \leq Ax\ & 0 \leq x \leq 1 \ & \mathbf{1}^Tx = k, \end{array} \]
- Let $w_i = \alpha/(\tau + x_i) \forall i$

until we reach a Boolean solution. Here, $\alpha$ and $\tau$ are adjusted to promote a sparse solution. Typical choices would be $\alpha = 1$ and $\tau = 10^{-4}$.

Intuitively, the weight vector $w$ is incentivizing elements of $x$ which were close to zero in the last iteration towards zero in the next iteration.

We create a parade route in the unit square $[0,1] \times [0,1]$ by generating points along a connected sequence of line segments. Possible guard locations are a set of randomly placed points in the unit square. Guards' coverage of points in the parade route is a function of distance between the two points. We also add buildings in the unit square to obstruct the guards' view. The guard has no effect on a point if his line of sight is obstructed. We generate the $A$ matrix for this problem instance below.

In [1]:

```
import numpy as np
import matplotlib.pyplot as plt
import cvxpy as cvx
def form_path(points,n):
x, y = [], []
pold = points[0]
for p in points[1:]:
x += list(np.linspace(pold[0],p[0],n))
y += list(np.linspace(pold[1],p[1],n))
pold = p
path = np.array([x,y]).T
return path
def form_grid(k):
xs = list(np.linspace(0,1,k))
ys = list(np.linspace(0,1,k))
locations = []
for x in xs:
for y in ys:
locations.append(np.array((x,y)))
return np.array(locations).T
def guard_sets(k,num,noise):
guard_set = []
grid = form_grid(k)
for i in range(num):
pert = noise*np.random.randn(*grid.shape)
guard_set.append( grid+pert )
return np.hstack(guard_set)
def inRect(p,rect):
x,y,w,h = rect
return x <= p[0] and p[0] <= x + w and y <= p[1] and p[1] <= y + h
def remove_guards(guards,buildings):
'''Remove guards inside buildings and outside unit square.'''
outside = []
for i, guard in enumerate(guards.T):
inside = False
for build in buildings:
if inRect(guard,build):
inside = True
break
if not inRect(guard,(0,0,1,1)):
inside = True
break
if not inside:
outside.append(i)
return guards[:,outside]
def intersect(p1,p2,xmin,xmax,ymin,ymax):
'''determine if a rectangle given by xy limits blocks the line of sight between p1 and p2'''
block = False
# if either point inside block
for p in [p1,p1]:
if xmin <= p[0] and p[0] <= xmax and ymin <= p[1] and p[1] <= ymax:
return True
# if the two points are equal at this stage, then they are outside the block
if p1[0] == p2[0] and p1[1] == p2[1]:
return False
if p2[0] != p1[0]:
for x in [xmin,xmax]:
alpha = (x-p1[0])/(p2[0] - p1[0])
y = p1[1] + alpha*(p2[1] - p1[1])
if 0 <= alpha and alpha <= 1 and ymin <= y and y <= ymax:
return True
if p2[1] != p1[1]:
for y in [ymin,ymax]:
alpha = (y-p1[1])/(p2[1] - p1[1])
x = p1[0] + alpha*(p2[0] - p1[0])
if 0 <= alpha and alpha <= 1 and xmin <= x and x <= xmax:
return True
return False
def p_evade(x,y,r=.5,minval=.1):
d = np.linalg.norm(x-y)
if d > r:
return 1
return (1-minval)*d/r + minval
def get_guard_effects(path, guards, buildings, evade_func):
guard_effects = []
for guard in guards.T:
guard_effect = []
for p in path:
prob = 1
if not np.any([intersect(p,guard,x,x+w,y,y+h) for x,y,w,h in buildings]):
prob = evade_func(p,guard)
guard_effect.append(prob)
guard_effects.append(guard_effect)
return np.array(guard_effects).T
locations = []
for x in xs:
for y in ys:
point = np.array((x,y))
detect_p = []
for r in path:
detect_p.append(p_evade(point,r,r=.5,m=0))
locations.append((point,np.array(detect_p)))
```

In [3]:

```
np.random.seed(0)
buildings = [(.1,.1,.4,.1),
(.6,.1,.1,.4),
(.1,.3,.4,.1),
(.1,.5,.4,.1),
(.4,.7,.4,.1),
(.8,.1,.1,.3),
(.8,.5,.2,.1),
(.2,.7,.1,.3),
(.0,.7,.1,.1),
(.6,.9,.1,.1),
(.9,.7,.1,.2)]
n = 10
points = [(.05,0),(.05,.25),(.55,.25),(.55,.6),(.75,.6),(.75,.05),(.95,.05), (.95,.45),(.75,.45), (.75,.65),(.85,.65),
(.85,.85),(.35,.85),(.35,.65),(.15,.65),(.15,1)]
path = form_path(points,n)
g = guard_sets(12,4,.02)
g = remove_guards(g,buildings)
guard_effects = get_guard_effects(path, g, buildings, p_evade)
A = 1 - np.log(guard_effects)
fig = plt.figure(figsize=(10,10))
ax = plt.subplot(111,aspect='equal')
for x,y,w,h in buildings:
rect = plt.Rectangle((x,y),w,h,fc='y',alpha=.3)
ax.add_patch(rect)
ax.plot(path[:,0],path[:,1],'o')
ax.plot(g[0,:],g[1,:],'ro',alpha=.3)
```

Out[3]:

In [14]:

```
num_guards = 12
tau = 1e-2
m,n = A.shape
w = np.zeros(n)
for i in range(3):
x = cvx.Variable(n)
t = cvx.Variable(1)
objective = cvx.Maximize(t - x.T*w)
constr = [0 <=x, x <= 1, t <= A*x,cvx.sum_entries(x) == num_guards]
cvx.Problem(objective, constr).solve(verbose=False)
x = np.array(x.value).flatten()
w = 2/(tau+np.abs(x))
fig = plt.figure(figsize=(5,5))
ax = fig.add_subplot(111)
ax.plot(x,'o')
xsol = x
```

In [10]:

```
fig = plt.figure(figsize=(10,10))
ax = plt.subplot(111,aspect='equal')
for x,y,w,h in buildings:
rect = plt.Rectangle((x,y),w,h,fc='y',alpha=.3)
ax.add_patch(rect)
ax.plot(path[:,0],path[:,1],'o')
ax.plot(g[0,:],g[1,:],'ro',alpha=.3)
ax.plot(g[0,xsol > .5],g[1,xsol > .5],'go',markersize=20,alpha=.5)
```

Out[10]:

In [ ]:

```
```