The optimal packing of eight points in the real projective plane
Dustin G. Mixon Hans Parshall




How can we arrange n lines through the origin in three-dimensional Euclidean space in a way that maximizes the minimum interior angle between pairs of lines?  Conway, Hardin and Sloane (1996) produced line packings for n ≤ 55 that they conjectured to be within numerical precision of optimal in this sense, but until now only the cases n ≤ 7 have been solved. In this paper, we resolve the case n = 8.  Drawing inspiration from recent work on the Tammes problem, we enumerate contact graph candidates for an optimal configuration and eliminate those that violate various combinatorial and geometric necessary conditions.  The contact graph of the putatively optimal numerical packing of Conway, Hardin and Sloane is the only graph that survives, and we recover from this graph an exact expression for the minimum distance of eight optimally packed points in the real projective plane.



Convergence in uncertain linear systems
    [
    Received: date / Accepted: date
=======================================





§ INTRODUCTION


Consider the fundamental line packing problem of packing n points in 𝐑𝐏^d-1 or 𝐂𝐏^d-1 so that the minimum distance is maximized.
We refer to such n-point sets as projective n-packings, and we take the distance between two points in projective space to be the interior angle of the corresponding lines in affine space.
The line packing problem has received considerable attention over the last century.
The case of 𝐂𝐏^1 is equivalent to packing points in S^2, namely, the Tammes problem tammes30.
The general line packing problem was originally formulated in fejes65 and subsequently studied in welch74,delsarte75,levenshtein82.
Several putatively optimal packings in 𝐑𝐏^d-1 were later provided in conway96, each conjectured to be within numerical precision of an optimal packing.
Most of these conjectures remain open, even despite a surge of progress over the last decade (see fickus18 and references therein) that has been motivated in part by emerging applications in compressed sensing bandeira13, digital fingerprinting mixon13, quantum state tomography renes04, and multiple description coding strohmer03.

In the present paper, we consider the special case of projective n-packings in 𝐑𝐏^2, that is, sets of n lines through the origin in 𝐑^3.  It is convenient to represent a projective n-packing by a set of n unit vectors Φ⊆ S^2 spanning the corresponding lines, where two such sets represent the same packing if and only if one can be obtained from the other by negating some of the vectors.  We denote the resulting equivalence class by [Φ].  The coherence of the projective n-packing Φ is defined as

    μ(Φ) := max_x,y∈Φ
    x≠ y |⟨ x,y ⟩|,

which satisfies μ(Φ')=μ(Φ) for every Φ'∈[Φ].  Setting

    μ_n := inf{μ(Φ) : Φ⊆ S^2, |Φ| = n},

we say that a projective n-packing Φ is optimal if μ(Φ) = μ_n.  The existence of optimal projective n-packings is guaranteed for every n by a compactness argument, but they have only been classified for n ≤ 7 fejes65,conway96,cohn12; see <ref> for an illustration.  In this paper, we compute μ_8 and identify an optimal projective 8-packing, which is unique up to isometry.




The projective 8-packing with minimum coherence μ_8 is unique up to isometry.  Furthermore, μ_8 is the third-largest root of

    1 + 5x - 8x^2 - 80x^3 - 78x^4 + 146x^5 - 80x^6 - 584x^7 + 677x^8 + 1537x^9,

which is given numerically by μ_8 ≈ 0.6475889787.


This improves upon the prior state of the art, which is summarized by the following lemma.
The lower bound follows from the Levenshtein bound levenshtein82, while the upper bound follows from an explicit numerical packing given in conway96, which was recently given exact coordinates in mixon19.


The coherence μ_8 of an optimal projective 8-packing satisfies 

    0.6 ≤μ_8 ≤ 0.647588979.



A set X ⊆ S^2 is said to be antipodal if X = -X, and we say X is an antipodal 2n-packing if X is antipodal with |X|=2n.
We may identify any projective n-packing Φ with an antipodal 2n-packing, namely

    X(Φ) := {x ∈ S^2 : x ∈Φ or -x∈Φ}.

For x,y ∈ S^2, we define their geodesic distance to be

    d(x,y) := arccos(⟨ x,y ⟩),

and we define the minimum distance of any finite X ⊆ S^2 to be

    ψ(X) := min_x,y ∈ X 
     x ≠ y d(x,y).

Of course, minimizing μ(Φ) over projective n-packings Φ is equivalent to maximizing ψ(X) over antipodal 2n-packings X.  Both formulations play a role in our proof of the main result.

The second formulation bears resemblance to the Tammes problem, originating in tammes30, of maximizing ψ(X) over arbitrary packings of n points X ⊆ S^2.  The earliest solutions to the Tammes problem, for n = 3, 4, 6, 12, are given by regular triangulations of the sphere that achieve equality in a general inequality on the optimal minimum distance fejes49.  For sufficiently large values of n, this inequality was sharpened in robinson61, leading to the solution for n = 24.  Further progress has been made by studying the contact graph of X ⊆ S^2, denoted G(X), which is drawn on the sphere with vertices X and with a geodesic edge between x,y ∈ X exactly when d(x,y) = ψ(X).  Classifying contact graphs led to solutions to the Tammes problem for n = 5, 7, 8, 9 schutte51 and n = 10, 11 danzer86.  More recently, computer enumeration and optimization of over a billion contact graph candidates led to the solutions for n = 13, 14 musin12,musin15.  The authors write in related work musinSteklov15 that “the direct approach to the solution of the Tammes problem based on computer enumeration of irreducible contact graphs … has actually been exhausted.”  All other cases of the Tammes problem remain open.  For instance, the n=20 case is open, though it has long been known that the vertices of the dodecahedron are suboptimal vdw52.

To prove our main result, we adapt the approach of contact graph enumeration and elimination.  To make the enumeration feasible, we leverage the antipodal constraint to significantly reduce the number of graphs that we need to generate and consider.  
To this end, Section <ref> identifies a set Γ_1 of 547 graphs with the property that every maximal connected set of optimal antipodal 16-packings has a representative packing with contact graph in Γ_1.
In Section <ref>, we consider all planar embeddings of members of Γ_1 modulo homeomorphism and reflection.
These equivalence classes can be represented as combinatorial embeddings, which we identify using the SPQR-tree data structure developed in dibattista89.
After removing embeddings that violate various necessary conditions on the contact graph's facial structure, we arrive at a set ℰ_1 of 217 combinatorial embeddings of contact graph candidates.
Section <ref> then imposes geometric constraints in order to eliminate all but four of these candidates, three of which are subgraphs of the fourth.
We then demonstrate in Section <ref> that the three subgraphs violate the Karush–Kuhn–Tucker optimality conditions, thereby isolating a single contact graph candidate G. For the sake of reproducibility, we offer a SageMath worksheet code for our computations.
By [Theorem II.2]mixon19, there is a unique antipodal 2n-packing (up to isometry) with coherence satisfying Lemma <ref> and whose contact graph is isomorphic to G.
The singleton set containing this packing is a maximal connected set of optimal packings, thereby implying uniqueness.
We conclude in Section <ref> with a discussion of possible directions for future work.



§ ENUMERATING CONTACT GRAPH CANDIDATES


In this section, we discuss various combinatorial features of the contact graph of an antipodal 2n-packing, and then we apply these features to enumerate contact graph candidates for the case where n=8.
First, for any X⊆ S^2, notice that no two edges of G(X) may cross, since otherwise, a simple application of the triangle inequality produces two points in X that are closer than ψ(X) apart.  This implies the following.


	For any X ⊆ S^2, G(X) is a planar graph.


We define the projective contact graph G([Φ]) of a projective n-packing Φ to have vertex set Φ and x↔ y precisely when |⟨ x,y⟩| = μ(Φ).  The contact graph G(X(Φ)) of the corresponding antipodal packing provides a canonical double cover of G([Φ]), where the covering map amounts to identifying antipodal pairs.  This covering map naturally provides a projective planar embedding of G([Φ]).


	For any projective packing Φ, G([Φ]) is a projective planar graph.


For finite X ⊆ S^2, we say x ∈ X can be shifted if every open neighborhood U ⊆ S^2 of x contains a point y ∈ U such that d(x',y) > ψ(X) for all x' ∈ X ∖{x}.  If the isolated vertices of G(X) are the only ones that can be shifted, then we call G(X) irreducible.  Irreducible contact graphs were introduced in schutte51 in order to study the Tammes problem.  They made the following straightforward observations.

 Suppose X ⊆ S^2 with |X| > 6 and G(X) irreducible.
	
		
  (i) If x ∈ X is not isolated in G(X), then 3 ≤(x) ≤ 5.
		
  (ii) Each face of G(X) is a convex spherical polygon.
		
  (iii) The length of each face of G(X) is at most ⌊ 2π / ψ(X) ⌋.
	


The proof of Lemma <ref>(iii) hinges on the following useful lemma.


For every convex spherical B⊆ S^2 and convex spherical polygon A⊆ B, the perimeter of A is at most the perimeter of B.



There exist finitely many hemispheres {H_i} such that A=⋂_i H_i.
Put S_0:=B and S_k+1:=S_k∩ H_k+1 for each k≥ 0.
Then by the triangle inequality, the perimeter of S_k+1 is at most that of S_k for every k, implying the result.



If (x) = 1, then x can be shifted.  Similarly, if (x) = 2 and ψ(X) < π/2, which holds for |X| > 6, then x can be shifted.
If (x) > 5, then the smallest angle between two edges emanating from x is at most π/3.  This forces the corresponding vertices adjacent to x to lie closer to one another than ψ(X), which is not possible.  This establishes (i).  For (ii), if any face of G(X) is not convex, then the non-convexity is witnessed by one of its vertices having interior angle greater than π. This vertex can then be shifted into the face, contradicting irreducibility.  Finally, (iii) follows from the fact that each face has perimeter at most 2π, which in turn follows from Lemma <ref> by taking A to be the face in question, which is allowed by (ii), and B to be any hemisphere that contains that face.


Suppose X ⊆ S contains a pair of antipodal points x, -x ∈ X.  As long as |X| > 6 so that ψ(X) < π/2, x can be shifted if and only if -x can be shifted as well.  It follows that by shifting antipodal vertices simultaneously, every antipodal X ⊆ S^2 with |X| > 6 and G(X) reducible can be transformed to an antipodal X' ⊆ S^2 with ψ(X') ≥ψ(X) such that either G(X') is irreducible or ψ(X') > ψ(X).  Hence, when classifying the contact graphs of optimal antipodal 2n-packings X, it suffices to consider irreducible candidates for G(X) and possible locations for isolated vertices.  In what follows, we make a few additional observations specific to antipodal sets, including a strengthening of <ref>(iii).

 Suppose X ⊆ S^2 is antipodal with |X| > 6 and G(X) irreducible.
	
		
  (i) The shortest path between an antipodal pair of vertices in G(X) has at least ⌈π/ψ(X)⌉ edges.
		
  (ii) Any two vertices incident to a common face of G(X) are not antipodal.
		
  (iii)  Each face of G(X) has length at most
	
    ⌊2π/ψ(X)sin(π - ψ(X)/2) ⌋.

	



	For (i), note that we need at least π / ψ(X) edges of length ψ(X) to form a path between x, -x ∈ X.  For (ii), suppose that x and -x were incident to a common face of G(X).  This face must then consist of two disjoint paths of at least ⌈π / ψ(X)⌉ edges.  Together with <ref>(iii), we see that each path must form a geodesic of length π and the face must form a spherical lune.  Since ψ(X) < π/2, we could then shift a point distinct from x, -x into the face, demonstrating the reducibility of G(X).  Hence, (ii) follows.  For (iii), observe that if a pair of points x,y ∈ X are not an antipodal pair, then d(x,y) ≤π - ψ(X); otherwise, the distance from x to -y would be less than ψ(X).  Then every face of G(X) lies in a spherical cap of geodesic diameter at most π - ψ(X), which has perimeter at most 2πsin(π - ψ(X)/2).  <ref> then gives an upper bound on the perimeter of each convex face of G(X), from which (iii) follows.


A short argument involving <ref> forbids the isolated vertices of a contact graph from residing within a triangular face.  Additional restrictions for isolated vertices were observed in [Lemma 2, Lemma 9]boroczky03.

 Suppose X ⊆ S^2 and |X| ≥ 13.

	
  (i) Every face of G(X) containing an isolated vertex has length at least 6.
	
  (ii) No hexagonal face of G(X) contains two isolated vertices.





 §.§ Generating Contact Graph Candidates



We now construct a finite set of graphs which contains, up to isomorphism, the contact graph of an optimal antipodal 16-packing.  For the remainder of this section, let Φ denote an optimal projective 8-packing with corresponding optimal antipodal 16-packing X = X(Φ), where the discussion preceding <ref> allows us to insist that G(X) be irreducible without loss of generality.  By <ref>(i), each non-isolated vertex in X has degree 3, 4 or 5 in G(X).  As G(X) provides a double cover for G([Φ]), each non-isolated vertex in G([Φ]) has degree 3, 4 or 5 as well.

By deleting isolated vertices from a graph, we obtain its essential part.  Observe that the essential part of G([Φ]) is necessarily connected, since otherwise, G(X) exhibits a face that is not simply connected, and therefore not a convex spherical polygon, violating <ref>(ii). Moreover, the essential part of G([Φ]) must have at least 6 vertices.  Otherwise, G([Φ]) would have 3 isolated vertices, which would imply that the essential part of G(X) consists of at most 5 antipodal pairs.  By <ref>(i), it would also contain a face of length 6, which is impossible since <ref>(ii) forbids antipodal pairs from lying on the same face.

We use  nauty within SageMath sagemath to generate connected graphs with order 6, 7, or 8, with minimum degree at least 3, and with maximum degree at most 5.  This produces a set Γ_0 of 1007 graphs, one of which must be isomorphic to the essential part of G([Φ]).  To obtain candidates for G(X), we generate all planar double covers for the graphs in Γ_0.  Given a graph G on vertices V with edges E, every double cover G̃ of G can be represented on the vertex set V ×{0,1} in the following way.  For each edge {v,w}∈ E, either include in G̃ the edges {(v,0), (w,0)}, {(v,1), (w,1)} or the edges {(v,0), (w,1)}, {(v,1), (w,0)}.  In this way, a graph with m edges has 2^m double covers, many of which are often isomorphic.

We need only retain double covers that are connected, planar, and consistent with the constraints of this section.  It will be helpful to set ψ_8 := arccos(μ_8); observe that ψ(X) = ψ_8 since X is optimal by assumption.  <ref> provides the bounds

    0.86638 ≤ψ_8 ≤ 0.9273.

The upper bound together with <ref>(i) ensures that every path between antipodal vertices (v,0) and (v,1) must travel along at least 4 edges; we eliminate graphs that violate this condition.  We also eliminate any graphs that contain a wheel subgraph, or a cycle whose vertices all share a common neighbor.  Indeed, since the edge lengths of G(X) are all equal, it could only contain a wheel graph on n + 1 vertices if the interior angle of each triangular face was 2π / n; this is only possible on S^2 for n = 3, 4, 5, but each of these would require an edge length of ψ(X) ≥ 1.107, contradicting (<ref>).  Retaining only one isomorphic copy of each double cover that meets these requirements, we obtain a set Γ_1 of 547 planar graphs that contains an isomorphic copy of the essential part of G(X).

Since we are considering reasonably small graphs,  generates Γ_0 very quickly.  However, extracting Γ_1 from Γ_0 is rather slow due to repeated graph isomorphism queries, taking several hours.  <ref> indicates that we could first trim Γ_0 by insisting that each member be projective planar.  In principle, one could implement this step with the help of an existing fast algorithm for detecting projective planarity [15.8]kocay17.  Unfortunately, no such algorithm is currently available in SageMath, and so our computer-assisted proof does not incorporate this potential speedup.



§ ENUMERATING COMBINATORIAL GRAPH EMBEDDINGS



At this point, we have a sizable collection of contact graph candidates, and later, we will eliminate most of these by geometric considerations.
Before we can do this, we must identify the different planar embeddings of our contact graph candidates up to homeomorphism.
Such equivalence classes can be encoded combinatorially by associating with each vertex a counterclockwise cyclic ordering of its neighbors.
Given a planar embedding of a graph modulo homeomorphism, we call this encoding the corresponding combinatorial embedding of that graph.
In this section, we describe a practical method to enumerate all combinatorial embeddings of a planar graph, which will allow us to recover all combinatorial embeddings for the graphs of Γ_1.

We first summarize some standard language from graph theory.  A cut vertex for G is a vertex whose deletion disconnects G; we say G is separable if it has a cut vertex or biconnected otherwise.  The blocks of a separable graph G are its maximal biconnected subgraphs.  A separation pair for G is a pair of vertices whose simultaneous deletion disconnects G, and we say G is triconnected if it has no separation pair.

It is well known whitney33 that the triconnected planar graphs are exactly those with a unique combinatorial embedding, up to reversing the cyclic ordering of each vertex simultaneously.  As observed in maclane37, it is “natural to try to reduce a graph G which is not triply connected to triply connected constituents.”  For our purposes, since the combinatorial embeddings of the triconnected constituents are unique modulo reflection, such a decomposition would enable one to enumerate all combinatorial embeddings of the original graph.  Similar decompositions are crucial to efficient graph algorithms such as linear-time planarity tests hopcroft74.

We will encode the triconnected constituents of graphs with the SPQR-tree data structure.  This was introduced in dibattista89 and provides an efficient representation of the set of all combinatorial embeddings of a biconnected planar graph dibattista96.  The main idea is to decompose a biconnected graph G into edge-disjoint subgraphs {G_i}; the combinatorial embeddings of G can be obtained from those of {G_i} by considering their possible relative arrangements.  If any of {G_i} is triconnected, we can report their two combinatorial embeddings.  Otherwise, we can iterate the procedure for each biconnected G_i and for each block of each separable G_i.  An example of such a decomposition is depicted in <ref>.
While the formal definition of SPQR-trees is quite technical, we found the presentation in gutwenger10 to be particularly readable.























 §.§ Generating Combinatorial Embeddings


In <ref>, we generated a set Γ_1 consisting of 547 planar graphs, one of which is guaranteed to be the essential part of G(X) for some optimal antipodal 16-packing X.  We now generate combinatorial embeddings of these graphs, and one of these necessarily corresponds to the planar embedding of G(X) modulo homeomorphism and reflection.  Of the 547 graphs in Γ_1, 263 are triconnected, each producing a unique combinatorial embedding modulo reflection.  The remaining 284 graphs are biconnected, so we can produce reflection representatives of their combinatorial embeddings by constructing their SPQR-trees.  All together, we generate a set ℰ_0 of 2540 combinatorial embeddings, but only a fraction of these exhibit the necessary facial structure identified in the previous section.
In particular, by <ref>(iii), all faces must have length at most six, and by <ref>, there must be a hexagonal face for every isolated vertex.
The members of ℰ_0 meeting these requirements form a set ℰ_1 consisting of only 217 combinatorial embeddings.  Five of these have order 14, while the remaining 212 have order 16.  One of these must represent the combinatorial embedding of G(X).



§ GEOMETRIC CONSTRAINTS


At this point, we have identified a set ℰ_1 comprised of 217 combinatorial embeddings of contact graph candidates, and at least one of these corresponds to an optimal antipodal 16-packing.
In this section, we leverage geometric constraints to eliminate almost all of these embeddings from consideration.
In principle, determining ψ_8 amounts to maximizing ψ such that an embedding from ℰ_1 can be drawn on S^2 with equilateral edge length ψ.  To help describe the feasibility region of this optimization problem, we take the interior angles {θ_i} of the faces of a given embedding to be the decision variables.
These angles must satisfy certain relationships imposed by spherical geometry, which results in linear and nonlinear constraints.  We will relax each of these constraints to produce a linear program that must be feasible if the embedding corresponds to an optimal antipodal 16-packing.
Note that this same approach was taken in musin12, musin15 to solve instances of the Tammes problem.

A fundamental tool here is the spherical law of cosines.


If a,b,c are geodesic side lengths for a spherical triangle with interior angle θ opposite a, then

    cos(a) = cos(b)cos(c) + sin(b)sin(c)cos(θ).



Let α_n(d) denote the interior angle of a spherical regular n-gon of side length d.  We record some straightforward consequences of <ref>.


An equilateral spherical triangle with side length ψ has interior angle

    α_3(ψ) = arccos(cos(ψ)/1 + cos(ψ)).

The opposite interior angles within a spherical rhombus are equal.  If a rhombus has side length ψ, then its adjacent interior angles θ, θ' satisfy

    tan(θ/2)tan(θ'/2) = (ψ).

In general, for an equilateral convex spherical polygon with n sides of length ψ, we can express each of its interior angles θ_1, …, θ_n as a function of ψ and any n - 3 of its interior angles.


We have the following additional restrictions on the interior angles of a contact graph.


Let X ⊆ S^2 and set ψ = ψ(X).

	
  (i) Every interior angle θ of G(X) satisfies α_3(ψ) ≤θ≤π.
	
  (ii) The sum of the interior angles that meet at a vertex in G(X) equals 2π.
	
  (iii) The sum of the interior angles of a face of length n in G(X) is at most nα_n(ψ).
	
  (iv) Suppose α,β,γ,δ,ϵ are the interior angles of a pentagonal face of G(X) ordered cyclically.  Then α≥ϵ only if δ≥β.
	
  (v) Suppose α,β,γ,δ,ϵ,ζ are the interior angles of a hexagonal face of G(X) ordered cyclically. Then α≥δ only if both γ≥ζ and ϵ≥β.









For (i), note that every interior angle θ must be at least α(ψ) since vertices of G(X) must be of geodesic distance ψ apart, while the convexity of the faces guarantees θ≤π.  Next, (ii) is obvious.  Since the area of a spherical polygon is determined by the sum of its interior angles, (iii) follows from a version of the isoperimetric inequality that compares the face's area to that of a regular spherical n-gon of the same perimeter.  The proofs of (iv) and (v) can be found in [Proposition 5]danzer86 and [Proposition 12]boroczky03, respectively.


Antipodal packings satisfy additional constraints that do not appear in the relevant Tammes literature.


Let X ⊆ S^2 be antipodal and set ψ = ψ(X).  Suppose G(X) contains a path from x∈ X to -x along four edges e_1, e_2, e_3, e_4.  Denote the angle between e_1 and e_2 by α, the angle between e_2 and e_3 by β, and the angle between e_3 and e_4 by γ.  Then

    cos(α) + cos(γ) = -2^2(ψ)      and    cos(β) ≤ 1 - 2^2(ψ).




Suppose the path along e_1, e_2, e_3, e_4 passes through the vertices x, x_1, x_2, x_3, -x.  Connecting x to x_2 by a geodesic segment results in a spherical triangle with side lengths ψ, ψ, d(x,x_2), where d(x,x_2) has opposite angle α. Similarly, connecting x_2 to -x by a geodesic segment results in a spherical triangle with side lengths ψ, ψ, d(x_2, -x), where d(x_2,-x) has opposite angle γ.  Moreover, we have d(x,x_2) + d(x_2,-x) = π since x and -x are antipodal.  Applying <ref>,

    cos(α) + cos(γ) = cos(d(x,x_2)) - cos^2(ψ)/sin^2(ψ) + cos(d(x_2,-x)) - cos^2(ψ)/sin^2(ψ) = -2^2(ψ).

Similarly, connecting x_1 to x_3 results in a spherical triangle with side lengths ψ,ψ,d(x_1,x_3), where d(x_1,x_3) has opposite angle β.  Moreover, we know that d(x_1,x_3) ≥π - 2ψ since we are traveling between antipodal points x and -x.  Applying <ref>,

    cos(β) = cos(d(x_1,x_3)) - cos^2(ψ)/sin^2(ψ)≤-cos(2ψ) - cos^2(ψ)/sin^2(ψ) = 1 - 2^2(ψ).





 §.§ Eliminating Combinatorial Embeddings


From <ref>, we have a set ℰ_1 of 217 combinatorial embeddings.  Fix an embedding E in ℰ_1, let {θ_i}_i≥1 denote the interior angles of each of its faces, and set

    m_i := infθ_i      and    	M_i := supθ_i,

where the infimum and supremum are taken over all drawings of E on S^2 with equilateral geodesic edge length ψ_8.  Supposing that such a drawing exists, then for any index i, a conclusion of the form m_i > M_i produces a contradiction, and so we may eliminate E from consideration.

Recall that <ref> implies 0.86638 ≤ψ_8 ≤ 0.9273, and so <ref> gives 1.1668 ≤α_3(ψ_8) ≤ 1.1864.  Moreover, by <ref>(i), each interior angle θ_i satisfies 1.1668 ≤θ_i ≤π.  This provides us with initial lower bounds on each m_i and initial upper bounds on each M_i.  Next, we relax by treating the edge length ψ as another decision variable.
Denote

    m_0 := infψ    and     M_0 := supψ,

where the infimum and supremum are taken over all drawings of the embedding with edge length ψ satisfying 0.86638≤ψ≤0.9273 and with interior angles θ_i∈[m_i, M_i].  Of course, if the embedding corresponds to an optimal contact graph, then M_0 = ψ_8.  We can use these bounds and simple interval arithmetic to relax each of the nonlinear constraints from the previous section to linear constraints and use linear programming to determine improved bounds on m_i and M_i.
For instance, suppose we currently have bounds of the form m_0≥ a and M_0≤ b. Then <ref> and <ref>(i) together imply

    m_i
    ≥α_3(ψ)
    ≥arccos(cos(a)/1+cos(b))

for every i≥1.
In a similar way, the remainder of <ref> and also <ref> can be used to improve estimates of various m_i and M_i from existing estimates of other m_j and M_j.
In addition, every part of <ref> converts existing bounds on m_0 and M_0 into inequality constraints on all of {θ_i}_i≥1, allowing one to improve existing estimates on each m_i and M_i by linear programming.
As one might expect, this process benefits from iteration.  If our estimates of any m_i and M_i ever imply m_i>M_i, then we may eliminate the corresponding embedding from consideration.  For each member of ℰ_1, we follow this general procedure for several iterations, and as a result, only 17 combinatorial embeddings remain.

For the remaining 17 embeddings, we improve our bounds on each m_i and M_i as follows.  Suppose that we know m_0 ≥ a and M_0 ≤ b.  Then we can partition [a,b] into several subintervals [a,b] = [a_1, b_1] ∪⋯∪ [a_k, b_k] and repeat the same procedure as in the previous paragraph under the assumption ψ∈ [a_j,b_j].  Each linear program involves relaxations of several nonlinear constraints on the interior angles θ_i involving ψ, so it comes as no surprise that we have stronger bounds on each θ_i when we restrict ψ to a small subinterval.  For several embeddings, each subinterval [a_j,b_j] results in an infeasible program, allowing us to eliminate that embedding.  In other cases, we instead obtain improved bounds on each m_i and M_i.  After this stage, only seven combinatorial embeddings remain.

We remove a few more embeddings from consideration by simultaneously bisecting feasiblility regions for ψ and each interior angle for a given face.  For a face of length n, this results in 2^n + 1 regions that can be checked for feasibility individually.  Two of our embeddings lead to infeasible programs in all 2^7 regions corresponding to one of their hexagonal faces.  For others, this again leads to improved upper bounds on ψ, and iterating our procedure finally shows that one last embedding leads to an infeasible program.





At this point, we have reduced our search to four combinatorial embeddings, which we denote by G, H_1, H_2, and H_3; see <ref> for an illustration.
One may verify that G corresponds to the numerical packing obtained by Conway, Hardin and Sloane conway96, which we made exact in mixon19.
Each H_ℓ can be obtained from G by deleting appropriate edges or vertices.
Since deleting an edge or isolating a vertex from a contact graph may correspond to an infinitesimal shift of its vertices, we should not expect the coarse approximations in the linear program described so far to be precise enough to rule out any of H_1, H_2, or H_3.
In the following section, we close this gap by passing to a dual optimization program.




§ OPTIMALITY CONDITIONS AND UNIQUENESS



To complete the proof of the main result, we first apply optimization theory to analyze each of the remaining four contact graph candidates (eliminating all but one), and then we leverage ideas from real algebraic geometry to demonstrate uniqueness.
For the first step, we follow musin15 by testing whether each candidate is consistent with the necessary Karush–Kuhn–Tucker optimality conditions [5.5.3]boyd04.
In words, these conditions imply that optimal packings are infinitesimally rigid and therefore must have an associated equilibrium stress matrix, which we denote by (ω_ij).


	Let X = {x_1, …, x_2n} be an optimal antipodal 2n-packing.  Set
	
    J(i) := {j : x_i ↔ x_j  in  G(X)},

and for each j ∈ J(i), let v_ij denote the unit vector at x_i tangent to the geodesic from x_i to x_j. There exist coefficients ω_ij = ω_ji∈𝐑 for 1 ≤ i, j ≤2n with i≠ j such that


  (i) ω_ij≥ 0 for all i,j with i≠ j,

  (ii) ω_ij = 0 whenever x_i and x_j are not adjacent in G(X),

  (iii) ∑_i,j:i≠ jω_ij = 1, and

  (iv) for each i, ∑_j ∈ J(i)ω_ij v_ij = 0.




For any antipodal 2n-packing Z={z_1,…,z_2n}, we may take z_j + n = -z_j for 1 ≤ j ≤ n without loss of generality by reordering if necessary.
As such, the antipodal 2n-packing problem in S^2 is equivalent to the following program, whose decision variables are z_1,…,z_2n∈𝐑^3 and δ∈𝐑:

    [                                         maximize                                                δ;                                       subject to δ - |z_i - z_j|^2 ≤ 0  for all  i,j  with  i≠ j,;                                                                  |z_i|^2 - 1 = 0 for all  i,  and;                                                            z_i + n + z_i = 0  for all  1 ≤ i ≤ n. ]

For the optimal packing X, we denote the square of its minimum Euclidean distance by Δ:=2 - 2cos(ψ(X)).
Then Δ is the maximum value of the above program with maximizer (x_1,…,x_2n,Δ).


For what follows, it is convenient to put z:=(z_1,…,z_2n) and define
	
    f(z,δ)    = δ,
    
    	g_ij(z,δ)    = δ - |z_i - z_j|^2   for  i,j  with  i≠ j,
    
    	h_i(z,δ)    = |z_i|^2 - 1  for  1 ≤ i ≤ 2n,  and
    
    	h'_i(z,δ)    = z_i + n + z_i  for  1 ≤ i ≤ n.

	The Karush–Kuhn–Tucker conditions guarantee coefficients ω_ij, λ_i, and λ_i' with
	
    ∇ f (x,Δ) = ∑_i=1^2n∑_j=1
    j≠ i^2nω_ij∇ g_ij(x,Δ) + ∑_i = 1^2nλ_i ∇ h_i(x,Δ) + ∑_i = 1^n λ_i' ∇ h_i'(x,Δ).

By symmetry, we may replace both ω_ij and ω_ji with their average, thereby ensuring ω_ij=ω_ji.
Dual feasibility ensures ω_ij≥ 0 for all i,j with i≠ j, implying (i).   Complementary slackness ensures ω_ij g_ij(x,Δ) = 0. As such, when j ∉J(i) and j≠ i, then g_ij(x,Δ) ≠ 0, and so ω_ij = 0.  This gives (ii). Notice that (iii) follows by considering the δ-component of the gradient vector in (<ref>), while its z_i-component yields 
	
    0 = ∑_j ∈ J(i) 2ω_ij(x_j - x_i) + 2(λ_i + λ_i') x_i.

To obtain (iv), project the above identity onto the 2-dimensional subspace that is parallel to the tangent plane to S^2 at x_i, and then observe that for j∈ J(i), the vectors x_j-x_i have the same norm, and therefore the same component in the x_i direction since x_i,x_j∈ S^2.




 §.§ Eliminating Subgraphs
  <ref> reports that every optimal antipodal packing X must have a corresponding stress matrix (ω_ij).
Note that we can verify whether a given (ω_ij) is a stress matrix once we have the graph G(X) and the tension directions (v_ij), both of which are determined by X.
In this subsection, we argue that for each combinatorial embedding H_ℓ illustrated in <ref>, then for every choice of tension directions (v_ij) that are consistent with the bounds derived in the previous section, there is no stress matrix (ω_ij) that is consistent with both H_ℓ and (v_ij).

Given H_ℓ, suppose for the sake of contradiction that there exists an optimal antipodal 16-packing X such that the essential part of G(X) belongs to the homeomorphism equivalence class of H_ℓ.
We start by putting <ref>(iv) in a standard form.
For each i, we apply an isometry Q_i that maps the tangent plane of S^2 at x_i to 𝐑^2 in such a way that one of the tension vectors v_ij is mapped to (1,0).
Put u_ij:=Q_i(v_ij).
Then <ref>(iv) is equivalent to ∑_j∈ J(i)ω_iju_ij=0.
Importantly, each u_ij takes the form u_ij=(cos(θ_ij),sin(θ_ij)), where each θ_ij is a sum of interior angles θ_ij=∑_k∈ K_ijθ_k over an index set K_ij that is determined by the combinatorial embedding of H_ℓ.
For each k, the previous section computed bounds on m_k and M_k such that θ_k∈[m_k,M_k], and as we will see, we can apply the relationships in <ref> to further improve these bounds.





Notice that bounds on the interior angles {θ_k} imply bounds on the angles θ_ij=∑_k∈ K_ijθ_k, which in turn imply bounds on the coordinates of u_ij that take the form cos(θ_ij) ∈ [c_ij,C_ij] and sin(θ_ij) ∈ [s_ij,S_ij].  If (ω_ij) satisfies <ref>(iv), then it must be the case that

    ∑_j ∈ J(i)ω_ij c_ij   ≤ 0 ≤∑_j ∈ J(i)ω_ij C_ij, and
    ∑_j ∈ J(i)ω_ij s_ij   ≤ 0 ≤∑_j ∈ J(i)ω_ij S_ij.

By linear programming, we can quickly determine whether the set Ω consisting of all (ω_ij) that satisfy both <ref>(i)–(iii) and (<ref>) is empty.
Overall, we can sharpen the bounds from <ref> by imposing the additional constraint that Ω is nonempty.
With this approach, we repeatedly partition the admissible intervals for each variable from <ref>, improving the bounds on each m_k and M_k until we demonstrate m_k>M_k for some k, thereby delivering the desired contradiction.




 §.§ The Optimal Configuration


By the process of elimination, every maximal connected set of optimal antipodal 16-packings in S^2 has a representative packing X whose contact graph belongs to the homeomorphism equivalence class of either G or its reflection.
Furthermore, no other irreducible graph is the contact graph of an optimal antipodal 16-packing.
Since G is irreducible with no isolated vertices, we may conclude that the contact graph of every optimal antipodal 16-packing is isomorphic to G.

By selecting one vector from each antipodal pair in X, we obtain an optimal projective 8-packing Φ = {φ_1, …, φ_8}.  If we abuse notation by treating Φ as a matrix with columns φ_1, …, φ_8, then its Gram matrix is given by Φ^TΦ.  Notice that (Φ^TΦ)_ij = μ_8 precisely when φ_i and φ_j are adjacent in G(X), while (Φ^TΦ)_ij = -μ_8 precisely when φ_i and -φ_j are adjacent in G(X).  Every other entry of Φ^TΦ is strictly less than μ_8 in absolute value.
We choose φ_1, …, φ_8∈ X so that the resulting Gram matrix takes the form

    Φ^TΦ = ([             1           a_1  lightgrayμ_8           a_2           a_3 lightgray-μ_8           a_4 lightgray-μ_8;           a_1             1 lightgray-μ_8           a_5 lightgray-μ_8           a_6 lightgray-μ_8           a_7;  lightgrayμ_8 lightgray-μ_8             1 lightgray-μ_8           a_8 lightgray-μ_8           a_9          a_10;           a_2           a_5  lightgrayμ_8             1          a_11          a_12  lightgrayμ_8  lightgrayμ_8;           a_3 lightgray-μ_8           a_8          a_11             1  lightgrayμ_8  lightgrayμ_8 lightgray-μ_8; lightgray-μ_8           a_6 lightgray-μ_8          a_12  lightgrayμ_8             1  lightgrayμ_8          a_13;           a_4 lightgray-μ_8           a_9  lightgrayμ_8  lightgrayμ_8  lightgrayμ_8             1          a_14; lightgray-μ_8           a_7          a_10  lightgrayμ_8 lightgray-μ_8          a_13          a_14             1 ]).

Overall, Φ^TΦ is rank-3 and positive semidefinite, a_j satisfies |a_j| < μ_8 for every 1 ≤ j ≤ 14, and μ_8 satisfies <ref>.
In mixon19, we used cylindrical algebraic decomposition to establish that there is a unique real matrix of this form.  Moreover, we computed μ_8 exactly as the third-largest root of

    1 + 5x - 8x^2 - 80x^3 - 78x^4 + 146x^5 - 80x^6 - 584x^7 + 677x^8 + 1537x^9,

which is given numerically by μ_8 ≈ 0.6475889787.
This completes the proof of the main result.



§ FUTURE WORK



In this paper, we identified an optimal projective 8-packing in 𝐑𝐏^2, and we established that it is unique up to isometry.
The unique optimality of this packing was predicted in conway96, where numerical solutions for putatively optimal projective n-packings in 𝐑𝐏^2 were presented for n ≤ 55.
For the n=9 case, they conjecture that there are two distinct isometry classes of optimal projective packings.
It would be natural to extend our work to study these larger packings, and such extensions would benefit from a more computationally efficient procedure.
The primary bottleneck appears in <ref>, where, given a projective contact graph candidate in Γ_0 on m edges, we naively produced all of its 2^m double covers, and then we retained a single isomorphic copy of each double cover that happened to be planar and consistent with other combinatorial constraints.
As indicated by <ref>, we only need to consider members of Γ_0 that are projective planar, but we did not leverage this information in our computation; as such, we suspect that a fast projective planarity algorithm [15.8]kocay17 could speed up our procedure.
Also, does there exist a faster-than-naive algorithm for generating planar double covers of a given projective planar graph?
Finally, it would be interesting to investigate comparable approaches for efficiently packing points in 𝐅𝐏^d-1 for 𝐅∈{𝐑,𝐂} with d ≥ 3 and, more generally, packing higher-dimensional subspaces in 𝐅^d.



§ ACKNOWLEDGMENTS


DGM was partially supported by AFOSR FA9550-18-1-0107, NSF DMS 1829955, and the Simons Institute of the Theory of Computing.



