On the Complexity of Minimum Riesz s-Energy Subset Selection in Euclidean and Ultrametric Spaces
May 06, 2026
Abstract
We study the computational complexity of exact cardinality-constrained minimum Riesz \(s\)-energy subset selection in finite metric spaces: given \(n\) points, select \(k<n\) points of minimum Riesz \(s\)-energy. The objective sums inverse-power pair interactions and therefore promotes well-separated subsets; as \(s\) becomes large, it increasingly approaches a bottleneck criterion governed by the closest selected pair, linking it to minimum pairwise distance (MPD). Building on the general-metric NP-hardness result of Pereverdieva et al. (2025), we prove that NP-hardness persists for point sets in the Euclidean plane when \(s\) is part of the input. In contrast, finite ultrametric spaces form an exact tractable regime: on rooted binary ultrametric trees with \(n\) leaves, an optimal size-\(k\) subset can be computed by dynamic programming in \(O(nk^2)\) time. We also discuss the ordered one-dimensional Euclidean case, where the classical MPD objective admits simple dynamic programming, but the additive Riesz energy does not appear to allow the same state compression. Finally, we explain why one natural route to fixed-\(s\) Euclidean hardness does not close: Fowler-style 3SAT gadgets, together with zeta-function bounds for far-field interactions, show why this approach still requires an exponent depending on \(k\). Together, these results provide a compact complexity landscape for a natural diversity or dispersion objective, distinguishing Euclidean hardness, ultrametric tractability, and the ordered one-dimensional case.
Keywords: Riesz \(s\)-Energy; Subset Selection; Computational Complexity; Ultrametrics; Minimum Pairwise Distance
1 Introduction↩︎
Riesz energies are a classical family of pairwise interaction objectives. For exponent \(s>0\), nearby points interact strongly through the term \(d(x,y)^{-s}\), and minimizing the total interaction encourages well-separated configurations. This viewpoint is central in potential theory and in the study of energy-minimizing point distributions on spheres and manifolds, where it is shown to lead to uniformly spaced configurations; see, for example, the survey of Brauchart and Grabner [1] and the reference work [2]. Beyond this classical setting, minimum-subset selection for discrete Riesz \(s\)-energy has recently received broader attention in applied multiobjective optimization and diversity maximization, where it is used to extract well-distributed finite Pareto-front approximations and to compare diversity indicators, and it has also appeared in closely related geometric optimization problems in theoretical computer science [3]–[5].
A useful perspective is that, on every fixed finite candidate set, minimum Riesz \(s\)-energy interpolates between pairwise repulsion and minimum pairwise distance dispersion. Indeed, as \(s\to\infty\), the contribution of the closest selected pair dominates the sum, so minimizing Riesz \(s\)-energy tends to maximizing the minimum pairwise distance. In this sense, the MPD optimization problem appears as the large-\(s\) limit of minimum Riesz \(s\)-energy subset selection; see Chapter 3 of Borodachov, Hardin, and Saff [2] and also Atta [3].
Max–min diversity has been studied for decades, not only in abstract metric spaces but also in geometric settings, where hardness arguments are closely related to geometric independent set and packing constructions [6]–[8]. Our Euclidean-plane hardness proof belongs to this line of work. At the same time, the Riesz objective is structurally subtler: unlike minimum pairwise distance (MPD), it depends on the full collection of pairwise interactions rather than only on the closest selected pair.
In this paper we study the cardinality-constrained subset-selection problem for the discrete Riesz \(s\)-energy. If \(S\) is a \(k\)-subset of a metric space, then the objective aggregates all pairwise interactions \(d(x,y)^{-s}\) over distinct selected points. Thus it differs from the maximum-sum dispersion objective, which directly rewards large distances, and from the classical MPD objective, which is determined solely by the closest selected pair. At the same time, for large \(s\) the Riesz objective increasingly concentrates on the smallest selected distance. In that sense it forms a natural bridge between pairwise-sum repulsion and minimum pairwise distance (MPD).
Our main theme is a sharp contrast between hard general metric cases and a tractable hierarchical regime. A general-metric NP-hardness result was already established by Pereverdieva et al. [5]. Building on that hard side of the landscape, we show that NP-hardness persists in the Euclidean plane when \(s\) is part of the input. By contrast, ultrametric spaces admit an exact dynamic programming algorithm. This tractable case is still nontrivial: one must distribute a fixed budget across competing branches of a hierarchy under a global objective, but the ultrametric structure collapses every cross-subtree interaction to a simple cardinality-dependent term.
The paper is deliberately selective. Rather than attempting a full survey of all geometric and parameterized variants, we isolate one coherent picture: hardness in general metrics, exact tractability on ultrametrics, and hardness in Euclidean spaces. Throughout the study, we relate Riesz energy to minimum pairwise distance (MPD), treating it as a limiting case.. This makes the theorem flow short and leaves the unresolved directions visible rather than buried.
1.0.0.1 Contributions.
The paper has five principal contributions, one of which is contextual.
For general finite metric spaces, we reproduce and refine the NP-hardness result of Pereverdieva et al. [5]. For completeness and in order to keep the paper self-contained, we include a short hardness proof in our notation and make the reduction argument explicit.
We prove Euclidean NP-hardness in the plane when \(s\) is part of the input.
We identify finite ultrametric spaces as a tractable island by deriving an exact dynamic programming algorithm running in \(O(nk^2)\) time on rooted binary ultrametric trees.
We record the one-dimensional comparison with minimum pairwise distance (MPD): on the line, the classical MPD objective admits a Bellman-style dynamic program, whereas no analogous state compression is available for minimum Riesz \(s\)-energy. This clarifies why the one-dimensional Euclidean case remains structurally delicate.
For the minimum pairwise distance (MPD) comparison, we use the known finite-instance fact that for sufficiently large \(s\), every Riesz-energy minimizer is MPD, and we highlight the implications of this limiting case for hardness and tractability.
These results leave a short list of structurally meaningful open cases, most notably fixed-exponent Euclidean hardness.
2 Preliminaries↩︎
Definition 1 (Minimum Riesz \(s\)-energy subset selection). Let \((X,d)\) be a finite metric space and let \(s>0\) be fixed. For a subset \(S\subseteq X\) with \(|S|\ge 2\), define the unordered-pair Riesz \(s\)-energy by \[E_{s}\!\left(S\right)=\sum_{\{u,v\}\subseteq S} d(u,v)^{-s}.\] We also set \(E_s(S)=0\) when \(|S|\le 1\). The optimization problem asks for a \(k\)-subset minimizing \(E_s(S)\).
We will use the associated decision problem.
RSSP-Decision: Given \((X,d)\), a target size \(k\), an exponent \(s>0\), and a threshold \(T\), decide whether there exists a subset \(S\subseteq X\) with \(|S|=k\) and \(E_s(S)\le T\).
We also use the classical minimum pairwise distance (MPD) objective \[\operatorname{MPD}\!\left(S\right):=\min_{\{u,v\}\subseteq S} d(u,v),\] with the convention that \(\operatorname{MPD}\!\left(S\right)=+\infty\) for \(|S|\le 1\).
Definition 2 (Ultrametric). A metric \(d\) on a set \(X\) is called an ultrametric* if \[d(x,z)\le \max\{d(x,y),d(y,z)\} \qquad\text{for all }x,y,z\in X.\] Equivalently, for every triple of points, the two largest pairwise distances are equal.*
Every finite ultrametric admits a rooted-tree representation: the points are the leaves, internal nodes carry nondecreasing heights, and the distance between two leaves is twice the height of their least common ancestor. Ultrametrics arise naturally from rooted phylogenetic trees, where distances may represent species separation in evolutionary time, and from the cophenetic distances of hierarchical clustering dendrograms, including hierarchical analyses of cultural vectors [9], [10].
3 General metric spaces via \(k\)-Clique↩︎
We start with the purely metric setting. The corresponding NP-hardness statement already appears in the comparative analysis of Pereverdieva et al. [5]; we include a short self-contained proof because it serves as the clean baseline for the rest of the paper.
Definition 3 (\(k\)-Clique). Given a graph \(G=(V,E)\) and an integer \(k\), decide whether there exists a subset \(V'\subseteq V\) of size \(k\) such that every pair of distinct vertices in \(V'\) is joined by an edge.
Lemma 1. Let \(G=(V,E)\). Define \(d:V\times V\to\{0,1,2\}\) by \[d(u,v)=\begin{cases} 0,&u=v,\\ 2,&u\ne v\text{ and }\{u,v\}\in E,\\ 1,&u\ne v\text{ and }\{u,v\}\notin E. \end{cases}\] Then \((V,d)\) is a metric space.
Proof. Non-negativity, symmetry, and \(d(u,v)=0\Leftrightarrow u=v\) are immediate. For the triangle inequality, the only nontrivial case is \(d(u,v)=2\). If \(u\neq w\neq v\), then each of \(d(u,w)\) and \(d(w,v)\) belongs to \(\{1,2\}\), so \(2\le d(u,w)+d(w,v)\). All other cases are trivial. ◻
Theorem 1. For every fixed \(s>0\), RSSP-Decision is NP-hard for general finite metric spaces.
Proof. Given \((G,k)\), build \((V,d)\) as in 1. For any \(k\)-subset \(S\), \[\binom{k}{2}\cdot 2^{-s}\le E_s(S)\le \binom{k}{2}.\] Set \(T:=\binom{k}{2}\cdot 2^{-s}\). Then \(G\) has a \(k\)-clique if and only if there exists a \(k\)-subset \(S\) with all pairwise distances equal to \(2\), which holds if and only if \(E_s(S)=T\le T\). This is a polynomial reduction from \(k\)-Clique. The argument uses only the strict distance gap between \(1\) and \(2\), and therefore does not require \(s\ge 1\). For fixed positive integer \(s\), all numbers in this reduction are rational; for a fixed non-integer exponent, the constant \(2^{-s}\) is part of the fixed exact-real problem specification. ◻
Remark 1. Theorem 1 is not the main novelty of the present paper. Its role is to anchor the hard side of the complexity landscape in a self-contained way, while the general-metric NP-hardness result itself already appears in the comparative analysis of Pereverdieva et al. [5].
4 Euclidean hardness in the plane when \(s\) is part of the input↩︎
We now restrict to Euclidean geometry and reduce geometric independent set to RSSP using a global threshold based on the gap between forbidden and admissible distances; see Figure 1.
Definition 4 (Geometric Independent Set). Given a finite point set \(P\subset\mathbb{R}^2\), a threshold \(\delta>0\), and an integer \(k\), decide whether there exists a subset \(S\subseteq P\) with \(|S|=k\) such that \(\|x-y\|\ge \delta\) for all distinct \(x,y\in S\).
This problem is NP-complete in the plane [8], [11], [12]. For a finite set \(X\subseteq P\), define \[(X\times X)_{\ge \delta}:=\{(x,y):x,y\in X,\;x\neq y,\;\|x-y\|\ge \delta\},\] \[(X\times X)_{< \delta}:=\{(x,y):x,y\in X,\;x\neq y,\;\|x-y\|< \delta\}.\] We then set \[D_{\min}(X):= \begin{cases} \min\{\|x-y\|:(x,y)\in (X\times X)_{\ge\delta}\}, & \text{if }(X\times X)_{\ge\delta}\neq\emptyset,\\ +\infty, & \text{otherwise,} \end{cases}\] and \[\delta_{\max}(X):= \begin{cases} \max\{\|x-y\|:(x,y)\in (X\times X)_{<\delta}\}, & \text{if }(X\times X)_{<\delta}\neq\emptyset,\\ -\infty, & \text{otherwise.} \end{cases}\] For the ambient point set \(P\), we abbreviate \[D_{\min}:=D_{\min}(P),\qquad \delta_{\max}:=\delta_{\max}(P).\] Figure 1 illustrates these quantities. The key point is that we need a strict gap between the largest forbidden distance \(\delta_{\max}\) and the smallest admissible distance \(D_{\min}\).
Lemma 2. Let \(\emptyset\neq S\subseteq P\). Then \[\delta_{\max}(S)\le \delta_{\max}(P)<\delta\le D_{\min}(P)\le D_{\min}(S).\] In particular, if both forbidden and admissible pairs occur in \(P\), then \[\delta_{\max}<D_{\min}.\]
Proof. The inclusions \[(S\times S)_{<\delta}\subseteq (P\times P)_{<\delta}, \qquad (S\times S)_{\ge\delta}\subseteq (P\times P)_{\ge\delta}\] immediately imply \[\delta_{\max}(S)\le \delta_{\max}(P), \qquad D_{\min}(P)\le D_{\min}(S),\] with the stated conventions when one of the sets is empty. By definition, every element of \((P\times P)_{<\delta}\) has distance strictly smaller than \(\delta\), and every element of \((P\times P)_{\ge\delta}\) has distance at least \(\delta\). Hence, because \(P\) is finite, \[\delta_{\max}(P)<\delta\le D_{\min}(P).\] Combining the inequalities yields the claim. ◻
Lemma 3 (Auxiliary exponent choice). Let \(c>0\) and let \(0<a<b\). Then the inequality \[cb^{-u}<a^{-u}\] holds precisely for \[u>\frac{\log(c)}{\log(b)-\log(a)}.\] In particular, if \(u_0\) is any real number with \[u_0>\frac{\log(c)}{\log(b)-\log(a)},\] then \(cb^{-u}<a^{-u}\) for every \(u\ge u_0\).
Proof. Since \(0<a<b\), we have \(\log(b)-\log(a)>0\). Taking logarithms, the inequality \[cb^{-u}<a^{-u}\] is equivalent to \[\log(c)<u(\log(b)-\log(a)),\] which is exactly the stated condition. The second claim is immediate. ◻
Before stating the main NP-hardness theorem of this section, we briefly outline the coding conventions for the Euclidean case:
Remark 2 (Encoding convention for the Euclidean reduction). For the following reduction we use the usual exact encoding for Euclidean instances in computational geometry. The input point coordinates and the threshold \(\delta\) of the geometric independent-set instance are rational numbers encoded in binary. Squared Euclidean distances are therefore rational, and distances such as \(D_{\min}\) and \(\delta_{\max}\) are represented exactly by the corresponding realizing pairs and the square roots of their rational squared distances. In the constructed RSSP instance, the exponent \(s\) is a positive integer encoded in binary, and the threshold is represented by the straight-line expression \(T=\binom{k}{2}D_{\min}^{-s}.\). Exponentiation is stored succinctly, for instance by repeated squaring, rather than by expanding numerators and denominators. If \(L\) is the bit length of the geometric input, standard separation bounds for distinct rational squared distances imply that the integer \(s\) chosen below has polynomial bit length in \(L\). Thus the constructed RSSP instance has polynomial encoding length. No decimal approximation of \(D_{\min}\), \(\delta_{\max}\), \(s\), or \(T\) is used in the reduction.
We can now state the central theorem for the Euclidean plane:
Theorem 2. RSSP-Decision in the Euclidean plane is NP-hard when \(s\) is part of the input.
Proof. Let \((P,\delta,k)\) be an instance of geometric independent set. If either \((P\times P)_{<\delta}=\emptyset\) or \((P\times P)_{\ge\delta}=\emptyset\), then the instance is trivial and can be handled directly. We therefore assume that both forbidden and admissible pairs occur. By Lemma 2, we have \(\delta_{\max}<D_{\min}\).
By comparing squared distances, compute in polynomial time the realizing pairs for \(\delta_{\max}\) and \(D_{\min}\), and represent these distances exactly as described above. Choose the positive integer \[s:=1+\left\lceil\frac{\log \binom{k}{2}}{\log(D_{\min}/\delta_{\max})}\right\rceil, \qquad T:=\binom{k}{2}D_{\min}^{-s}.\] The exponent \(s\) is encoded in binary and \(T\) is encoded by the displayed expression. By Lemma 3 with \(c=\binom{k}{2}\), \(a=\delta_{\max}\), and \(b=D_{\min}\), this choice guarantees \[\binom{k}{2}D_{\min}^{-s}<\delta_{\max}^{-s}.\] We claim that for every \(k\)-subset \(S\subseteq P\):
if \(S\) is \(\delta\)-independent, then \(E_s(S)\le T\);
if \(S\) is not \(\delta\)-independent, then \(E_s(S)>T\).
If \(S\) is \(\delta\)-independent, then every selected pair lies in \((S\times S)_{\ge\delta}\), so every selected distance is at least \(D_{\min}(S)\ge D_{\min}\) by Lemma 2. Hence \[E_s(S)\le \binom{k}{2}D_{\min}^{-s}=T.\] Conversely, if \(S\) is not \(\delta\)-independent, then some selected pair has distance strictly smaller than \(\delta\), and therefore at most \(\delta_{\max}(S)\le \delta_{\max}\) by Lemma 2. Hence \[E_s(S)\ge \delta_{\max}^{-s}.\] Therefore \(T<\delta_{\max}^{-s}\). Hence every non-independent \(k\)-subset has energy strictly larger than \(T\).
Thus the GIS instance is a YES-instance if and only if the constructed RSSP instance is a YES-instance. ◻
Remark 3. The essential point in Theorem 2 is the global threshold: the ambient gap between \(\delta_{\max}\) and \(D_{\min}\), illustrated in Figure 1, is converted into an energy separation by choosing \(s\) large enough that a single forbidden distance dominates the total contribution of all admissible pairs. This is also the mechanism by which the reduction points toward the MPD limit.
Remark 4. Since \(\mathbb{R}^2\) embeds isometrically into \(\mathbb{R}^d\) for every fixed \(d\ge 2\) via \((x_1,x_2)\mapsto (x_1,x_2,0,\dots,0)\), Theorem 2 immediately implies that RSSP-Decision is NP-hard in every fixed Euclidean dimension \(d\ge 2\) when \(s\) is part of the input.
We close this section by discussing three obstructions that arise when one tries to remove \(s\) from the input. These show that several standard routes are insightful, but ultimately do not work.
Remark 5. We discuss a direct reduction to 3SAT via Fowler’s gadget [8] in Appendix B. This route is instructive because it shows that one need not use a completely naive separation argument in which a single forbidden pair must dominate the total contribution of all \(\binom{k}{2}\) admissible-pair contributions in a \(\delta\)-independent \(k\)-subset. Owing to the constant geometric gap in that construction — for instance, a separation of at least \(2/3\) between below-threshold and above-threshold distances — the exponent \(s\) can be chosen to grow only linearly with \(k\), rather than proportionally to \(\binom{k}{2}\). Thus the construction shows that the exponent need not compensate for the full quadratic number of pair contributions in an admissible subset. However, this still does not yield a fixed-\(s\) hardness result. Although the required exponent grows only linearly with \(k\), it remains tied to the subset size and hence to the input instance. Therefore this route improves substantially over the naive quadratic estimate, but still requires choosing \(s\) as part of the construction rather than keeping it uniformly constant.
Remark 6. One might want to use a scaled distance matrix to achieve separation between admissible and non-admissible cases for fixed \(s\). The main obstruction is the homogeneity of the Riesz kernel: \(E_s(\lambda X)=\lambda^{-s}E_s(X)\). Thus, if all relevant distances are scaled by the same factor \(\lambda\), then the energies of admissible and non-admissible subsets are both multiplied by the same factor \(\lambda^{-s}\). Hence uniform scaling does not improve their relative separation. In particular, a distance pattern \((\lambda,2\lambda)\) yields only the constant factor \(2^s\) between one forbidden pair and one admissible pair. More precisely, if in an admissible subset all pairs have distance at least \(2\lambda\), whereas a non-admissible subset contains at least one pair at distance at most \(\lambda\), then one forbidden pair contributes at least \(\lambda^{-s}\), while an admissible \(k\)-subset may still have total energy as large as \(\binom{k}{2}(2\lambda)^{-s}\). Thus one would need \(2^s>\binom{k}{2}\) for a single forbidden pair to dominate the full energy of every admissible \(k\)-subset. For fixed \(s\), this fails once \(k\) is sufficiently large. Therefore such a simple \((\lambda,2\lambda)\) scaling cannot yield the required separation for fixed \(s\).
Remark 7 (Heuristic obstruction to broadcast amplification). A plausible route to fixed-\(s\) hardness would be to replicate “check” gadgets in additional geometric layers, with the aim that one incompatible logical pair should trigger many short distances. However, if the replicated check points are fixed in the instance and are not themselves selected, then their contribution to the Riesz \(s\)-energy factorizes. Indeed, if the selected representatives are \(x_u\) and \(y_v\) and \(Z\) is any fixed auxiliary set, then \[E(u,v)=\sum_{z\in Z}\|x_u-z\|^{-s} +\sum_{z\in Z}\|y_v-z\|^{-s} +\|x_u-y_v\|^{-s} +\mathrm{const}.\] Hence, the auxiliary layers contribute only unary terms depending on \(u\) and \(v\) separately; they do not create new pair-specific penalties beyond the direct interaction \(\|x_u-y_v\|^{-s}\). Therefore, merely stacking fixed check layers in a third dimension does not by itself amplify a single bad logical pair into many independent bad terms. Any successful fixed-\(s\) construction would need an additional mechanism that forces the co-selection of companion points, or else must encode the amplification directly in the distances among the selected points themselves.
Taken together, these three obstructions indicate that eliminating \(s\) from the input is a genuinely nontrivial challenge: while each route is conceptually informative, none presently leads to a fixed-\(s\) hardness construction in Euclidean space.
5 Ultrametrics and exact dynamic programming↩︎
Ultrametrics are metrics induced by rooted hierarchies. Every finite ultrametric space admits a rooted-tree representation: the points are the leaves, internal nodes carry nondecreasing heights, and the distance between two leaves is determined by the height of their least common ancestor. Such spaces arise naturally from rooted phylogenetic trees in biodiversity applications, where distances encode species separation in evolutionary time, and from the cophenetic distances of hierarchical clustering dendrograms, including hierarchical analyses of cultural vectors [9], [10]. They therefore form a natural structured class in which one may ask whether minimum Riesz \(s\)-energy subset selection remains hard or becomes tractable.
A recent exact algorithmic contribution is the pruning algorithm of [13], which computes a size-\(k\) subset of maximum phylogenetic diversity (PD) on a phylogenetic tree. In the specifically ultrametric setting, exact algorithmic results appear to be comparatively sparse; among the verified references included here, [14] gives polynomial-time methods for counting all size-\(k\) maximum-PD sets on rooted ultrametric trees and for optimizing a linear score over those maximum-PD sets. A broader exact ultrametric framework is provided by [15], who prove polynomial-time solvability for selecting a size-\(k\) subset maximizing any weakly subset-monotone diversity function extending an ultrametric; this applies in particular to phylogenetic diversity in the tree setting, as well as to the pairwise-sum diversity \(\delta_{\mathrm{sum}}\), the minimum-pairwise-distance diversity \(\delta_{\min}\), and Weitzman diversity \(\delta_W\).
Also for Riesz \(s\)-energy subset selection, we also show that finite ultrametric spaces form an exact, tractable island. The key structural fact is that, at any internal node of the rooted tree, all cross-subtree distances are equal. Hence the full interaction between two child subproblems collapses to a single cardinality-dependent cross term.
Theorem 3 (Ultrametric recurrence). Let \(F_u(t)\) denote the minimum Riesz \(s\)-energy among all subsets of size \(t\) chosen from the leaves below node \(u\) in a rooted binary ultrametric tree. Then for a leaf \(\ell\), \[F_\ell(0)=0,\qquad F_\ell(1)=0,\] and for an internal binary node \(u\) with children \(v\) and \(w\), \[F_u(t)=\min_{t_v+t_w=t}\bigl(F_v(t_v)+F_w(t_w)+t_v t_w\,\Delta_u^{-s}\bigr),\] where \(\Delta_u\) is the common distance between any leaf in the subtree of \(v\) and any leaf in the subtree of \(w\).
Corollary 1. On a rooted binary ultrametric tree with \(n\) leaves, minimum Riesz \(s\)-energy subset selection can be solved exactly in \(O(nk^2)\) time and \(O(nk)\) space.
Remark 8. Theorems 3 and Corollary 1 identify finite ultrametric spaces as the main exact tractable regime isolated in this paper. A fuller derivation of the recurrence, an explicit dynamic-programming algorithm, a biological worked example, and a brief numerical validation against brute force are deferred to Appendix 10.
6 The large-\(s\) limit and minimum pairwise distance (MPD)↩︎
For the comparison with minimum pairwise distance (MPD), we recall the following known finite-instance observation; see, for example, Atta [3].
Proposition 1 (Eventual minimum pairwise distance equivalence). Let \((X,d)\) be a finite metric space and fix \(k\). Then there exists \(s_0>0\) such that for every \(s\ge s_0\), every minimizer of \(E_s\) over all \(k\)-subsets of \(X\) is also a maximizer of the MPD objective \(\operatorname{MPD}\!\left(S\right)\).
Proof. We use a similar reasoning than in the proof of theorem 2. Let \(\mathcal{S}\) be the finite family of all \(k\)-subsets of \(X\), and let \[D^*:=\max\{\operatorname{MPD}\!\left(S\right):S\in\mathcal{S}\}.\] Choose any \(T\in\mathcal{S}\) with \(\operatorname{MPD}\!\left(T\right)=D^*\). Then \[E_s(T)\le \binom{k}{2}(D^*)^{-s}.\] Now let \[R:=\max\{\operatorname{MPD}\!\left(S\right):S\in\mathcal{S},\;\operatorname{MPD}\!\left(S\right)<D^*\}.\] Because \(\mathcal{S}\) is finite, we have \(R<D^*\). For every \(S\in\mathcal{S}\) with \(\operatorname{MPD}\!\left(S\right)<D^*\), some selected pair has distance at most \(R\), hence \[E_s(S)\ge R^{-s}.\] Choose \(s_0\) large enough that \[R^{-s}>\binom{k}{2}(D^*)^{-s} \qquad\text{for all }s\ge s_0.\] Then for all such \(s\), every subset with \(\operatorname{MPD}\!\left(S\right)<D^*\) has strictly larger Riesz energy than \(T\). Therefore every minimizer of \(E_s\) must satisfy \(\operatorname{MPD}\!\left(S\right)=D^*\), that is, it must be MPD optimal. ◻
Remark 9. This known finite-instance observation is central to our comparison with minimum pairwise distance (MPD). It is interpretive rather than a new uniform complexity result: it explains why many large-\(s\) phenomena naturally converge to the classical MPD picture.
7 Open problems and boundary of the tractable regime↩︎
The preceding theorems draw a sharp boundary, but not a complete classification. The remaining open cases are structurally meaningful rather than incidental.
7.0.0.1 Fixed-exponent Euclidean hardness.
Our Euclidean hardness proof uses \(s\) as part of the input. Extending NP-hardness to every fixed exponent \(s>0\) in low-dimensional Euclidean spaces remains open. This is arguably the most natural missing theorem on the hard side. Appendix 11 records a Fowler-style fixed-\(s\) route and explains precisely why the present admissible-versus-forbidden threshold estimates do not yet close such a reduction.
7.0.0.2 One-dimensional and ordered Euclidean inputs.
The line (and also \(\ell_1 ordered staircases\), see [16]) appears significantly more structured than arbitrary Euclidean instances, but one should be careful not to overinterpret this. The contrast with minimum pairwise distance (MPD) is already sharp in low Euclidean dimension. For the MPD objective on sorted points \(x_1<\dots<x_n\), an exact Bellman-style recurrence can be written in terms of the last chosen point: \[M(i,t)=\max_{j<i}\min\{M(j,t-1),\,x_i-x_j\},\] where \(M(i,t)\) denotes the best possible bottleneck value of a \(t\)-subset ending at \(x_i\). Thus the MPD optimization problem is polynomial-time solvable on the line by dynamic programming, whereas in dimension two it is already NP-hard via the same geometric independent-set machinery that underlies our plane hardness discussion [17]. By contrast, no analogous left/right state compression works for minimum Riesz \(s\)-energy: when a new point \(x_i\) is added, the extra cost is \[\sum_{x\in S}|x_i-x|^{-s},\] which depends on the full internal configuration of the already selected subset \(S\), not merely on its cardinality or its extreme points. In particular, the naive Bellman principle that suffices for MPD on the line breaks down for Riesz energy (see [18], where also a concrete counterexample is provided). Thus the one-dimensional Euclidean case remains open here for a substantive structural reason, rather than because of a missing implementation detail.
Notably, however, the Riesz \(s\)-energy on ordered point sets exhibits considerable positive structure. For example, it is straightforward to show that it satisfies a Monge condition, which can be exploited to obtain algorithmic speedups and certificates for multilevel dynamic programming. This may be useful, for instance, when considering higher-order neighborhoods in the Bellman recursion. A detailed discussion of these aspects, however, exceeds the scope of the present paper and is left for future work.
In the limit \(s \rightarrow \infty\), the picture becomes clearer, because the limiting maximum minimum pairwise distance (MPD) problem on the line admits a polynomial time algorithm based on a simple greedy feasibility test. To test whether there exists a \(k\)-subset with MPD at least \(\tau\), start with the leftmost point and then scan the points from left to right, always selecting the first point whose distance from the last selected point is at least \(\tau\). At the end of the scan, one checks whether at least \(k\) points have been selected. This greedy test is exact: choosing each feasible point as far to the left as possible cannot reduce the remaining feasible space to its right. The optimal MPD value is one of the pairwise distances \(x_j-x_i\), \(1\le i<j\le n\). Hence an exact implementation can explicitly form and sort all pairwise distances and then apply binary search with the greedy feasibility test. This gives an \(O(n^2\log n)\) time bound and \(O(n^2)\) storage. Alternatively, one may apply the Frederickson–Johnson selection and ranking algorithm for sorted matrices to the implicit sorted matrix of pairwise differences [19]. This avoids materializing the \(O(n^2)\) pairwise distances and yields an \(O(n\log n)\)-type exact bound after the input points have been sorted, but requires the nontrivial sorted-matrix-search machinery.
7.0.0.3 Parameterized complexity.
The problem also invites parameterized analysis in \(k\), dimension, or structural parameters of thresholded interaction graphs. We leave this aside here in order to keep the paper focused on the primary hard-versus-tractable dichotomy. Please refer to [3] for a first result in the direction of PTAS constructions.
8 Complexity landscape↩︎
The results of this paper isolate a compact but informative landscape for minimum Riesz \(s\)-energy subset selection. On the hard side, the problem is NP-hard already on finite metric spaces for every fixed \(s>0\), and it remains NP-hard in the Euclidean plane — hence in every fixed Euclidean dimension \(d\ge 2\) — when the exponent is part of the input. Thus neither unrestricted metric structure nor low-dimensional Euclidean geometry is, by itself, enough to make the problem easy.
On the tractable side, finite ultrametric spaces admit exact dynamic programming because the interaction between two child subproblems depends only on the selected cardinalities and the common cross-subtree distance. This identifies ultrametrics as a genuine tractable island rather than a degenerate case. For the comparison with minimum pairwise distance (MPD), a central ingredient is the known finite-instance large-\(s\) observation that every Riesz \(s\)-energy minimizer is eventually MPD optimal. This explains why large-\(s\) behavior increasingly resembles bottleneck dispersion even though the finite-\(s\) objective is globally additive over all pairs.
The low-dimensional Euclidean comparison is particularly revealing. For minimum pairwise distance (MPD) on the line, the ordered structure supports a Bellman-style recurrence, whereas for Riesz \(s\)-energy the same state compression fails because adding a new point contributes a sum over all previously selected points. Thus the one-dimensional Riesz case remains open for a structural reason, not merely because the right implementation has not yet been found. Table 1 summarizes this picture.
| Setting | Minimum Riesz \(s\)-energy subset selection | Large-\(s\) / minimum-pairwise-distance viewpoint |
|---|---|---|
| General finite metric spaces | NP-hard for every fixed \(s>0\) (Theorem 1); general-metric hardness already appears in [5] | NP-hard [17] |
| Euclidean plane, \(s\) in input | NP-hard (Theorem 2) | Consistent with the fact that the MPD optimization problem is already NP-hard in the plane |
| Euclidean plane, fixed exponent \(s\) | Open in this paper; for \(s>2\), Appendix 11 gives an \(O(k)\) far-field bound for Fowler-style gadgets. | NP-hard in 2-D [17] |
| Finite ultrametric spaces | Exact DP in \(O(nk^2)\) (Theorem 3, Corollary 1) | Exact DP in \(O(nk^2)\); comparison uses the known finite-instance large-\(s\) link to MPD optimality (Proposition 1) |
| One-dimensional Euclidean inputs | Open in this paper; left–right DP is not exact for Riesz energy. | Polynomial-time solvable for MPD, e.g., by a Bellman-style left–right dynamic program; |
9 Conclusion↩︎
Minimum Riesz \(s\)-energy subset selection exhibits a concise and, to our knowledge, previously unrecorded complexity picture. Building on the earlier general-metric hardness context established in the comparative-analysis work of Pereverdieva et al. [5], we showed that the problem remains hard in the Euclidean plane when the exponent \(s\) is part of the input. At the same time, finite ultrametric spaces admit exact dynamic programming, showing that hierarchical metrics recover nontrivial algorithmic structure. Together with the known finite-instance large-\(s\) link to minimum pairwise distance (MPD), this shows that the problem naturally interpolates between pairwise repulsion and bottleneck-style separation. This makes the boundary between general geometry, the one-dimensional ordered setting, and hierarchical structure especially interesting: ultrametrics are exactly solvable, the Euclidean plane is hard, and the one-dimensional case remains structurally subtle because the natural state compression that works for MPD does not directly transfer to Riesz \(s\)-energy. Beyond exact complexity, recent work also points toward the need to explore practical algorithmic directions, including fast Riesz-based subset-selection heuristics in multiobjective optimization and approximation results for related geometric facility-location variants [3], [18], [20].
10 Appendix A: Dynamic programming on finite ultrametric spaces↩︎
This appendix gives a fuller discussion of the ultrametric dynamic programming result stated in the main text. We first derive the recurrence, then formulate the corresponding algorithm explicitly, and finally illustrate it on a small abstract tree example and report a brief numerical validation.
10.1 Derivation of the recurrence↩︎
Let \(u\) be an internal node of a rooted binary ultrametric tree with children \(v\) and \(w\). Denote by \(L(u)\) the leaves below \(u\), and let \(\Delta_u\) be the common distance between every leaf in \(L(v)\) and every leaf in \(L(w)\). If \(A\subseteq L(v)\) and \(B\subseteq L(w)\), then every cross pair contributes the same amount \(\Delta_u^{-s}\) to the Riesz energy. Therefore \[E_s(A\cup B)=E_s(A)+E_s(B)+|A||B|\,\Delta_u^{-s}.\] This identity is the only nontrivial structural ingredient. Once it is available, the Bellman recurrence follows immediately by minimizing over all splits of a target cardinality.
Proof of Theorem 3. Any subset \(S\subseteq L(u)\) of size \(t\) decomposes uniquely as \(S=A\cup B\) with \(A\subseteq L(v)\) and \(B\subseteq L(w)\). Every pair inside \(A\) contributes to \(E_s(A)\), every pair inside \(B\) contributes to \(E_s(B)\), and every cross pair contributes the same amount \(\Delta_u^{-s}\). Since there are \(|A||B|\) cross pairs, we obtain \[E_s(S)=E_s(A)+E_s(B)+|A||B|\,\Delta_u^{-s}.\] Minimizing over all splits \(t=|A|+|B|\) yields the stated recurrence. The leaf initialization is immediate from the definition of Riesz \(s\)-energy. ◻
Proof of Corollary 1. The recurrence is evaluated bottom-up. For each internal node and each target cardinality \(t\le k\), all splits \(t_v+t_w=t\) are tried. This is the same counting as in a knapsack convolution, hence \(O(k^2)\) work per internal node and \(O(nk^2)\) total work over the tree. Storing one table of size \(O(k)\) per node gives \(O(nk)\) space. ◻
10.2 Algorithmic form↩︎
Algorithm 2 records the resulting exact dynamic program. A fully documented reference implementation and validation suite for this algorithm is available at https://github.com/emmerichmtm/RieszSEnergyUltrametricsDP.
10.3 Worked abstract tree example↩︎
We illustrate the recurrence on a small abstract rooted tree with six taxa. The tree consists of two cherries, \(\{a,b\}\) and \(\{c,d\}\), together with two further taxa \(e\) and \(f\). The first cherry merges at height \(3.5\), the second at height \(4\), and the two resulting groups merge at height \(5.5\). Figure 3 shows the corresponding stylized ultrametric tree. Such numbers may be interpreted abstractly as times in a phylogenetic tree, as in biodiversity applications, or as cophenetic heights in a hierarchical clustering of cultural vectors [9], [10].
If leaf-to-leaf distances are taken as twice the least-common-ancestor heights, then \[d(a,b)=7, \qquad d(c,d)=8, \qquad d(x,y)=11\] for all other distinct pairs \(x,y\). Take \(s=1\). Then the pair \(\{a,b\}\) has energy \(1/7\), the pair \(\{c,d\}\) has energy \(1/8\), and every pair drawn from different top-level branches has energy \(1/11\). For three-taxon subsets, one obtains for instance \[E_1(\{a,c,e\})=\frac{1}{11}+\frac{1}{11}+\frac{1}{11}=\frac{3}{11}\approx 0.272727,\] whereas \[E_1(\{a,b,e\})=\frac{1}{7}+\frac{1}{11}+\frac{1}{11}=\frac{25}{77}\approx 0.324675.\] Thus the energy-minimizing three-taxon subsets avoid the closer cherry pair \(\{a,b\}\) whenever possible, exactly as the recurrence predicts.
10.4 Numerical validation and implementation↩︎
To validate the recurrence computationally, we implemented the ultrametric DP and compared its output against brute-force enumeration on a deterministic set of small random rooted binary ultrametric trees. On all tested instances the dynamic-programming optimum matched the brute-force optimum up to floating-point roundoff. A small reference Python implementation accompanies the present work See: https://github.com/emmerichmtm/RieszSEnergyUltrametricsDP.
11 Appendix B: A Fowler-style fixed-\(s\) route in the Euclidean plane and its limitation↩︎
This appendix records an alternative route toward Euclidean hardness at a fixed exponent \(s\), inspired by the classical gadget construction of Fowler, Paterson, and Tanimoto [8] and translated here from unit squares to unit discs. The aim is not to claim a complete fixed-\(s\) hardness proof, but to document a geometrically natural approach and to isolate precisely where the present estimates stop short. At the same time, the analysis yields a useful structural fact: in the relevant gadget class, the admissible far-field Riesz energy is only linear in the number of selected points, rather than quadratic as a naive pair-count estimate would suggest.
11.1 Geometric idea↩︎
The intended reduction route is directly from 3SAT. As in Fowler et al., one arranges wire gadgets with two states, crossover gadgets preserving the transmitted state, and clause gadgets that allow one additional choice if and only if a literal is satisfied. In the present context, the choices are represented by points in the Euclidean plane, or equivalently by equal-radius discs centred at those points.
The geometric idea is simple. A forbidden local configuration forces two selected discs to come too close, creating a very short centre-to-centre distance and therefore a large Riesz contribution. One may therefore hope to separate satisfying from unsatisfying assignments by making every inconsistent choice incur such an overlap penalty, while every consistent choice incurs only smaller contributions from more distant interactions.
Figures [fig:boolean-circuit-app]–4 reproduce the figure structure from the longer manuscript and illustrate this idea.
![Figure 4: Disk-based gadget layout in the Euclidean plane following the structure of Fowler, Paterson, and Tanimoto [8]: a Boolean circuit encoding of a 3SAT instance, wire gadgets with two states, crossover gadgets preserving the transmitted state, and clause gadgets in which one extra disk fits if and only if a literal is satisfied. In a full fixed-s proof attempt, admissible selections would correspond to non-overlapping choices consistent with a satisfying assignment.](/article/2605.04715/_pdflatex/yruihctd.png "fig:")
We arrange the disc centres so that every admissible selected pair is at distance at least \(2r\), while any forced conflict creates a pair at distance at most \(\tfrac{3}{2}r\). Thus a forbidden overlap contributes at least \[\left(\tfrac{3}{2}r\right)^{-s},\] whereas an admissible pair contributes at most \[(2r)^{-s}.\] This produces a constant local gap.
Lemma 4 (Local overlap gap). At exponent \(s>0\), any forced-overlap pair of selected discs whose centres have distance at most \(\tfrac{3}{2}r\) contributes at least \[\left(\tfrac{3}{2}r\right)^{-s}\] to the Riesz energy, while any admissible non-overlap pair whose centres have distance at least \(2r\) contributes at most \[(2r)^{-s}.\]
Proof. This is immediate from monotonicity of the map \(d\mapsto d^{-s}\) for \(s>0\). ◻
Lemma 4 explains why this route is attractive: forbidden local configurations are penalized more heavily than admissible ones. To obtain a reduction, however, one must compare a single such overlap penalty with the entire energy of an admissible \(k\)-point selection. This is where the global budget enters.
11.2 Packing-based far-field budget in the Euclidean plane↩︎
To bound the admissible side globally, we use the standard packing estimate for equal-radius discs in the Euclidean plane. Figure 5 is the same packing picture used in the longer draft.
None
Figure 5: Packing of equal-radius discs in the Euclidean plane. Around a fixed selected disc, the number of non-overlapping discs in the \(i\)th layer grows at most linearly in \(i\). This leads to a convergent sum for every fixed \(s>2\)..
Fix one selected disc with centre \(p\). Around this disc, consider the concentric annuli \[A_i(p)=\bigl\{x\in\mathbb{R}^2 : 2ir \le \|x-p\| < 2(i+1)r\bigr\}, \qquad i\ge 1.\] Because the selected discs are non-overlapping, only linearly many other disc centres can lie in the \(i\)th annulus.
Lemma 5 (Pointwise far-field budget). Fix \(s>2\). Let \(p\) be the centre of a selected disc in a non-overlapping configuration of equal-radius discs with radius \(r\), and let \(\zeta\) denote the Riemann zeta function. Then \[\sum_{q\in S\setminus\{p\}} \|p-q\|^{-s} \le 6(2r)^{-s}\sum_{i\ge 1} i^{1-s} = 6(2r)^{-s}\zeta(s-1).\] In particular, for every fixed \(s>2\), the total contribution of all other selected points to one selected point is bounded by a constant depending only on \(s\) and on the geometric scale \(r\).
Proof. In the Euclidean plane, at most \(6i\) non-overlapping equal-radius discs can have centres in the \(i\)th annulus around the centre \(p\) of a fixed selected disc. Every such centre is at distance at least \(2ir\) from \(p\), hence each such point contributes at most \[(2ir)^{-s}.\] Therefore the total contribution from the \(i\)th annulus is at most \[6i\,(2ir)^{-s}=6(2r)^{-s}i^{1-s}.\] Summing over all annuli gives \[\sum_{q\in S\setminus\{p\}} \|p-q\|^{-s} \le 6(2r)^{-s}\sum_{i\ge 1} i^{1-s} = 6(2r)^{-s}\zeta(s-1),\] which converges for every fixed \(s>2\). ◻
Remark 10 (Basel problem special case). For the concrete case \(s=3\), Lemma 5 becomes \[\sum_{q\in S\setminus\{p\}} \|p-q\|^{-3} \le 6(2r)^{-3}\sum_{i\ge 1} i^{-2} = 6(2r)^{-3}\zeta(2) = 6(2r)^{-3}\frac{\pi^2}{6}.\] Thus the classical Basel value \(\zeta(2)=\pi^2/6\) appears naturally as the packing-based bound on the total contribution to one selected point.
The pointwise estimate immediately implies a linear global bound.
Corollary 2 (Linear admissible-energy budget). For every fixed \(s>2\), there exists a constant \(C_s>0\) such that every admissible non-overlapping \(k\)-point selection \(S\) in this gadget class satisfies \[E_s(S)\le C_s\,k.\]
Proof. Write the unordered Riesz energy as \[E_s(S) = \frac{1}{2}\sum_{p\in S}\sum_{q\in S\setminus\{p\}}\|p-q\|^{-s}.\] By Lemma 5, each inner sum is bounded by the same constant \[6(2r)^{-s}\zeta(s-1).\] Hence \[E_s(S)\le \frac{1}{2}\,k\cdot 6(2r)^{-s}\zeta(s-1),\] which proves the claim. ◻
Remark 11. Corollary 2 is the key structural message of this appendix. A naive global estimate would bound the admissible energy by \(O(k^2)\) merely because there are \(\binom{k}{2}\) pairs. The packing argument shows that this is far too crude: in the present Euclidean-plane gadget class, the admissible far-field energy is only \(O(k)\).
11.3 Why this still does not close a fixed-\(s\) reduction↩︎
Suppose now that a selected set contains at least one forbidden overlap pair. Then by Lemma 4 its energy is at least \[\left(\tfrac{3}{2}r\right)^{-s},\] plus the contributions of all remaining pairs. On the admissible side, Corollary 2 gives the upper bound \[E_s^{\mathrm{adm}}(k)\le C_s\,k.\] Thus the overlap penalty must be compared not with a quadratic admissible budget, but with a linear one.
Proposition 2 (Interpretation of the obstruction). The Fowler-style route in the Euclidean plane shows that the obstruction to a fixed-\(s\) reduction is not a quadratic explosion of admissible pair contributions. Rather, the admissible-energy budget is only linear in \(k\), but still grows with \(k\). Therefore a uniform reduction at fixed \(s\) would require overlap penalties strong enough to dominate that linear growth.
Proof. By Corollary 2, every admissible \(k\)-point selection has energy at most \(C_s k\). A single forbidden overlap contributes at least \(\left(\tfrac{3}{2}r\right)^{-s}\) by Lemma 4. For a reduction based solely on this gap, one would need the latter quantity to dominate the former uniformly in \(k\). Since the admissible budget still grows linearly with \(k\), the present estimates do not by themselves yield a fixed-\(s\) reduction. ◻
This is the main point of the appendix. The failure of the fixed-\(s\) route is subtler than a naive pair-count argument would suggest. The obstruction is not that admissible long-range interactions explode quadratically. On the contrary, packing controls the per-point contribution by a constant and hence the total admissible energy linearly in \(k\). The real issue is that even this sharpened linear budget must still be dominated by the penalty of a forbidden overlap.
This also explains why the successful Euclidean hardness proof lets the exponent depend on \(k\). Increasing \(s\) amplifies the ratio between the forbidden short distance \(\tfrac{3}{2}r\) and the admissible distance \(2r\), so that one forbidden pair can dominate the entire linear admissible-energy budget.
Remark 12. The value of the Fowler-style route is therefore twofold. First, it provides a natural geometric mechanism for encoding inconsistency by overlap penalties in the Euclidean plane. Second, it clarifies why the successful reduction uses an exponent depending on \(k\): not because the admissible far-field energy is uncontrolled, but because even its sharp linear bound still grows with \(k\).