Informatica

In the most general case, we assume that the constrained set ${\mathcal{R}^{\cup }}$ is a (finite) union of arbitrary closed and bounded subsets (“region-units”) of the problem space $\mathcal{S}$: The bounded property is not restrictive for the problem, since the set $\{{p_{1}},\dots ,{p_{N}}\}$, as a finite set of points, is always bounded, and it can be shown that this implies that all optimal centres ${\mu _{k}}$ for the problem (1) always lie within some bounded subset of $\mathcal{S}$. In addition, the closed property for each region-unit ${\mathcal{R}_{m}}$ ensures that for any $p\in \mathcal{S}$, we can always find the nearest and the furthest points ${\mu _{m}^{\inf }}\in {\mathcal{R}_{m}}$ and ${\mu _{m}^{\sup }}\in {\mathcal{R}_{m}}$, such that

As a less general example (for the multi-locations with polyhedral barriers problem we study in the numerical part of this paper), we define ${\mathcal{R}^{\cup }}$ to be a union of (arbitrary) polygons. We note that any such set ${\mathcal{R}^{\cup }}$ of (possibly, overlapping) polygons can be decomposed into a union of (non-overlapping) triangles by some triangulation method,3 where by “non-overlapping” we mean that all the triangles will “touch” only along the sides. Such decomposition is also possible if some of the region-units ${\mathcal{R}_{m}},\hspace{0.2778em}m=1,\dots ,M$ are segments; then ${\mathcal{R}^{\cup }}$ can be decomposed into a union of “non-overlapping” segments and triangles.

Going to higher dimensions, e.g. for $\mathcal{S}={\mathbb{R}^{3}}$, ${\mathcal{R}^{\cup }}$ could be decomposed into a union of segments, triangles and tetrahedrons (Paulavičius and Žilinskas, 2013), etc., we can continue to higher dimensional Euclidean spaces. This “decomposition into simplices” is a separate research topic and we do not consider this question in our paper. Therefore, in our numerical part for the multi-Weber with barriers problem we assume that ${\mathcal{R}^{\cup }}={\textstyle\bigcup _{m=1}^{M}}{\mathcal{R}_{m}}$ is provided as a union of triangles ${\mathcal{R}_{m}}$, $m=1,\dots ,M$.

For an arbitrary cluster ${\mathcal{C}_{k}}\subseteq \{1,2,\dots ,N\}$ and an arbitrary centre ${\mu _{k}}\in \mathcal{S}$, define the related loss for this cluster-centre pair $({\mathcal{C}_{k}},{\mu _{k}})$: With this definition, the loss in (1) can be written this way: As can be seen from the above equation, the optimal centres given a predefined partition of the clients into clusters can be found independently. Therefore, (1) can be rewritten this way: Formulation (3) can be noticed to lead us to the classical location-allocation (Cooper, 1964) (k-means-type) algorithm: Repeating the procedure above, the loss (3) decreases with every step, and at the termination we are guaranteed to arrive at a local solution of problem (1), i.e. at a solution which cannot be improved by switching cluster memberships or changing cluster centres. For the multi-Weber with polyhedral barriers problem, see, for example, such a location-allocation algorithm by Bischoff et al. (2009).

To find the global solution of problem (1), we could theoretically enumerate all the possible partitions of the index set ${\mathbb{N}_{1}^{N}}=\{1,\dots ,N\}$ into K clusters ${\mathcal{C}_{1}},\dots ,{\mathcal{C}_{K}}$, find the optimal centres for each possible clustering, and pick the best solution.

Such an approach would be brutal and would allow to solve only very small problems, as there are ${S_{2}}(N,K)$ of such partitions into clusters.4 Nevertheless, we have shown in our recent article (Kepalas and Žilinskas, 2024) on a related problem that having a good pruning criteria, problems of size $50/7$ (i.e. $N=50$ points into $K=7$ clusters) become possible to solve using this approach. We will discuss the pruning criteria in Section 3.2. For now, consider the tree $\mathcal{T}(N,K)$ in Fig. 1, which visualizes such an enumeration procedure.

In the tree, the n-th level (distance from the root ${v_{1}}$) of the node corresponds to the n-th element in the set $\{{p_{1}},\dots ,{p_{N}}\}$, and the label λ of the corresponding node represents the cluster (index) to which this element is assigned along the path. As an example, the highlighted path ${v_{1}}{v_{2}}{v_{3}}{v_{4}}{v_{5}}$ in Fig. 1 results in a sequence which corresponds to clustering In general, the enumeration procedure is as follows: It can be shown by simple arguments that this way the nodes ${v_{n}}$ at the level n enumerate all the possible partitions of the set $\{1,2,\dots ,n\}$ into no more than K groups. We notate such partition associated with any vertex ${v_{n}}\in \mathcal{T}(N,K)$ with $\mathcal{P}({v_{n}})$, i.e. $\mathcal{P}({v_{n}})=\{{\mathcal{C}_{1}}({v_{n}}),\dots ,{\mathcal{C}_{K}}({v_{n}})\}$.5 When $n\lt N$, to emphasize that not all N points are yet assigned to some cluster, we sometimes call such a partition “partial” and add the symbol ∂ to our notation, i.e. $\partial \mathcal{P}({v_{n}})=\mathcal{P}({v_{n}})$ for $n\lt N$. Also, we notate the set of all the vertices at the level N with $\mathcal{L}(N,K)$ (this is the leaf-node set of the tree).

The pruning criterion is based on the following simple observation.

Label Consider some ${v_{N}}\in \mathcal{L}(N,K)$ and the corresponding partition Define Now suppose we are at some vertex ${v_{n}}\in \mathcal{T}(N,K)$ at the level $n\lt N$, i.e., we are in the path ${v_{1}}{v_{2}}\dots {v_{n}}$ which corresponds to a partial clustering $\partial {\mathcal{C}_{1}}({v_{n}}),\dots ,\partial {\mathcal{C}_{K}}({v_{n}})$. Consider the following problem, where ${v_{N}}$ is a leaf-node of the tree $\mathcal{T}(N,K)$, i.e. ${v_{N}}\in \mathcal{L}(N,K)$: Looking at the enumeration tree in Fig. 1, in problem (7) we only consider clusterings which correspond to the branch starting at the vertex ${v_{n}}$.

Consider the unique path ${v_{1}}{v_{2}}\dots {v_{n}}$ from ${v_{1}}$ to ${v_{n}}$ and label the corresponding partial partition with $\partial {\mathcal{P}_{n}}=\partial \mathcal{P}({v_{n}})$. For short, we notate problem (7) with $[{\mathbb{P}_{1}^{N}}/K]|[{\mathcal{R}^{\cup }},\partial {\mathcal{P}_{n}}]$, and the optimal solution with ${\textbf{Loss}^{\text{opt}}}([{\mathbb{P}_{1}^{N}}/K]|[{\mathcal{R}^{\cup }},\partial {\mathcal{P}_{n}}])$.

We now state the following theorem. An idea for pruning criteria. Suppose we have an upper bound ${L^{\mathcal{UB}}}$ for the problem $[{\mathbb{P}_{1}^{N}}/K]|[{\mathcal{R}^{\cup }}]$, i.e. suppose we know that Such a bound is given by the optimal solution of any (fixed) particular clustering ${\mathcal{C}_{1}},{\mathcal{C}_{2}},\dots ,{\mathcal{C}_{K}}$. A “good quality” lower bound can be obtained by a location-allocation algorithm as discussed in Section 2.2. Now, if we know that Theorem 2 implies that ${\textbf{Loss}^{\text{opt}}}([{\mathbb{P}_{1}^{N}}/K]|[{\mathcal{R}^{\cup }},\partial {\mathcal{P}_{n}}])\geqslant {L^{\mathcal{UB}}}$ and the branch corresponding to $\partial {\mathcal{P}_{n}}$ can be cut, i.e. there is no need to explore it, as we will not find a better solution there.

Take some arbitrary cluster ${\mathcal{C}_{k}}\subseteq \{1,\dots ,N\}$ and consider the problem i.e. the problem of finding the optimal centre for a given set of points. For a distance metric d induced by a set of barriers, this is sufficiently complicated problem in itself (Klamroth, 2002). The book suggests solving (9) by formulating a mixed integer programming problem. Other suggested methods in the literature use heuristic algorithms (Aneja and Parlar, 1994; Bischoff and Klamroth, 2007).

Instead of trying to solve (9) directly, we will concentrate on finding a lower bound for the problem, i.e. we will define ${\textbf{Loss}^{\mathcal{LB}}}(\{{p_{i}},i\in {\mathcal{C}_{k}}\}|{\mathcal{R}^{\cup }})$ (here, the superscript $\mathcal{LB}$ stands for the “lower bound”) which will satisfy The main thing we seek for this definition is that given the bounds for some cluster ${\mathcal{C}_{k}}$, it would be “easy” to find the bounds for cluster ${\mathcal{C}_{k}}\cup \{j\}$, e.g. if we “increment” the cluster with “client” ${p_{j}}$.

Define the “diameter” $\text{diam}$ of an arbitrary region-unit ${\mathcal{R}_{m}}$ as In addition, suppose that distance metric d satisfies the triangle inequality: given any three points $x,y,z\in \mathcal{S}$, it is always true that In particular, this property holds for the distance metric d induced by a set of barriers, i.e. when $d(x,y)$ is equal to the length of the shortest path between x and y (and the length of an arbitrary path from x to y is defined by the standard length of a curve in the Euclidean geometry).

Now take an arbitrary point ${p_{n}}\in \{{p_{1}},{p_{2}},\dots ,{p_{N}}\}$. Consider the nearest and the furthest point for ${p_{n}}$ in region-unit ${\mathcal{R}_{m}}$, as was defined in (11). By the triangle inequality, we have Thus, if region-unit ${\mathcal{R}_{m}}$ is sufficiently small ($\text{diam}({\mathcal{R}_{m}})\approx 0$), we can expect that ${d_{nm}^{\hspace{0.1667em}\inf }}\approx {d_{nm}^{\hspace{0.1667em}\sup }}$. (If $\text{diam}({\mathcal{R}_{m}})$ is not sufficiently small, we could divide ${\mathcal{R}_{m}}$ into smaller parts and redefine ${\mathcal{R}^{\cup }}$.)

Now consider some cluster ${\mathcal{C}_{k}}$ and suppose that $\text{diam}({\mathcal{R}_{m}})\leqslant \varepsilon $ for any ${\mathcal{R}_{m}}\in \{{\mathcal{R}_{1}},{\mathcal{R}_{2}},\dots ,{\mathcal{R}_{M}}\}$ for some $\varepsilon \gt 0$. By combining (16) and (25), where $|{\mathcal{C}_{k}}|$ is the number of items in cluster ${\mathcal{C}_{k}}$. Therefore, we have the following bounds: Similarly, for ${\textbf{Loss}^{\text{opt}}}(\{{p_{i}},i\in {\mathcal{C}_{k}}\}|{\mathcal{R}^{\cup }})$ we obtain and noting the definition of ${\textbf{Loss}^{\mathcal{LB}}}(\{{p_{i}},i\in {\mathcal{C}_{k}}\}|{\mathcal{R}^{\cup }})$ in (18), inequalities (27) can be rewritten as Summing over all the clusters of a particular partition $\mathcal{P}=\{{\mathcal{C}_{1}},\dots ,{\mathcal{C}_{K}}\}$, we obtain (use the definitions in (6), (21) and inequality above) In particular, the gap between the lower bound and the optimal solution for any partition $\mathcal{P}$ is at most $N\varepsilon $. Thus, we can control the level of this gap by subdividing ${\mathcal{R}^{\cup }}={\textstyle\bigcup _{m=1}^{M}}{\mathcal{R}_{m}}$ into sufficiently small region-units ${\mathcal{R}_{m}}$. These bounds also hold for the globally optimal solution:

We have defined ${\textbf{Loss}^{\mathcal{LB}}}(\{{p_{i}},i\in {\mathcal{C}_{k}}\}|{\mathcal{R}^{\cup }})$ in (18) and have proven in Theorem 4 that Also, we have defined ${\textbf{Loss}^{\mathcal{LB}}}([{\mathbb{P}_{1}^{N}}/K]|[{\mathcal{R}^{\cup }}])$ in (22) which satisfies We have also proven in Theorem 5 the pruning criteria: where $\partial {\mathcal{P}_{n}}=\{\partial {\mathcal{C}_{1}},\dots ,\partial {\mathcal{C}_{K}}\}$ is a (given) partition of index-set ${\mathbb{N}_{1}^{n}}$.

To compute the distances $d({p_{i}},{\mu _{k}})$ between client-points and facility-points, we need to find the so-called visibility-graph $\mathcal{G}(V,E)$ (Viegas and Hansen, 1985). Each point $p\in {\mathcal{R}^{\cup }}$ defines a visibility-polygon as illustrated in Fig. 3. This polygon contains the set of points “seen” from the corresponding point directly through a ray; to reach the points outside this polygon, we must bypass some barriers. It is easy to convince oneself that (see Viegas and Hansen, 1985), in case point $q\in {\mathcal{R}^{\cup }}$ is not seen directly from point p, the shortest path $d(p,q)$ consists of a set of segments with the end-points at the corners of the barrier-regions.6 Therefore, the visibility-polygons for the barrier corner-points are particularly interesting (as in Fig. 3b). Barrier corner-points, together with client-points, define the set of vertices V for the visibility-graph $\mathcal{G}(V,E)$. The (directed) edges E are inserted as follows: if vertex $v\in V$ is “seen” from a barrier-corner-vertex w (i.e. v belongs to the visibility polygon of w), we add the edge $w\mapsto v$ to graph (with edge weight equal to the Euclidean-length from point $p(v)$ to point $p(w)$).

To compute the shortest paths from a given facility-point μ to each client-vertex $v({p_{1}}),\dots ,v({p_{N}})$, we could insert an additional vertex $v(\mu )$ to the graph, find the set of vertices $v\in V$ “seen” from μ, and insert the corresponding edges. Now, the shortest-paths can be found with Dijkstra’s algorithm (Dijkstra, 1959).

Visibility graph for our example problem is shown in Fig. 4.

To apply the ideas from the previous sections, we have to compute the lower-bound distances ${d_{nm}^{\hspace{0.1667em}\inf }},\hspace{0.1667em}n\in {\mathbb{N}_{1}^{N}}$, $m\in {\mathbb{N}_{1}^{M}}$, i.e. the shortest possible distance ${d_{nm}^{\hspace{0.1667em}\inf }}$ from the client-point ${p_{n}}$ to the region-unit (triangle-facility in our case) ${\mathcal{R}_{m}}$.

For this task, we use the visibility-graph $\mathcal{G}(V,E)$ described in the previous subsection.

To find the lower-bound distances ${d_{nm}^{\hspace{0.1667em}\inf }},\hspace{0.1667em}n\in {\mathbb{N}_{1}^{N}}$, we firstly insert an additional vertex $v({\mathcal{R}_{m}})$ into the graph. Now, we have to take care of the edges. This can be done as follows: for each $v\in V$, intersect the visibility polygon $\mathcal{VP}(v)$ of v with triangle ${\mathcal{R}_{m}}$: $\mathcal{I}(v,{\mathcal{R}_{m}})=\mathcal{VP}(v)\cap {\mathcal{R}_{m}}$.7 If the intersection $\mathcal{I}(v,{\mathcal{R}_{m}})$ is empty, do nothing. Otherwise, find the nearest point $\textbf{proj}(v,{\mathcal{R}_{m}})$ from $p(v)$ to $\mathcal{I}(v,{\mathcal{R}_{m}})$.8 Now, add to $\mathcal{G}(V,E)$ a directed-edge $v({\mathcal{R}_{m}})\mapsto v$ with weight $\| p(v)-\textbf{proj}(v,{\mathcal{R}_{m}}){\| _{2}}$. The distances ${d_{nm}^{\hspace{0.1667em}\inf }},\hspace{0.1667em}n\in {\mathbb{N}_{1}^{N}}$ for the triangle-facility ${\mathcal{R}_{m}}$ can now again be found with Dijkstra’s algorithm (Dijkstra, 1959).

The output of this procedure is illustrated in Fig. 5.

To compute an upper bound for the sum of distances from triangle-facility ${\mathcal{R}_{m}}$ to all the clients, we could pick any point $\mu \in {\mathcal{R}_{m}}$ and calculate the distances as outlined in Section 5.1, i.e. using the visibility graph. However, in our case where all the region-units are triangles, for finding an upper-bound faster, we can use the following simple idea.

To find the globally-optimal centre for a set of points, i.e. optimal centree for some cluster ${\mathcal{C}_{k}}$, we use the triangle-partitioning algorithm in the style of Paulavičius and Žilinskas (2013). For a start, we place all the initial triangles (region-units ${\mathcal{R}_{m}}$) in a MinPQ9 (Sedgewick and Wayne, 2011), ordered by the lower-bound value Initially, the MinPQ contains M elements. Now, pick from MinPQ the triangle-facility ${\mathcal{R}_{\min }}$ with the smallest lower-bound value (32). Delete this triangle from MinPQ. Next, partition ${\mathcal{R}_{\min }}$ into two smaller triangles ${\mathcal{R}_{\text{left}}}$ and ${\mathcal{R}_{\text{right}}}$ along the longest edge (see the partitions of the triangle-facility from Fig. 5b in Fig. 6a, Fig. 6b). Both these triangles have an improved (increased) lower-bound estimates because they were obtained by subdividing ${\mathcal{R}_{\min }}$ into two smaller pieces. Since the triangles have been shrunk, the gap between the lower bound and the upper bound for these triangles improves, i.e. we have shown that inequality (28) holds with ${\varepsilon _{\text{left}}}=\text{diam}({\mathcal{R}_{\text{left}}})$ and ${\varepsilon _{\text{right}}}=\text{diam}({\mathcal{R}_{\text{right}}})$. Insert ${\mathcal{R}_{\text{left}}}$ and ${\mathcal{R}_{\text{right}}}$ into MinPQ. Continue deleting/partitioning triangles in MinPQ until the distance between the best known lower-bound and the best known upper-bound is within a predefined accuracy ϵ, i.e. until condition holds for some “small” ϵ (take for example $\epsilon =0.001$). The procedure is illustrated in Fig. 6.

Other steps in the algorithm of Section 4 given a routine to calculate ${d_{nm}^{\hspace{0.1667em}\inf }}$ and a routine to find the optimal centre for a cluster (as described in Section 5.3) are rather straightforward. We just emphasize that during enumeration procedure at any vertex ${v_{n}}\in \mathcal{T}(N,K)$, we use formula (17a) to find in $O(M)$ time.

Some examples of problems solved within ϵ-accuracy are shown in Figs. 9 and 10. These are optimal solutions for largest problem instances we can solve in acceptable time. The results suggest themselves to be compared with the work of Rosing (1992), where the author had studied the unconstrained multi-Weber problem without barriers. This paper was published more than 30 years ago, but the results of their work are comparable to ours in the sense that, just like us (see Fig. 10), they reported global solutions for up to 35/5 and 30/6 problems. Our contribution is that we have managed to solve constrained multi-Weber problems with barriers of the same size, i.e. the ideas presented in this paper allow us to solve much more general class of uncapacitated multi-location problems.