# Disentangling Hype from Practicality: On Realistically Achieving Quantum
Advantage

Torsten Hoefler  Microsoft CorporationOne Microsoft
WayRedmondWashingtonUSA98052 ETH ZurichUniversitaestsstrasse
6ZurichZurichSwitzerland8092<EMAIL_ADDRESS>,  Thomas Häner
<EMAIL_ADDRESS>and  Matthias Troyer<EMAIL_ADDRESS>Microsoft
CorporationOne Microsoft WayRedmondWashingtonUSA98052

###### Abstract.

Quantum computers offer a new paradigm of computing with the potential to
vastly outperform any imagineable classical computer. This has caused a gold
rush towards new quantum algorithms and hardware. In light of the growing
expectations and hype surrounding quantum computing we ask the question which
are the promising applications to realize quantum advantage. We argue that
small data problems and quantum algorithms with super-quadratic speedups are
essential to make quantum computers useful in practice. With these guidelines
one can separate promising applications for quantum computing from those where
classical solutions should be pursued. While most of the proposed quantum
algorithms and applications do not achieve the necessary speedups to be
considered practical, we already see a huge potential in material science and
chemistry. We expect further applications to be developed based on our
guidelines.

††copyright: none

Operating on fundamentally different principles than conventional computers,
quantum computers promise to solve a variety of important problems that seemed
forever intractable on classical computers. Leveraging the quantum foundations
of nature, the time to solve certain problems on quantum computers grows more
slowly with the size of the problem than on classical computers—this is called
_quantum speedup_. Going beyond quantum supremacy (arute19, ), which was the
demonstration of a quantum computer outperforming a classical one for an
artificial problem, an important question is finding meaningful applications
(of academic or commercial interest) that can realistically be solved faster
on a quantum computer than on a classical one. We call this a practical
quantum advantage, or _quantum practicality_ for short.

Figure 1. Quantum speedup: The time needed to solve certain problems with
quantum algorithms increases more slowly than that of any known classical
algorithm as the problem size N increases. To be practical, however, we need
more than an asymptotic speedup: the crossover time where quantum advantage
gets realized needs to be reasonably short and the crossover problem size not
too large. (For illustration, the time axis is scaled such that the quantum
algorithm is a straight line.)

There is a maze of hard problems that have been suggested to profit from
quantum acceleration: from cryptanalysis, chemistry and materials science, to
optimization, big data, machine learning, database search, drug design and
protein folding, fluid dynamics and weather prediction. But which of these
applications realistically offer a potential quantum advantage in practice?
For this, we cannot only rely on asymptotic speedups but must consider the
constants involved. Being optimistic in our outlook for quantum computers, we
will identify clear guidelines for quantum practicality and use them to
classify which of the many proposed applications for quantum computing show
promise and which ones would require significant algorithmic improvements to
become practically-relevant.

To establish reliable guidelines, or lower bounds for the required speedup of
a quantum computer, we err on the side of being optimistic for quantum and
overly pessimistic for classical computing. Despite our overly-optimistic
assumptions, our analysis will show that a wide range of often-cited
applications is unlikely to result in a practical quantum advantage without
_significant_ algorithmic improvements. We compare the performance of only a
single classical chip that is fabricated today similar to the one used in the
NVIDIA A100 GPU which fits around 54 billion transistors (nvidia2020, ) with
an optimistic assumption for a hypothetical quantum computer that may be
available in the next decades with 10,000 error-corrected logical qubits, 10
$\mu s$ gate time for logical operations, the ability to simultaneously
perform gate operations on all qubits and all-to-all connectivity for fault
tolerant two-qubit gates.111Note that no quantum error correction scheme
exists today that allows simultaneous execution of gates and all-to-all
connectivity without at least a $O(\sqrt{N})$ slowdown for $N$ qubits.

| GPU | ASIC | Future Quantum
---|---|---|---
I/O bandwidth | 10,000 Gbit/s | 10,000 Gbit/s | 1 Gbit/s
Operation throughput |  |  |
16-bit floating point | 195 Top/s | 550 Top/s | 10.5 kop/s
32-bit integer | 9.75 Top/s | 215 Top/s | 0.83 kop/s
binary (boolean logical) | 4,992 Top/s | 77,000 Top/s | 235 kop/s
Table 1. Performance comparison. We compare the peak performance of a single
classical chip that can be manufactured today (similar to an NVIDIA A100 GPU,
or an ASIC with a similar number of transistors) with a future quantum
computer with 10,000 error-corrected logical qubits, 10$\mu s$ gate time for
logical operations and all-to-all connectivity. We consider an estimate of the
I/O bandwidth (namely the number of operations per second) and three types of
operations: logical binary operations, 16-bit floating point, 32-bit integer
or fixed-point arithmetic multiply add operations.

#### I/O bandwidth

We first consider the fundamental I/O bottleneck that limits quantum computers
in their interaction with the classical world, which determines bounds for
data input and output bandwidths. Scalable implementations of quantum random
access memory (QRAM (Giovannetti_2008, ; Giovannetti_20082, )) demand a fault-
tolerant error corrected implementation and the bandwidth is then
fundamentally limited by the number of quantum gate operations or measurements
that can be performed per unit time. We assume only a single gate operation
per input bit. For our optimistic future quantum computer the resulting rate
is 10,000 times smaller than for an existing classical chip (see Table 1). We
immediately see that any problem that is limited by accessing classical data,
such as search problems in databases, will be solved faster by classical
computers. Similarly, a potentially exponential quantum speedup in linear
algebra problems (hhl, ), vanishes when the matrix has to be loaded from
classical data, or when the full solution vector should be read out. More
generally, quantum computers will be practical for _“big compute” problems on
small data_ , not big data problems.

#### Crossover scale

With quantum speedup, asymptotically fewer operations will be needed on a
quantum computer than on a classical computer. Due to the high operational
complexity and slower gate operations, however, each operation on a quantum
computer will be slower than a corresponding classical one. As sketched in
Figure 1, classical computers will thus always be faster for small problems
and quantum advantage is realized beyond a problem-dependent crossover scale
where the gain due to quantum speedup overcomes the constant slowdown of the
quantum computer. To have real practical impact, the crossover time needs to
be short, not more than weeks. Constants matter in determining the utility for
applications, as with any runtime estimate in computing.

#### Compute performance

To model performance, we employ the well-known work-depth model from classical
parallel computing to determine upper bounds of classical silicon-based
computations and an extension for quantum computations. In this model, the
work is the total number of operations and applies to both classical and
quantum executions. In Table 1 we provide concrete examples using three types
of operations: logical operations, 16-bit floating point, and 32-bit integer
or fixed-point arithmetic operations for numerical modeling. For the quantum
costs, we consider only the most expensive parts in our estimates, again
benefiting quantum computers: For arithmetic, we count just the dominant cost
of multiplications, assuming that additions are free. Furthermore, for
floating point multiplication, we consider only the cost of the multiplication
of the mantissa (10 bits in fp16). We ignore all further overheads incurred by
the quantum algorithm due to reversible computations, as well as the
significant cost of mapping to a specific hardware architecture with limited
qubit connectivity.

#### Crossover times for classical and quantum computation

To estimate lower bounds for the crossover times, we next consider that while
both classical and quantum computers have to evaluate the same functions
(usually called oracles) that describe a problem, quantum computers require
fewer evaluations thereof due to quantum speedup. At the root of many quantum
acceleration proposals lies a quadratic quantum speedup, including the well-
known _Grover algorithm_ (grover1, ; grover2, ). For such an algorithm, a
problem that needs $X$ function calls on a quantum computer requires
quadratically more, namely on the order of $X^{2}$ calls on a classical
computer. To overcome the large constant performance difference between a
quantum computer and a classical computer, which Table 1 shows to be more than
a factor of $10^{10}$, a large number of function calls $X\gg 10^{10}$ is
needed for the quantum speedup to deliver a practical advantage. In Table 2,
we estimate upper bounds for the complexity of the function that will lead to
a cross-over time of $10^{6}$ seconds, or roughly two weeks.

| Maximum number of operations for practical
---|---
Operation type |  quadratic speedup |  cubic speedup |  quartic speedup
16-bit floating point |  0.2 |  45,800 |  2,800,000
32-bit integer |  0.003 |  1,630 |  130,000
Binary (logical) |  68 |  12,500,000 |  712,000,000
Table 2. Crossover operation counts for quantum algorithms with quadratic,
cubic, and quartic speedups. We determine the number of operations that can be
afforded per function call (see Figure 1) for a quantum computer to show an
advantage over a classical computer using a quantum algorithm with quadratic,
cubic, and quartic quantum speedup. The number of oracle calls required to
reach the crossover point with a quadratic, cubic, and quartic speedup is
computed using the relative runtimes of a single oracle evaluation, and the
total runtime of $10^{6}$ seconds is then used to compute how many basic
operations can be afforded in each oracle call. Since we make optimistic
assumptions for a future quantum computer, we ignore overheads of reversible
arithmetic for quantum computing and limit the classical computer to a single
chip that can be manufactured today. The actual crossover operation counts
will be significantly smaller. A similar analysis for quantum algorithms with
exponential speedups yields promising operation budgets for all datatypes.

We see that with quadratic speedup even a single floating point or integer
operation leads to crossover times of several months. Furthermore, at most 68
binary logical operations can be afforded to stay within our desired crossover
time of two weeks, which is too low for any non-trivial application. Keeping
in mind that these estimates are pessimistic for classical computation (a
single of today’s classical chips) and overly optimistic for quantum computing
(only considering the multiplication of the mantissa and assuming all-to-all
qubit connectivity), we come to the clear conclusion that quadratic speedups
are insufficient for practical quantum advantage. The numbers look better for
cubic or quartic speedups where thousands or millions of operations may be
feasible, and we hence conclude, similarly to Babbush et al. (Babbush2021, ),
that at least cubic or quartic speedups are required for a practical quantum
advantage, and taking into account .

As a result of our overly-optimistic assumptions in favor of quantum
computing, these conclusions will remain valid even with significant advances
in quantum technology of multiple orders of magnitude.

#### Practical and impractical applications

We can now use the above considerations to discuss several classes of
applications where our fundamental bounds draw a line for quantum
practicality. The most likely problems to allow for a practical quantum
advantage are those with exponential quantum speedup. This includes the
simulation of quantum systems for problems in chemistry, materials science,
and quantum physics, as well as cryptanalysis using Shor’s algorithm (shor, ).
The solution of linear systems of equations for highly structured problems
(hhl, ) also has an exponential speedup, but the I/O limitations discussed
above will limit the practicality and undo this advantage if the matrix has to
be loaded from memory instead if being computed based on limited data or
knowledge of the full solution is required (as opposed to just some limited
information obtained by sampling the solution).

Equally importantly, we identify likely dead ends in the maze of applications.
A large range of problem areas with quadratic quantum speedups, such as many
current machine learning training approaches, accelerating drug design and
protein folding with Grover’s algorithm, speeding up Monte Carlo simulations
through quantum walks, as well as more traditional scientific computing
simulations including the solution of many non-linear systems of equations,
such as fluid dynamics in the turbulent regime, weather, and climate
simulations will not achieve quantum advantage with current quantum algorithms
in the foreseeable future. We also conclude that the identified I/O limits
constrain the performance of quantum computing for big data problems,
unstructured linear systems, and database search based on Grover’s algorithm
such that a speedup is unlikely in those cases. Furthermore, Aaronson et al.
(aaronson21, ) show that the achievable quantum speedup of unstructured black-
box algorithms is limited to $\mathcal{O}(N^{4})$. This implies that any
algorithm achieving higher speedup must exploit structure in the problem it
solves.

These considerations help with separating hype from practicality in the search
for quantum applications and can guide algorithmic developments. Specifically,
our analysis shows that 1) it is necessary for the community to focus on
super-quadratic speedups, ideally exponential speedups and 2) one needs to
carefully consider I/O bottlenecks when deriving algorithms to exploit quantum
computation best. Therefore, _the most promising candidates for quantum
practicality are small-data problems with exponential speedup_. Specific
examples where this is the case are quantum problems in chemistry and
materials science (feynman1981simulating, ), which we identify as the most
promising application. We recommend to use precise requirements models
(resest, ) to get more reliable and realistic (less optimistic) estimates in
cases where our rough guidelines indicate a potential practical quantum
advantage.

## Methods

Here we provide more details for how we obtained the numbers above. We compare
our quantum computer with a single microprocessor chip similar to the one used
in the NVIDIA A100 GPU (nvidia2020, ). The A100 chip is around $850mm^{2}$ in
size and manufactured in TSMC’s 7nm N7 silicon process. A100 shows that such a
chip fits around 54.2 billion transistors and can operator at a cycle time of
around 0.7ns.

### Determining peak operation throughputs

In Table 1, we provide concrete examples using three types of operations:
logical operations, 16-bit floating point, and 32-bit integer arithmetic
operations for numerical modeling. Other datatypes could be modeled using our
methodology as well.

#### Classical NVIDIA A100

According to its datasheet, NVIDIA’s A100 GPU, a SIMT-style von Neumann load
store architecture, delivers 312 tera-operations per second (Top/s) with half
precision floating point (fp16) through tensor cores and 78 Top/s through the
normal processing pipeline. NVIDIA assumes a 50/50 mix of addition and
multiplication operations and thus, we divide the number by two, yielding 195
Top/s fp16 performance. The datasheet states 19.5 Top/s for 32-bit integer
operations, again assuming a 50/50 mix of addition and multiplication, leading
to an effective 9.75 Top/s. The binary tensor core performance is listed as
4,992 Top/s with a limited set of instructions.

#### Classical Special Purpose ASIC

Our main analysis assumes that we build a special-purpose ASIC using a similar
technology. If we were to fill the equivalent chip-space of an A100 with a
specialized circuit, we would use existing execution units, for which the size
is typically measured in gate equivalents (GE). A 16-bit floating point unit
(FPU) with addition and multiplication functions requires approximately 7 kGE,
a 32-bit integer unit requires 18 kGE (mach2020, ), and we assume 50 GE for a
simple binary operation. All units include operand buffer registers and
support a set of programmable instructions. We note that simple addition or
multiplication circuits would be significantly cheaper. If we assume a
transistor-to-gate ratio of 10 (jones2018, ) and that 50% of the total chip
area is used for control logic of a dataflow ASIC with the required buffering,
we can fit $54.2B/(7k\cdot 10\cdot 2)=387k$ fp16 units. Similarly, we can fit
$54.2B/(18k\cdot 10\cdot 2)=151k$ int32, or $54.2B/(50\cdot 10\cdot 2)=54.2M$
bin2 units on our hypothetical chip. Assuming a cycle time of 0.7ns, this
leads to a total operation rate of 0.55 fp16, 0.22 int32, and 77.4 bin Pop/s
for an application-specific ASIC with the A100’s technology and budget. The
ASIC thus leads to a raw speedup between roughly 2 and 15x over a programmable
circuit. Thus, on classical silicon, the performance ranges roughly between
$10^{13}$ and $10^{16}$ op/s for binary, int32, and fp16 types.

#### Hypothetical future quantum computer

To determine the costs of N-bit multiplication on a quantum computer, we
choose the controlled adder from Gidney (Gidney2018halvingcostof, ) and
implement the multiplication using N single-bit controlled adders, each
requiring $2N$ CCZ magic states. These states are produced in so called “magic
state factories” that are implemented on the physical chip. While the
resulting multiplier is entirely sequential, we found that this construction
allows for more units to be placed on one chip than for a low-depth adder
and/or for a tree-like reduction of partial products since (1) the number of
CCZ states is lower (and thus fewer magic state factories are required) and
(2) the number of work-qubits is lower. The resulting multiplier has a CCZ-
depth and count of $2N^{2}$ using $5N-1$ qubits ($2N$ input, $2N-1$ output,
$N$ ancilla for the addition).

To compute the space overhead due to CCZ factories, we first use the analysis
of Gidney and Fowler (Gidney2019efficientmagicstate, ) to compute the number
of physical qubits per factory when aiming for circuits (programs) using
$\approx 10^{8}$ CCZ magic states with physical gate errors of $10^{-3}$. We
approximate the overhead in terms of logical qubits by dividing the physical
space overhead by $2d^{2}$, where we choose the error-correcting code distance
$d=31$ to be the same as the distance used for the second level of
distillation (Gidney2019efficientmagicstate, ). Thus we divide Gidney and
Fowler’s 147,904 physical qubits per factory (for details consult the
ancillary spreadsheet (field B40) of Gidney and Fowler) by $2d^{2}=2\cdot
31^{2}$ and get an equivalent space of 77 logical qubits per factory.

For the multiplier of the 10-bit mantissa of an fp16 floating point number, we
need $2\cdot 10^{2}=200$ CCZ states and $5\cdot 10=50$ qubits. Since each
factory takes 5.5 cycles (Gidney2019efficientmagicstate, ) and we can pipeline
the production of CCZ states, we assume 5.5 factories per multiplication unit
such that multipliers don’t wait for magic state production on average. Thus,
each multipler requires 200 cycles and $5N+5.5\cdot 77=50+5.5\cdot 77=473.5$
qubits. With a total of 10,000 logical qubits, we can implement $21$ 10-bit
multipliers on our hypothetical quantum chip. With 10$\mu s$ cycle time, the
200 cycle latency, we get the final rate of less than
$10^{5}cycle/s/(200cycle/op)\cdot 21=10.5kop/s$. For int32 (N=32), the
calculation is equivalent. For binary, we assume two input and one output
qubit for the (binary) adder (Toffoli gate) which does not need ancillas. The
final results are summarized in Table 1.

### A note on parallelism

We assumed massively parallel execution of the oracle on both the classical
and quantum computer (i.e., oracles with a depth of one). If the oracle does
not admit such parallelization, e.g., if depth = work in the worst case
scenario, then the comparison becomes more favorable towards the quantum
computer. One could model this scenario by allowing the classical computer to
only perform one operation per cycle. With a 2 GHz clock frequency, this would
mean a slowdown of about 100,000 times for fp16 on the GPU. In this _extremely
unrealistic_ algorithmic worst case, the oracle would still have to consist of
only several thousands of fp16 operations with a quadratic speedup. However,
we note that in practice, most oracles have low depth and parallelization
across a single chip is achievable, which is what we assumed in the main text.

### Determining maximum operation counts per oracle call

In Table 2, we list the maximum number of operations of a certain type that
can be run to achieve a quantum speedup within a runtime of $10^{6}$ seconds
(a little more than two weeks). The maximum number of classical operations
that can be performed with a single classical chip in $10^{6}$ seconds would
be: 0.55 fp16, 0.22 int32, and 77.4 bin Zop. Similarly, assuming the rates
from Table 1, for a quantum chip: 7, 4, 2,350 Gop, respectively.

We now assume that all calculations are used in oracle calls on the quantum
computer and we ignore all further costs on the quantum machine. We start by
modeling algorithms that provide polynomial $X^{k}$ speedup, for small
constants $k$. For example, for Grover’s algorithms (grover1, ), $k=2$. It is
clear that quantum computers are asymptotically faster (in the number of
oracle queries) for any $k>1$. However, we are interested to find the oracle
complexity (i.e., the number of operations required to evaluate it) for which
a quantum computer is faster than a classical computer within the time-window
of $10^{6}$ seconds.

Let the number of operations required to evaluate a single oracle call be $M$
and let the number of required invocations be $N$. It takes a classical
computer time $T_{c}=N^{k}\cdot M\cdot t_{c}$, whereas a quantum computer
solves the same problem in time $T_{q}=N\cdot M\cdot t_{q}$, where $t_{c}$ and
$t_{q}$ denote the time to evaluate an operation on a classical and on a
quantum computer, respectively. By demanding that the quantum computer should
solve the problem faster than the classical computer and within $10^{6}$
seconds, we find

$\sqrt[k-1]{\frac{t_{q}}{t_{c}}}\leq N\leq\frac{10^{6}}{t_{q}\cdot M},$

which allows us to compute the maximal number of basic operations per oracle
evaluation such that the quantum computer still achieves a practical speedup:

$M\leq 10^{6}\cdot\sqrt[k-1]{\frac{t_{c}}{t^{k}_{q}}}.$

### Determining I/O bandwidth

We use the I/O bandwidth specified in NVIDIA’s A100 datasheet for our
classical chips. For the quantum computer, we assume that one quantum gate is
required per bit of I/O. Using all 10,000 qubits for reading/writing, this
yields an estimate of the I/O bandwidth $B\approx\frac{10,000}{10^{-5}}=1$
Gbit/s.

### Acknowledgments

We thank Luca Benini for helpful discussions about ASIC and processor design
and related overheads and Wim van Dam and all anonymous reviewers for comments
that improved an earlier draft of this work.

## References

  * [1] Scott Aaronson, Shalev Ben-David, Robin Kothari, Shravas Rao, and Avishay Tal.  Degree vs. approximate degree and quantum implications of huang’s sensitivity theorem.  In Proceedings of the 53rd Annual ACM SIGACT Symposium on Theory of Computing, STOC 2021, page 1330–1342, New York, NY, USA, 2021. Association for Computing Machinery.
  * [2] Frank Arute, Kunal Arya, Ryan Babbush, Dave Bacon, Joseph Bardin, Rami Barends, Rupak Biswas, Sergio Boixo, Fernando Brandao, David Buell, Brian Burkett, Yu Chen, Jimmy Chen, Ben Chiaro, Roberto Collins, William Courtney, Andrew Dunsworth, Edward Farhi, Brooks Foxen, Austin Fowler, Craig Michael Gidney, Marissa Giustina, Rob Graff, Keith Guerin, Steve Habegger, Matthew Harrigan, Michael Hartmann, Alan Ho, Markus Rudolf Hoffmann, Trent Huang, Travis Humble, Sergei Isakov, Evan Jeffrey, Zhang Jiang, Dvir Kafri, Kostyantyn Kechedzhi, Julian Kelly, Paul Klimov, Sergey Knysh, Alexander Korotkov, Fedor Kostritsa, Dave Landhuis, Mike Lindmark, Erik Lucero, Dmitry Lyakh, Salvatore Mandrà, Jarrod Ryan McClean, Matthew McEwen, Anthony Megrant, Xiao Mi, Kristel Michielsen, Masoud Mohseni, Josh Mutus, Ofer Naaman, Matthew Neeley, Charles Neill, Murphy Yuezhen Niu, Eric Ostby, Andre Petukhov, John Platt, Chris Quintana, Eleanor G. Rieffel, Pedram Roushan, Nicholas Rubin, Daniel Sank, Kevin J. Satzinger, Vadim Smelyanskiy, Kevin Jeffery Sung, Matt Trevithick, Amit Vainsencher, Benjamin Villalonga, Ted White, Z. Jamie Yao, Ping Yeh, Adam Zalcman, Hartmut Neven, and John Martinis.  Quantum supremacy using a programmable superconducting processor.  Nature, 574:505–510, 2019.
  * [3] Ryan Babbush, Jarrod R. McClean, Michael Newman, Craig Gidney, Sergio Boixo, and Hartmut Neven.  Focus beyond quadratic speedups for error-corrected quantum advantage.  PRX Quantum, 2:010103, Mar 2021.
  * [4] Michael E. Beverland, Prakash Murali, Matthias Troyer, Krysta M. Svore, Torsten Hoefler, Vadym Kliuchnikov, Guang Hao Low, Mathias Soeken, Aarthi Sundaram, and Alexander Vaschillo.  Assessing requirements to scale to practical quantum advantage, 2022.
  * [5] NVIDIA Corporation.  NVIDIA A100 Tensor Core GPU Architecture, 2020.  https://www.nvidia.com/content/dam/en-zz/Solutions/Data-Center/nvidia-ampere-architecture-whitepaper.pdf.
  * [6] Richard P Feynman.  Simulating physics with computers, 1981.  International Journal of Theoretical Physics, 21(6/7), 1981.
  * [7] Craig Gidney.  Halving the cost of quantum addition.  Quantum, 2:74, June 2018.
  * [8] Craig Gidney and Austin G. Fowler.  Efficient magic state factories with a catalyzed $|CCZ\rangle$ to $2|T\rangle$ transformation.  Quantum, 3:135, April 2019.
  * [9] Vittorio Giovannetti, Seth Lloyd, and Lorenzo Maccone.  Architectures for a quantum random access memory.  Physical Review A, 78(5), nov 2008.
  * [10] Vittorio Giovannetti, Seth Lloyd, and Lorenzo Maccone.  Quantum random access memory.  Physical Review Letters, 100(16), apr 2008.
  * [11] L K Grover.  Quantum mechanics helps in searching for a needle in a haystack.  Physical Review Letters, 79(2), 7 1997.
  * [12] Lov K. Grover.  A fast quantum mechanical algorithm for database search.  In Proceedings of the Twenty-Eighth Annual ACM Symposium on Theory of Computing, STOC ’96, page 212–219, New York, NY, USA, 1996. Association for Computing Machinery.
  * [13] Aram W. Harrow, Avinatan Hassidim, and Seth Lloyd.  Quantum algorithm for linear systems of equations.  Phys. Rev. Lett., 103:150502, Oct 2009.
  * [14] Scott Jones.  GlobalFoundries: 7nm, 5nm and 3nm Logic, current and projected processes, 2018.  https://semiwiki.com/semiconductor-manufacturers/intel/7544-7nm-5nm-and-3nm-logic-current-and-projected-processes/.
  * [15] Stefan Mach, Fabian Schuiki, Florian Zaruba, and Luca Benini.  Fpnew: An open-source multiformat floating-point unit architecture for energy-proportional transprecision computing.  IEEE Transactions on Very Large Scale Integration (VLSI) Systems, 29(4):774–787, 2021.
  * [16] Peter W. Shor.  Polynomial-time algorithms for prime factorization and discrete logarithms on a quantum computer.  SIAM J. Comput., 26(5):1484–1509, oct 1997.


We consider ternary disc packings of the plane, i.e. the
packings using discs of three different radii. Packings in which
each “hole” is bounded by three pairwise tangent discs are called
triangulated. There are 164 pairs $(r,s)$, $1{>}r{>}s$,
allowing triangulated packings by discs of radii 1, $r$ and $s$.  In
this paper, we enhance existing methods of dealing with
maximal-density packings in order to find ternary triangulated
packings which maximize the density among all the packings with the
same disc radii.  We showed for 16 pairs that the density is
maximized by a triangulated ternary packing; for 15 other pairs, we
proved the density to be maximized by a triangulated packing using
only two sizes of discs; for 40 pairs, we found non-triangulated
packings strictly denser than any triangulated one; finally, we
classified the remaining cases where our methods are not

§ DENSITY OF DISC AND SPHERE PACKINGS

Given a finite set $S$ of discs, a packing of the plane by $S$ is a
collection of translated copies of discs from $S$ with disjoint

Given a packing $P$, its density $\delta(P)$ is the
proportion of the plane covered by the discs. More formally,

$$ \delta(P) := \limsup_{n\rightarrow\infty}\frac{\text{area}([-n,n]^2\cap P)}{\text{area}([-n,n]^2)} .$$

Nowadays, the density of disc packings is widely studied in different
contexts.  The worst-case optimal density of packings in triangular
and circular containers is found in [14, 12]. In
computer science, there are various connections between sphere
packings and error-correcting codes [4]. Researchers in
chemical physics used Monte Carlo simulations on 2-disc packings and,
among others, obtained lower bounds on the maximal density of packings
with particular disc sizes [11]. Two other groups of
physicists found lower bounds on maximal densities of packings in
$\RR^3$ with 2 sizes of spheres [26, 31]. Upper bounds on
the density are usually much harder to obtain.

The main problem we are interested in is the following: given a finite
set of ball sizes in $\RR^2$ (or $\RR^3$), find a packing of the plane
(or of the space) maximizing the density.

Disc packings self-assembled from colloidal nanodiscs and
nanorods in [32] (at the top) which very accurately
correspond to triangulated packings (on the bottom).

Answering this question has a few practical applications. Chemists,
for example, are interested in the disc and sphere sizes maximizing
the density in order to eventually design compact materials using
spherical nanoparticles of given
sizes [32, 11, 26].  Figure <ref>
gives an illustration of experimental results from

The first known studies of the densest packings go back to Kepler. Many
advances in this area have been made since then.

§.§ Spheres

In a Kepler manuscript dated by 1611, we find a description of the
“cannonball” packing followed by an assertion that it is a densest
1-sphere packings (i.e. packings by equally sized spheres) of
the three-dimensional Euclidean space. This assertion is widely known
by name of the Kepler conjecture.  The “cannonball” packing, also
called face-centered-cubic (FCC) packing, belongs to a family of
packings formed by stacking layers
of spheres centered in the vertices of a triangular lattice, like it is
shown in Figure <ref>. After placing the first two layers,
at each step, there are two choices of how to place the next
layer. This gives us an uncountable set of packings having the same
density. These packings are called close-packings of equal

The density $\delta(P)$ of packing $P$ of $\RR^3$ by unit spheres
never exceeds the density of a close-packing:

\begin{equation}
\delta(P)\leq\frac{\pi}{3\sqrt{2}}
\,.
\label{eq:kepler}
\end{equation}

The first advancement in a proof of the Kepler conjecture was made by
Gauss who, in 1831, showed that close sphere packings maximize the
density among all possible lattice packings, i.e. those where
the disc centers form a lattice [19]. However, the proof of
the whole conjecture took four centuries to be found. Hilbert included
this conjecture, also named “the sphere packing problem”, in his
famous list of 23 problems published in 1900.

The Kepler conjecture was finally proved in a series of 6 papers
submitted by Hales and Furgeson in 1998 [24, 21].  Their
computer-assisted proof took 8 years to be fully reviewed. In 2003,
Hales founded a project called Flyspeck in order to fully verify his
proof by an automated theorem prover. Flyspeck was completed in 2014
including the proof of the Kepler conjecture in the list of computer
verified proofs [20].

The rough idea of the proof consists of locally redistributing (or
weighting) the density function and showing
inequality (<ref>) for this redistributed (weighted)
density. Lagarias calls this approach
“localization” [29]. All in all, in our work, we use
the same general ideas which are discussed in detail in
Section <ref>.

First step of construction of a 3D close-packing.

2D hexagonal packing.

§.§ Discs

§.§.§ 1-disc packings

The two-dimensional variant of the Kepler conjecture claims the 2D
hexagonal packing on the plane (see Fig. <ref>) to have the
highest density among all planar packings by identical discs.

In 1772, Lagrange proved it to be a densest among lattice packings. The
general result was first shown by Thue in 1910 [33].  His
proof was however considered incomplete, a reliable proof was given by
Fejes-Tóth in 1942 [16].

A packing by a set of discs is called saturated if no more
discs from this set can be added to the packing without intersecting
already placed discs. In our setup, we always assume packings to be
saturated since we are interested in the upper bounds on the density
and adding discs to a packing augments it.

The proof of the two-dimensional Kepler conjecture contains the basics
of the strategy used to prove similar results for packings with
several disc sizes, like binary packings (discussed in the next
section) and ternary packings which are studied in this paper. We thus
find it useful to provide the idea of this proof, following its
version given in [5].

Let $P^*$ denote the hexagonal packing of identical discs of radius
1. Our aim is to show for any saturated two-dimensional packing $P$
using discs if radius 1, that its density does not exceed the density
$\delta^*:=\frac{\pi}{2\sqrt{3}}$ of $P^*$.

First, let us consider the Delaunay
triangulations[See [7], especially chapter 23, for
the definition and properties of Delaunay triangulations.] of disc
centers of $P^*$ and $P$ (see
Fig. <ref>, <ref>).  Notice that in the
triangulation of $P^*$, all triangles are equilateral triangles of
side $2$. We define the density of a packing in a given triangle to be
equal to the proportion of this triangle covered by discs of the
packing. The density in any triangle of the triangulation of $P^*$
equals $\frac{\pi}{2\sqrt{3}}=\delta^*$.

It turns out that the density of any triangle in the triangulation of
$P$ is less or equal to $\delta^*$, as proven
in [5]. This allows us to conclude.

This proof is rather simple due to its “locality”: instead of
showing that the density of the whole packing $P$ is bounded by
$\delta^*$, we show it for each triangle of its triangulation (which
is a stronger assertion). Intuitively, the smaller are the units we
work on, the more “local” the proof is.

The Delaunay triangulation of the hexagonal packing $P^*$.

An example of a packing $P$ with its Delaunay triangulation.

Packings of the plane where, as in the hexagonal one, each “hole” is
bounded by three pairwise tangent discs are called
triangulated. More formally,

A packing is called triangulated if the graph formed by
connecting the centers of every pair of tangent discs is
a triangulation.

Fejes Tóth in [17] called such packings “compact”:
since saturated triangulated packings have no “huge holes”, they
intuitively look the most compact. Moreover, around each disc, its
neighbors form a corona of tangent discs which looks like a locally
“optimal” way to pack. For these reasons, triangulated packings
appear to be the best candidates to maximize the density on the whole

Notice, that for a fixed $n$, there exists only a finite number of
$n$-tuples of disc radii $(r_1,\cdots,r_n)$
s.t. $1{=}r_1{>}\cdots{>}r_n{>}0$ allowing a triangulated packing
where all $n$ disc sizes are present [30].

§.§.§ 2-disc packings

Let us now consider binary packings of the plane. We study the
following question: given two discs of radii $1$ and $r<1$, what is
the maximal density of a packing by copies of these discs? We can
always obtain $\frac{\pi}{2\sqrt{3}}$, the density of the hexagonal
packing, by using only one of these discs which gives as a lower bound
on the maximal density. Florian in [13] derived an
upper bound on the density which is equal to the density in the
triangle formed by 2 small and one big pairwise tangent
discs. [9] gives tighter lower and upper bounds of
maximal density of binary packings of the plane, for all values of

9 triangulated periodic binary packings maximizing the
density among packings with the respective disc sizes.

There are 9 values of $r$ allowing triangulated binary packings where
the both disc sizes are
present [28].
Such packings are shown in Fig. <ref>. Each of the
depicted packings is periodic, i.e. if $P$ is a packing in
question, there are two non co-linear vectors $u$ and $v$, called
periods, such that $P+u=P+v=P$. Notice that in this paper, we always
consider packings of the whole plane, and since the triangulated
packings we show here and below are all periodic, it is enough to
represent their fundamental domain (a parallelogram formed by the
period vectors, marked in black in Fig. <ref>) to
see how the whole plane is packed.

Notice that for each of these values of $r$, there is actually an
infinite number of packings having the same density as the one
depicted in Fig. <ref>.  First, changing a finite
portion of a packings does not affect its density. Moreover, for
$b_1$, $b_3$, and $b_7$, there exist non-periodic triangulated
packings with a different global structure, having the same density as
the ones from Fig. <ref> [28]. For the sake of
simplicity, we choose to depict the periodic ones.

It turns out that for these 9 radii $r$, the density is maximized by a
triangulated binary packing – namely, the ones shown in
Figure <ref> [22, 23, 27, 1].

This result suggests the following conjecture [3].

If a finite set of discs allows a triangulated saturated packing,
then the density of packings by these discs is maximized on a
triangulated packing.

This holds for 1-disc packings and 2-disc (binary)
packings. To study this conjecture, the next step is to verify it for
3-disc (ternary) packings which was the main motivation of
our work.

§ RESULT, PLAN OF THE PAPER

Let us turn to the ternary packings. To begin with, we need to find
the sizes of discs allowing triangulated ternary packings. This
problem was solved in [10]: there are 164 pairs $(r,s)$
featuring triangulated packings with discs of radii $1,r,s$. In this
paper, the triplet of discs with radii associated to each of such
pairs is called a case.

The ternary cases are indexed by positive integers from $1$ to $164$,
like in [10].  To avoid confusion, the binary cases (pairs of
disc radii allowing binary triangulated packings) are denoted by
$b_1,\dots,b_9$ which respectively correspond to the cases 1–9
in [1].

The Connelly conjecture (Conjecture <ref>) is applicable
only to the cases having triangulated saturated packings. This
eliminates 15 cases where no triangulated packing is saturated and
leaves us with 149 cases.

Our main contribution is a classification of 71 cases formulated in
the following theorem:

*  For the 16 following cases: 53, 54, 55, 56, 66, 76, 77, 79, 93, 108,
115, 116, 118, 129, 131, 146, the density is maximized by a
triangulated ternary packing.

*  For the cases 1–15, the density is maximized by
triangulated binary packings. For cases 1–5, it is the triangulated
packing of $b_8$; for case 6 — $b_4$; for cases 7–9 — $b_7$; for
cases 10–16 — $b_9$.

*  For the 40 following cases: 19, 20, 25, 47, 51, 60, 63, 64, 70,
73, 80, 92, 95, 97, 98, 99, 100, 104, 110, 111, 117, 119, 126, 132,
133, 135, 136, 137, 138, 139, 141, 142, 151, 152, 154, 159, 161,
162, 163, 164, there exists a non-triangulated packing denser than
any triangulated one.

The values of radii corresponding to the cases from
Theorem <ref> are given in [10]. The triangulated packings
maximizing the density for the cases from Th. <ref> are
depicted in Fig. <ref>. For Th. <ref>, the binary
triangulated packings which maximize the density are present in
Fig. <ref> while the ternary triangulated packings
are given in Fig. <ref>. Triangulated ternary packings and
non-triangulated binary denser packings for Th. <ref> are given
in Fig. <ref> and in the Appendix <ref>.

All in all, we proved the Connelly conjecture to be false and
classified the 149 cases where it was applicable in several groups: 16
cases for which the conjecture holds (Th. <ref>), 15 cases
where the density is maximized on a triangulated packing using only
two discs out of three (Th. <ref>), 40 (periodic) counter
examples to the initial conjecture (Th. <ref>), and the other
cases where our proof strategy does not work. Figure <ref>
represents each case as a point with coordinates $(r,\frac{s}{r})$ and
its number from [10]. The color of the point and the number
corresponds to the class we assigned to the case.

The 16 triangulated ternary packings proved to maximize the
density (the numbers correspond to the numbering in [10]).

The “map” of the 164 cases with triangulated ternary
packings. The cases where no triangulated packing is saturated are
marked in  grey. The cases with a ternary
triangulated packing proved to maximize the density are marked by

with larger case numbers.  The cases where we proved a
triangulated binary packing to maximize the density are marked by
dark green
The cases with
counter examples are  red
The cases featuring two coronas (find the details in
Section <ref>) are  orange. The
cases with empty polytopes (see Section <ref>) are
blue. The remaining cases are marked in black
(Section <ref>).

Section <ref> is dedicated to the cases where a ternary
triangulated packing is proved to maximize the density. We start by
proving Th. <ref>
(Sections <ref>–<ref>). We explain the
approach used in the similar proof for binary packings from [1]
and how we enhance it to make it work in our context. The first
improvement was the generalization of the code universal to all the
cases (instead of treating them one by one as in [1]). The
second necessary generalization we made was leaving a bunch of
parameters as free variables instead of fixing them arbitrary in the

Our proof, as quite a few recent results in the domain, like [21, 12, 9], is based on computer calculations. The
main details of the implementation are provided in
Section <ref> (the complete code is given in the
supplementary materials). We prove Th. <ref> in
Section <ref> by adjusting the proof of Th. <ref>.

Cases from Th. <ref> are treated in Section <ref>. We
obtain a counter-example for each of these cases by applying the
flip-and-flow method [2] on the triangulated binary
packings with disc radii ratio close to the radii ratios of pairs of
discs of this case.

Section <ref> is dedicated to the remaining cases.
Section <ref> presents the 22 cases where one of the
discs appears with at least two different neighborhoods. These cases
are analogous to the case b5 from Fig. <ref> where
a small disc is either surrounded by 6 small discs or by two small
discs and two big ones. Our proof technique is not sufficient to treat
such cases. Handling them requires a less local approach, like it was
done in [1] for b5.

Section <ref> treats the 52 cases where we did not find a
set of constants satisfying all required inequalities needed in our
proof. Even though after several attempts with higher and higher
precision, we concluded that the existence of valid constants is quite
unlikely, it cannot be rigorously proved for the moment. We thus leave
this as an open problem.

Finally, Section <ref> is dedicated to the 4 cases where the
existence of such set of constants is more probable since we could
find the parameters satisfying the majority of constraints, but a few
of them were still not be satisfied. Whether the density is maximized
in these cases is also an open problem.

§ WHEN A TRIANGULATED PACKINGS IS THE DENSEST

In this section, we prove the first two parts of <ref>.

§.§ Proof strategy of Th. <ref>

We follow almost the same steps of the proof as in [1] where the
same result is proven for binary triangulated packings and
in [8] which treats computationally the “simplest” case
among the ternary triangulated packings (case 53).

From the theoretical point of view, the transition from binary
packings to ternary ones seem to be straightforward. In practice,
however, we have much more cases to treat (149 instead of 9) and for
each of them, the problem is much more complex due to the high number
of local combinatorial configurations in possible packings. This
requires a more refined and sensitive choice of parameters than
in [1].

This section is strongly based on [1]: we use the idea of the
proof and quite a few intermediate results. Thus, for the sake of
simplicity, we preserve the same notations.

Let us describe the theoretical background of the proof which is
common for all cases, the only difference being the choice of the
parameters described in Section <ref>.

We are given 3 discs of radii $1,r$ and $s$, $1{>}r{>}s$ and a ternary
triangulated packing of the plane by copies of these discs conjectured
to maximize the density, let us denote it by $P^*$. Our aim is to
prove that for any other packing $P$ using the same discs, its density
$\delta(P)$ does not exceed the density $\delta^*$ of $P^*$.

The main idea common to all the results about the maximal density of
triangulated packings was called “cell balancing” by
Heppes [23] and it perfectly matches this title.  It
consists of two steps: first we locally “redistribute” the density
among some well-defined cells (triangles of the triangulation
in [1, 27, 23] and a mixture of Delaunay simplices
and Voronoi cells, both encoded in so-called decomposition stars,
in [25])
preserving the global density value. Then we prove that the
redistributed density of any cell of $P$ never exceeds $\delta^*$ .

First, let us define triangulations for packings by several sizes of
discs. The FM-triangulation of a packing was introduced
in [18] (it is a particular case of weighted Delaunay
triangulations [7]). Some of its useful properties are given
in [1] (Section 4). The vertices of the FM-triangulation are the
disc centers. There is an edge between two disc centers iff there is a
point $p\in\RR^2$ and a distance $d>0$ such that $p$ is at distance
$d$ from the both discs and at least $d$ from any other

Let $\Tau$ and $\Tau^*$ respectively denote the FM-triangulations
$P$ and $P^*$. The cells we are interested in are triangles of these
triangulations. Instead of working with densities, we introduce an
additive function $E$, called emptiness, which, for a triangle
$T$ in $\Tau$, is defined by

$$E(T):=area(T)\times\delta^* - area(T\cap P)\,.$$

This function was used in [27] by the name of
“excess”. It was inspired by “surplus area” introduced
in [23] defined as
$area(T) - \frac{area(T\cap P)}{\delta^*}$, identical to emptiness up
to multiplication by $\delta^*$. A similar but more complex function
called “score” is used in the proof of the Kepler
conjecture [29].

The emptiness function reflects how “empty” the triangle is compared
to $\delta^*$. Indeed, $E(T)$ is positive if the density of $T$ is
less than $\delta^*$, negative if $T$ is denser, and equals zero if
$\delta(T)=\delta^*$. We use it rather than the density because of its
additivity: the emptiness of a union of two triangles equals the sum
of their emptiness values. This property does not hold for the

To prove that $\delta\leq\delta^*$, it is enough to show that
$ \sum\limits_{T\in\Tau}E(T)\geq 0 $ [1]. This intuitively means
that $P$ is globally more empty and less dense than $P^*$.

Instead of working directly with the emptiness, we define a so-called
potential which plays the role of density redistribution mentioned
above. We do it since this function constructed explicitly is easier
to manipulate. We will construct a potential $U$ such that for any
triangle $T\in\Tau$, its potential does not exceed its emptiness:

\begin{equation}
\label{ineq:tri}
E(T)\geq U(T)
\end{equation}

and the sum of potentials of all triangles in $\Tau$ is
non negative:

\begin{equation}
\label{ineq:U}
\sum_{T\in\Tau}U(T) \geq 0
\end{equation}

If, for $P^*$, there exists $U$ satisfying (<ref>) and
(<ref>) for any packing $P$, then $P^*$ maximizes the density
among packings using the same disc radii:

$$(\text{\ref{ineq:tri}),(\ref{ineq:U})}\;\Longrightarrow\;\sum_{T\in\Tau}E(T)\geq 0 \; \Longrightarrow\; \delta^*\geq\delta\,.$$

The rest of the proof consists in construction of potential $U$
satisfying both (<ref>) and (<ref>) for any
packing $P$.

§.§ How we choose the potential

To construct the potential, we follow the idea first used by Kennedy
in [27] who introduced the “localizing potential” being
inspired by a statistical mechanics notion of
We define the total potential of a triangle as a sum of potentials of
“smaller” units. It will consist of three vertex potentials
$\Uv{}$ defined in the next section and three edge potentials
$\Ue{}$ (Section <ref>). Emptiness redistribution takes place
through vertex and edge potentials: the sum of vertex potentials
around each vertex in the triangulation will be non negative as well
as the sum of two edge potentials of triangles sharing this
edge. These two conditions guarantee us the “global” inequality

All in all, the next four sections follow Sections 5.1-5.2, 6 and 7
of [1]: we are constructing the potential in a way that the
global inequality $(\ref{ineq:U})$ is satisfied, at the same time
trying to minimize the potential of each triangle.  The last step is
to verify for any triangle that the local inequality
$(\ref{ineq:tri})$ is also satisfied, this part is discussed in
Section <ref>.

§.§.§ Vertex potentials

We denote the first part of the potential by $\Uv{}$, it is composed
of three vertex potentials: if $A$, $B$, and $C$ are the vertices of
triangle $T$ in the FM-triangulation $\Tau$ of packing $P$,

Notice that we use the
same notations as in [1] except that we add a dot above the
vertex potentials and a bar above the edge potentials to differentiate
them more easily.

We seek to choose the vertex potentials in a way to satisfy the
following inequality around any vertex $v\in\Tau$:

\begin{equation}
\label{eq:ver}
\sum_{T\in\Tau |v \in T}\Uv{v}(T) \geq 0\tag{$\bullet$}\,.
\end{equation}

Case 54, its triangulated packing (left) and its 10 tight
triangles (right).

This inequality implies global inequality (<ref>) for vertex

A triangle is called tight if it is formed by three pairwise
tangent discs. It is known that the properties of the FM-triangulation
imply that all triangles of a triangulated packing are tight. For each
set of three distinct disc radii, there are 10 tight triangles, one
for each triplet of discs (see Figure <ref>). Let
$E_{xyz}$ denote the emptiness of the tight triangle formed by discs
of radii $x,y,z$. We denote by $V_{xyz}$ the vertex potential of this
triangle in the vertex corresponding to the $y$-disc. Since we set
$V_{xyz}=V_{zyx}$, there are 18 of them.

Three coronas of packing 54.

Inequalities $(\ref{ineq:tri})$ and $(\ref{ineq:U})$, applied to
$\Tau^*$, imply the following constraints on the vertex potentials of
tight triangles appearing in $\Tau^*$.
For each tight triangle $T$ in $\Tau^*$, $\Uv{}(T)=E(T)$ and for each
vertex $v$ in $\Tau^*$,

\begin{equation}
\label{eq:tri_cor}
\sum\limits_{T\in\Tau^*|v\in T}{\Uv{v}(T)}=0.
\end{equation}

The latest means that
the sum of vertex potentials in $v$ of all triangles surrounding it
(the set of these triangles is called the corona of $v$) equals
zero. Let us illustrate it with the case $54$ whose triangulated
packing is given in Fig. <ref>. There are only
three types of tight triangles present in this packing, so the first
class of constraints produces the following equations:

$$3V_{111}=E_{111},\qquad V_{1s1}+2V_{11s}=E_{1s1},\qquad  V_{1rs}+V_{rs1}+V_{s1r}=E_{1rs}\,. $$

The equations for the three coronas given in Fig. <ref>
are respectively:

$$V_{111}+2V_{11s}+4V_{r1s}=0,\qquad 6V_{1rs}=0,\qquad  V_{1s1}+4V_{1sr}=0\,.$$

We thus have $18$ variables and $6$ equations with only $5$ of them
being independent due to periodicity of the packing: since the sum of
the emptiness of tight triangles of the fundamental domain of the
packing equals zero, there is a linear combination of tight triangle
potentials equal to zero.

We also set the potentials of all other tight triangles equal to their
emptiness. In our example, this gives us $7$ more equations that are
still independent and leaves us with $6$ free variables.

Until now we followed the strategy used in [1] which does not
work for the cases proved here, except the case $53$ proved earlier
in [8]. The problem is that in [1, 8], the
free variables mentioned above are all directly set to zero to
simplify further computations. The first generalization we made was to
keep these $6$ variables as free variables and fix them later, in
order to satisfy all inequalities $(\ref{ineq:tri})$ around the
vertices. The $6$ tight vertex potentials left free are the vertex
potentials of the isosceles (but not equilateral) tight triangles in
the vertex adjacent to the two equal edges. More formally, the free
variables are the following:

$$\{V_{xyx}\;|\;x,y\in\{1,r,s\},x\neq y\}.$$
These variables are always
independent and independent of the equations.

The next step is to choose the vertex potentials of all the other (non
tight) triangles in a way that inequalities (<ref>) and
(<ref>) hold.

As in Section 5.2 of [1], we define them as follows:

where $x,y,z$ are the disc sizes of $T$, $v$ is the vertex
corresponding to the $y$-disc, $\hat{v}$ is the angle of $T$ in $v$,
$\widehat{xyz}$ is the angle in the vertex of the $y$-disc in the
tight triangle with discs of radii $x,y,z$, and $m_y$ is a constant
defined below. The difference between the vertex potential of $T$ and
the vertex potential of the tight triangle with the same discs is
proportional to their angle difference, or “how different triangle
$T$ is from the tight one”. The constants $m_1,m_r,m_s$ reflecting
the “importance” of this angle deviation in vertex potentials should
be chosen carefully, which is explained below.

Given a triangle $T$ with discs of radii $x,y,z$ in $\Tau$, let $T^*$
denote the tight triangle formed by the same discs.

Our aim now is to choose the tight vertex potentials and $m_1,m_r,m_s$
in a way that inequality (<ref>) holds around each vertex.
This is the case if, for each vertex $v$ with disc of radius $q$ in
$\Tau$, the following inequality is satisfied:

\begin{equation}
\label{eq:mq}
m_q \geq \frac{-\sum\limits_{j=1}^k{\Uv{}(T^*_j)}}{\left|2\pi - \sum\limits_{j=1}^k{\widehat{v}(T^*_j)}\right|}\,,
\end{equation}

where $T_1,\dots,T_{k}$ is the corona of $v$ in $\Tau$.

The proof is identical to the one for binary packings given
in [1]. Notice that the only case where the denominator equals
zero is when $T^*_1,\dots,T^*_{k}$ form a corona in the triangulated
packing. The sum of vertex potentials of triangles in such corona
equals zero by (<ref>).

We thus have to choose the tight vertex potentials and $m_1,m_r,m_s$
satisfying these inequalities. As explained in the proof of
Proposition 3 of [1], it is enough to perform an exhaustive
search on a finite number of configurations to assure this
inequality. It thus can be done by a computer.

The $6$ values of potentials $V_{xyz}$ of tight triangles which were
left as free variables, as well as $m_1,m_r,m_s$, are chosen to
satisfy inequality (<ref>). The solutions of this inequality (a
subset of $\RR^9$) are the appropriate combinations of values of
potentials and $m_1,m_r,m_s$. The details of computer implementation
of this part and how we choose the precise values are described in
Section <ref>.

Furthermore, using Proposition 4 from Section 6 of [1], we can
cap each vertex potential with a constant value depending only on the
disc type of the vertex. The proof of this proposition roughly
consists in showing that for any vertex, as soon as the vertex
potential of a triangle in its corona exceeds a certain value, the
remaining vertex potentials will never be “too negative” which
implies inequality (<ref>) in this vertex.

This allows us to diminish vertex potential
while still satisfying its non negativity around each vertex
(<ref>). It was not needed in [8] for case 53 but
turns out to be necessary for all the other cases we consider. Thus,
the vertex potential of triangle $T$ in vertex $v$ corresponding to a
$q$-disc will be rewritten as $Z(\Uv{v}(T))=\min{(\Uv{v}(T), Z_q)}$

$$Z_q:=2\pi\,\left|\min\limits_{T_{xqz}\in \Tau}{\frac{V_{xqz}}{\widehat{xqz}}}\right|.$$

§.§.§ Edge potentials

As in Section 7 of [1], in order to keep the potential lower
than the emptiness in all triangles (to satisfy
inequality (<ref>)), we introduce its second part, $\Ue{}$,
called the edge potential. We can now write down the total potential
$U(T)$ of a triangle $T$:

We need $\Ue{}$ to compensate the vertex potential of “stretched”
triangles: those where one angle is too large. Such triangles feature
high vertex potential and low emptiness. The edge potential of a
triangle is equal to the sum of edge potentials corresponding to its
three edges: if $e,d,f$ are the edges of triangle $T$,

We define the edge potential of $T$ in $e$ as

\begin{equation*}
\Ue{e}(T) :=
\begin{cases}
q_{xy} d_e  & \text{if } |e|>l_{xy} \\
0 & \text{otherwise}
\end{cases}
\end{equation*}

where $x,y$ are the radii of the discs with centers in the endpoins of
edge $e$ and $d_e$ is the signed distance from the center $X$ of the
circumscribed circle of $T$ to $e$ ($d_e$ is positive iff $X$ and $T$
are at the same side with respect to $e$). The choice of the constants
$q_{xy},l_{xy}$ is explained below.

With this definition, the sum of edge potentials of two triangles
sharing this edge is always non negative. This is shown in [1]
(Section 7, Proposition 3). Thus, adding the edge potential keeps the
global inequality (<ref>) valid. Meanwhile, it affects the
local distribution of potentials among neighbour triangles approaching
us to our aim which is to satisfy $E(T)\geq U(T)$ in all triangles.

We pick the pairs of constants $q_{xy}$ and $l_{xy}$ (there are 6 of
them, one for each pair $(x,y)$ of disc radii) in order to compensate
high vertex potential of stretched triangles. To find constants
allowing to do this, we compute the vertex potential and the emptiness
of the most “dangerous” triangles: those with two pairs of tangent
discs (see an example in Figure <ref>). We represent them
as curves in function of the length of the only “flexible” edge.

Case 54, triangles with discs $sss$, $srs$, $s1s$: behavior
of the emptiness (in  dark blue) and the potential
while stretching an edge. Initial vertex potential $\Uv{v}(T)$ is marked
by the dotted red line, capped vertex $Z(\Uv{v}(T))$
potential is the dashed red line, the total potential
$U(T)$ is the bold red line. The
dark green segment designates $d_e$ around the moment it
becomes negative, and the  pink asterisk
indicates the value of

To choose $q_{xy}$ and $l_{xy}$ for given $(x,y)$, we consider three
triplets of discs radii: $x1y, xry,xsy$ and for all of them we trace
the curves as ones shown in Figure <ref>. The aim of edge
potential is not to let the capped vertex potential $Z(\Uv{v}(T))$
(dashed red line) exceed the emptiness $E(T)$
(dark blue line). In all the cases considered here, the capped
vertex potential and the emptiness have at most one intersection
except the leftmost point corresponding to the tight triangle (the
neighborhood of this point is a special case treated in detail in
Section <ref>). This intersection is the side length such
that stretching the triangle even more causes vertex potential to be
greater than emptiness.

Let $l^*_{xqy}$ be the side length where this intersection occurs for
a triangle formed by discs with radii $x,q,y$. Notice that $d_e$ is a
decreasing monotonous function on the side length $|e|$. If
$d_{l^*_{xqy}}<0$ (which is the case for all the proved cases), then
it is enough to set $l_{xy}$ equal to at most $l^*_{xqy}-\alpha$ with
a small enough $\alpha$ (we used $\alpha=10^{-5}$). We thus set

As an illustration, when choosing $q_{ss}$ and $l_{ss}$ for case 54
(see Fig. <ref>), the triangle with discs of radii $s,s,s$
has the leftmost intersection of $E(T)$ and $Z(\Uv{v}(T))$
( on the first graph). This
means, $l_{ss}=l^*_{sss}-\alpha$.

Then we choose the coefficient $q_{xy}$ in a way that
$Z(\Uv{v}(T)) + d_eq_{xy}$ stays below $E$ starting from
$l^*_{xqy}-\alpha$ for all $q=1,r,s$ (which is always possible since $d_e$
is negative on this

These choices guarantee the total potential $U(T)$
(bold red line in Fig. <ref>) to be
below the emptiness all the way from tight to stretched triangles.

§.§ Verifying local inequality (<ref>)

In Section <ref>, we constructed the potential in a
way that the global inequality $(\ref{ineq:U})$ is satisfied. We
should now prove that for any triangle $T$, (<ref>) holds:
$E(T)\geq U(T)$.  This part widely uses computer calculations: we execute
them in SageMath [<https://www.sagemath.org/>].  We first treat
the triangles that are close to the tight ones, so-called
epsilon-tight triangles, in Section <ref>. For all
the other triangles, (<ref>) is verified in
Section <ref>.

§.§.§ Epsilon-tight triangles

As in Section 8 of [1], the definition of the vertex potential
allows us to directly show that it does not exceed the emptiness in
the triangles “close enough” to the tight ones.  We do it by
comparing the differentials of the two functions. We are interested in
this case since tight triangles are the limiting point where the
emptiness is equal to the potential.

Let $\Tau_{\epsilon}$ denote the set of epsilon-tight
triangles, those where each pair of discs is at distance at most
$\epsilon$. We can bound the variation $\Delta E$ of emptiness on

$$\Delta E \geq \sum\limits_{i=1,2,3}\min_{\Tau_\epsilon}{\frac{\partial E}{\partial x_i}\Delta x_i}, $$

where $x_1,x_2,x_3$ are the lengths of the sides.

We consider only the values of $\epsilon$ that are strictly smaller
than the $l_{xy}-$constants so that edge potentials do not play a role
in this case: since the smallest $l_{xy}$ is greater than $\epsilon$,
we can find the lower bound of the variation $\Delta U$ of potential:

$$\Delta U \leq \sum\limits_{i=1,2,3}\max_{\Tau_\epsilon}\frac{\partial U}{\partial x_i}\Delta x_i\,.$$

As we have set the potential equal to the emptiness on tight
triangles, choosing $\epsilon$ such that for $i=1,2,3$,

\begin{equation}
\label{ineq:eps}
\min_{\Tau_\epsilon}{\frac{\partial E}{\partial x_i}\Delta x_i} \geq \max_{\Tau_\epsilon}\frac{\partial U}{\partial x_i}\Delta x_i
\end{equation}

guarantees us that the local inequality (<ref>) is satisfied
by all triangles in $\Tau_\epsilon$.

To find the greatest value of $\epsilon$ satisfying this inequality,
we use a dichotomic search, aiming to maximize $\epsilon$ (indeed,
maximizing the set $\Tau_\epsilon$ we minimize the set of the
remaining triangles that should be treated afterwards). At each step
of the binary search, we compute both values in interval arithmetic
(explained in the next section and Section <ref>) and if the
intervals of these values do not intersect and satisfy
inequality (<ref>), the given $\epsilon$ is valid.

§.§.§ Other triangles

The final step is to verify the local inequality $(\ref{ineq:tri})$ on
the set of non epsilon-tight triangles that could appear in a
triangulation of a packing.  Thanks to the properties of
FM-triangulations of saturated packings, this set is compact. Indeed,
the length of the edge connecting centers of discs of radii $x$ and
$y$ in a triangle is at least $x+y$ and at most $x+y+2s$. The
support circle of a triangle in a packing is the circle tangent
to the three discs with centers in the vertices of the triangle.
Thanks to FM-triangulation properties, the support circle never
intersects other discs of the triangulation. The upper bound on the
side length comes from the fact that the radius of the support
circle is at most s (otherwise the packing would not be saturated).
Thus, if the centers of the two discs are respectively $X$ and $Y$ and
the support circle center is $O$, applying the triangle inequality to
$XOY$, we get $|XY|\leq x+y+2s$.

Working with intervals rather than precise values allows us to
decompose a compact continuum set of cases into a finite one.  Given a
triplet of disc radii, each side of each triangle with these 3 discs
in vertices is bounded. Therefore, to prove inequality $E(T)\geq U(T)$
for all triangles with given disc radii, it is enough to show it for
the triangle with sides represented by intervals depending on these
radii. More precisely, we should verify the following inequality:

$$E(\overline{T}_{xyz}) \geq U(\overline{T}_{xyz}),$$

where $\overline{T}_{xyz}$ is a triangle with discs of radii
$x,y,z$ and sides represented by intervals
$[y+z, y+z+2s],[x+z, x+z+2s],[x+y,x+y+2s]$. If it holds (the returned
value is ), then the inequality holds for all possible
triangles with discs $x,y,z$. If the returned value is ,
this means either that the inequality is false or that the intervals
in question intersect (and thus are incomparable). In this case, we
subdivide initial intervals augmenting
precision. Section <ref> describes in detail how it is done
in practice, using interval arithmetic tools of SageMath.

§.§ Computer implementation

As many proofs of the domain, notably the proof of the Kepler
Conjecture [21], the proofs of the maximal density for
triangulated packings, like ours and those
from [1, 8, 27], essentially rely on computer
calculations. In this section, we discuss the details of computer

The treatment of each case consists of two steps. We first choose all
the values necessary to define the potential: tight vertex potentials
$V_{xqy}$, constants $m_q$ and capping values $Z_q$ from
Section <ref>, the value of $\epsilon$
(Section <ref>), and the constants $l_{xy},q_{xy}$ of the
edge potentials (Section <ref>). We choose them in a way to
satisfy the “global” inequality (<ref>). The second step is to
verify the “local” inequality (<ref>) on all possible
triangles explained in Section <ref>.

§.§.§ Interval arithmetic

We use interval arithmetic in two completely different contexts: to
work with real numbers non representable in computer memory and to
verify inequalities on uncountable but compact sets of values. More
precisely, we use intervals to store the values of radii of discs
which are algebraic numbers obtained as roots of polynomials
in [10] as well as the value of $\pi$. The other situation where
we use intervals is to verify the local inequalities on a compact
continuum set of triangles in Section <ref>.

In interval arithmetic, each value is represented by an interval
containing this value and whose endpoints are exact values finitely
representable in computer memory (floating-point numbers). Performing
functions in interval arithmetic preserves both properties. More
precisely, if $x_1,\dots,x_n$ are intervals, and $f$ is an $n$-ary
function, the interval $f(x_1,\dots,x_n)$ must contain
$f(y_1,\dots,y_n)$ for all
$(y_1,\dots,y_n)\in x_1{\times}\dots{\times} x_n$ and its endpoints are
floating-point numbers.

To verify an inequality on two intervals $x_1 < x_2$, it is enough to
compare the right endpoint of $x_1$ and the left endpoint of
$x_2$. The returned value is  only if each pair of values
from these intervals satisfy the inequality. However, if the result is
, that does not mean that the inequality is false on the
numbers represented by $x_1$ and $x_2$, it might also mean that these
intervals overlap.

We worked with interval arithmetic implemented in
SageMath [6], called Arbitrary Precision Real
Intervals [<https://doc.sagemath.org/html/en/reference/rings_numerical/sage/rings/real_mpfi.html>]. The
intervals endpoints are floating-point numbers, the precision we use
in the majority of cases is the default precision of the library where
the mantissa encoding has 53 bits.

§.§.§ Polyhedra

In Section <ref>, we choose the values of vertex potentials
in tight triangles and constants $m_1,m_r,m_s$ in a way to satisfy all
the necessary constraints. These constraints together define a subset
of $\RR^9$ (where the variables are 6 tight vertex potentials
$V_{1r1},V_{1s1},V_{r1r},V_{rsr},V_{s1s},V_{srs}$ and 3 constants
$m_1,m_r,m_s$). We use the Polyhedra
of SageMath to work with them (it allows us to store the solutions of
a system of linear inequalities as a convex polyhedron).

Even more constraints are added in Section <ref>, since
there should exist a positive value of $\epsilon$ satisfying the
inequalities (<ref>). To guarantee that, we verify if
inequality (<ref>) holds for , in other
words, we make sure that this inequality holds for some non-negative
$\epsilon$. We do it in SageMath: to compute both parts of the
inequality, we use interval arithmetic and calculations of
derivatives. The obtained inequalities are intersected with the
polyhedron calculated above. For all the cases considered in this
section, this intersection is not empty (the cases where it was empty
are discussed in Section <ref>). Then we find the maximal
value of $\epsilon>0$ allowing the intersection not to be empty and
this permits us to fix $\epsilon$.

For all the cases we treat in this section, these constraints together
define a compact polyhedron in $\RR^9$ (where the variables are the 6
tight vertex potentials and $m_1,m_r,m_s$).

After we get a polyhedron of valid values, we are free to choose a
point inside to fix them. Our aim at that step is to minimize
potentials of all triangles in order to satisfy (<ref>). We
thus find the three vertices of the polyhedron minimizing $m_1$, $m_r$
and $m_s$ respectively, compute a linear combination of them (the
weights that worked well in practice were respectively 1,1 and
and take a point between this one and the center of the polyhedron in
order to avoid the approximations problems on the border which are
discussed in the next paragraph. Our method to choose the point
described above is a heuristic.

Implementing construction of polyhedra, we encounter the following
problem: the Polyhedra class does not allow coefficients of
constraints to be intervals, while some of the coefficients of our
inequalities are stored as such due to their dependency of $\pi$ and
disc radii. Polyhedra do not support intervals as a base ring for a
good reason: solutions of a system of linear inequalities with
interval coefficients might not form a convex polyhedron. We choose to
replace the intervals with their centers and work with an
approximation of the actual set of valid values for tight potentials
and $m_1,m_r,m_s$. Our polyhedron is stored in a field of rational
values, since this field is computationally quite efficient.

That means, after choosing a point inside this approximated polyhedra,
we can not know if this point actually satisfies all the
constraints. To make sure it does, we then rigorously verify that all
the inequalities with interval coefficients hold in this point.

§.§ Proof of Th. <ref>

Cases 1-18 are special: they are called large separated
in [10] since they do not contain pairs of adjacent medium and
small discs (see Fig. <ref> for the first 15). For each of
these cases, in addition to ternary triangulated packings, there are
other triangulated packings using only two discs out of three. It
happens because the radii of small and medium discs coincide with the
radii of small discs of two cases among $b_1$-$b_9$. It is thus
possible to assemble packings having the same density as the binary
packings of mentioned cases using only two of three discs. It turns
out that in all these cases, the density of one of the mentioned
binary packings exceeds the density of the ternary one. That means,
for each of cases 1-18, the densest packing among the triangulated
ones is a binary packing corresponding to a case from $b_1$-$b_9$
(Fig. <ref>).

Triangulated ternary packings for cases 1-15, where a
triangulated binary packing maximizes the density. For cases 1–5,
it is the triangulated packing of $b_8$; for case 6 — $b_4$; for
cases 7–9 — $b_7$; for cases 10–16 — $b_9$.

Indeed, each of these ternary packings is formed as a “combination”
of two binary packings one of which is denser than the other. Thus,
the densest of the binary packings will also be denser than its
combination with a less dense packing.

We were able to show that the denser triangulated binary packing
maximizes the density among all packings (not only triangulated ones)
for the cases from 1 to 15 (Fig. <ref>). The proof is almost
the same as in Section <ref>.

Let $i$ be the case number and $P_3$ denote its triangulated ternary
packing. Let $P^*_2$ denote the densest triangulated binary packing
using two discs of case $i$ and let $P_2$ denote the less dense
triangulated binary packing using two discs of case $i$. We already
know that $P^*_2$ is denser than the two others,
$\delta(P^*_2)>\delta(P_3)>\delta(P_2)$. Our aim is to show that
$P^*_2$ maximizes the density among all packings by the discs of case

The only difference with the strategy used for other cases concerns
vertex potentials from Section <ref>. Since $P^*_2$ uses
only two discs out of three, it features only 2 coronas instead of
3. Thus, these 2 coronas together with the 10 equations for tight
triangles, give us at most 11 independent equations instead of 12.

We now need to chose 7 free variables instead of 6. We can pick 6
tight potentials of isosceles triangles as before. There remains to
choose the last free variable. Vertex potentials of equilateral tight
triangles can not be picked because of the equations of type
$V_{xxx}=E_{xxx}$: they are already fixed. The remaining vertex
potentials of isosceles triangles ($V_{xxy},\; x\neq y$) can not be
used since they are dependent of the first 6 free variables and the
equations $2V_{xxy}+V_{xyx}=E_{xyx}$. The only candidates thus are
$V_{1rs},V_{1sr},V_{r1s}$; we add one of them.

For cases 16, 17, and 18, the densest binary packing is $b_5$ which
features two different coronas around the small disc, so our method is
not applicable to them as discussed in Section <ref>.

To summarize, for cases 1-18, among triangulated packings, the density is
maximized by a binary packing, not a ternary one as in the Connelly
conjecture. However, whether this packing maximizes the density among
all packings is still an open question for cases 16, 17 and 18.

§ COUNTER-EXAMPLES: PROOF OF TH. <REF>

Starting to work on the density of ternary saturated triangulated
packings, we believed the Connelly conjecture to hold, i.e. that for
all of the 149 cases, a triangulated packing would maximize the
density. Realization that our proof strategy failed for many of them
made us suspect the conjecture to be false. Knowing that the density
of binary triangulated packings (all of them are given in
Figure <ref>) often exceeds the density of ternary
triangulated packings in question gave us an idea to use them in order
to find counter examples.

The first result we obtained was for case 110 [15]. After
generalization, we ended up with 40 counter examples (19, 20, 25, 47,
51, 60, 63, 64, 70, 73, 80, 92, 95, 97, 98, 99, 100, 104, 110, 111,
117, 119, 126, 132, 133, 135, 136, 137, 138, 139, 141, 142, 151, 152,
154, 159, 161, 162, 163, 164). They are all non triangulated packings
using only two discs out of three which have greater densities than
triangulated packings using all three discs. We obtained each of them
deforming a triangulated binary packing with discs whose size ratio is
close to the one of a pair of discs in the triplet associated to the
case. Tiny deformations do not dramatically lower the density and
these packings are dense enough to outplay the ternary triangulated

Left: a triangulated binary packing of case $b_7$. Middle:
a deformation where the small discs are replaced with the small
discs of case 73. Right: a triangulated periodic packing of case
73, its fundamental domain and description are given
in [10].

Let us explain our method on an example.  Recall that the pairs of
discs allowing binary triangulated packings are denoted by
$b_1,\dots,b_9$ while the triplets with ternary triangulated packings
are indexed by positive integers from $1$ to $164$.  Let us consider
case 73, its triangulated ternary packing is given in
Figure <ref>, on the right. Notice that the radius of the small
disc ($s_{73}\approx 0.263$) of case 73 is close to the radius of the
small disc ($r_{b_7}\approx 0.281$) of case $b_7$. Let us deform the
triangulated binary packing of $b_7$ (Figure <ref>, on the
left) replacing the small disc of $b_7$ by the small disc from
$73$. We choose a deformation which breaks as few contacts between
discs as possible: the one given in the center of Figure <ref>.
Observe that the only broken contact is between the two small discs:
they are not tangent anymore. The density of this new non-triangulated
packing $\delta'\approx 0.9245$ is higher than the density of the
triangulated packing $73$ $\delta_{73}\approx 0.9206$
(Figure <ref>, on the right).

This method is called flip-and-flow [2]. The 40 counter
examples were found by computer search. First, for each case $b_i$, we
find the set of pairs of radii from the cases $1$-$164$ with radii
ratio “close enough” (we choose the distance heuristically) to the
ratio of the discs of $b_i$. Then we deform the triangulated packing
of $b_i$ to obtain packings with the found disc ratios. Our way to
deform packings was chosen in order to minimize the number of broken
contacts between discs since intuitively it is the best way to keep
the density high. Finally, the densities of 40 packings obtained by
our method were higher than the densities of the respective ternary
triangulated packings which leaves us with the counter examples given
in Appendix <ref>.

Our method is not universal: there might be other deformations for
certain cases to obtain even higher density and even more counter
examples. Besides that, there might be other cases with ternary
counter examples (notably, among the cases discussed in
Sections <ref>).

§ OTHER CASES

§.§ 2 coronas

Among the necessary conditions on vertex potentials in tight triangles
given in Section <ref>, we saw that the sum of potentials
in the corona around any vertex of triangulated packing $T^*$ must be
equal to zero. In all the proved cases, each disc has only one
possible corona in $T^*$. It is not always the case, more precisely,
among the cases where $T^*$ is saturated, and for which we did not
find counter examples, there are 22 cases where one of the discs
appears with at least two different coronas in $T^*$: 16, 17, 18, 36,
49, 52, 57, 58, 65, 78, 84, 90, 106, 114, 120, 148, 153, 155, 156,
157, 158, 160. Each of these cases features a supplementary corona
consisting of 6 discs of the same size as the central one. We thus
have to add a supplementary condition $6V_{xxx}=0$, where $x$ is the
radius of the disc with two coronas. This however contradicts the
condition $3V_{xxx}=E_{xxx}$ in all of these cases. Our density
redistribution would need to be less local to solve this problem. In
the context of binary triangulated packings, such a case ($b_5$, see
Figure <ref>) is treated in detail in Section 5.3
of [1].

§.§ Empty polyhedra

In Section <ref>, we construct a polyhedron in $\RR^9$ aiming
to contain all valid values of tight vertex potentials and
$m_1,m_r,m_s$. In this Section, we talk about the 52 cases where the
polyhedron obtained by our computations is empty: 21, 22, 23, 26, 27,
34, 35, 46, 48, 50, 59, 61, 67, 68, 69, 71, 72, 74, 81, 82, 83, 85,
86, 87, 88, 89, 91, 94, 96, 101, 102, 103, 105, 107, 109, 112, 113,
121, 122, 123, 124, 125, 127, 128, 130, 134, 140, 143, 145, 147, 149,

The polyhedron formed by the inequalities from
Section <ref>) and inequality (<ref>) for
$\epsilon=0$, represents the values satisfying (<ref>)
featuring a non-negative valid $\epsilon$. These constraints are
necessary for our proof to be correct. If this polyhedron is empty
there are no valid values of tight potentials and $m_1,m_r,m_s$ and
thus our strategy of proof is not applicable.

Nevertheless, our computations are limited by computer memory which
can represent only certain values. Normally, we avoid this problem by
using interval arithmetic (Section <ref>). However, we can
not apply this solution with polyhedra. First, as mentioned in
Section <ref>, in SageMath, the polyhedra module does not
support the interval field as a base ring. Implementing another way to
represent “interval polyhedra” would be unreasonable due to memory
and time constraints of calculations: the polyhedra are constructed
from thousands of inequalities, and performing computations in
interval field significantly increases time and memory costs. Instead,
we use the ring of rationals to store the inequalities
coefficients. Therefore, the polyhedron we work with is an
approximation of the actual polyhedron and may not contain all the
valid sets of values.

Yet, we believe that the polyhedra in question are probably actually
empty in these cases, so the precision issues are not the principal
obstacle.  All in all, some of the cases from this section might
actually maximize the density but we would need an essentially
different approach to be able to prove it. Looking forward, further
attempts to treat these cases would likely need to use a less local
density distribution.

§.§ The 4 mysterious cases

Triangulated ternary packings of the four mysterious

In the four remaining cases (45, 62, 75, and 144) the polyhedron from
Section <ref> is not empty, like for the cases from the
previous section. Nevertheless, we could not find a point in it to
guarantee the local inequality (<ref>) in all triangles: the
problematic triangles are always those close to one of the tight
ones. Minimizing $m_q$ and the tight potentials is an obvious strategy
to minimize the potentials and eventually satisfy (<ref>) but
the capping constants $Z_q$ also dramatically affect potentials.

Trying to find appropriate values of $V_{xyz}, m_q$ and $z_q$, we
represented them with all the constraints coming from (<ref>)
as a linear optimization problem.  This allowed us to encode
problematic triangles violating (<ref>) as constraints and
add them to the system, one by one, each time one appears during local
verification (Section <ref>), in hope to finally
“converge” to a solution which would satisfy (<ref>) on all
triangles. However, this method failed: no solutions were found.

The fact that we could not choose a set of appropriate constants in
these cases does not prove that they do not exist (due to the
approximation issues already discussed in the previous section as well
as the new ones coming from encoding our constraints into a rational
linear problem). We, however, believe that these cases, as well as
those from the previous section, just can not be treated by our proof
methods. They probably require a less local emptiness redistribution
then the one we use.

[1]
N. Bedaride and T. Fernique.
Density of Binary Disc Packings: The Nine Compact Packings.
Discrete and Computational Geometry, 67:1–24, 2022.

[2]
R. Connelly and S. J. Gortler.
Packing Disks by Flipping and Flowing.
Discret. Comput. Geom., 66:1262–1285, 2021.

[3]
R. Connelly, S. Gortler, E. Solomonides, and M. Yampolskaya.
Circle packings, triangulations, and rigidity, 2018.

[4]
J. Conway and N.J.A. Sloane.
Sphere Packings, Lattices and Groups.
Grundlehren der mathematischen Wissenschaften. Springer New York,

[5]
H. Chang and L. Wang.
A Simple Proof of Thue's Theorem on Circle Packing.
arXiv: Metric Geometry, 2010.

[6]
The Sage Developers.
Sage mathematics software (version 9.0).
<http://www.sagemath.org>, 2020.

[7]
S. L. Devadoss and J. O'Rourke.
Discrete and Computational Geometry.
Princeton University Press, 2011.

[8]
T. Fernique.
A densest ternary circle packing in the plane.
<https://arxiv.org/abs/1912.02297>, 2019.

[9]
T. Fernique.
Density of binary disc packings: Lower and upper bounds.
Experimental Mathematics, pages 1–12, 2022.

[10]
T. Fernique, A. Hashemi, and O. Sizova.
Compact packings of the plane with three sizes of discs.
Discret. Comput. Geom., 66(2):613–635, 2021.

[11]
E. Fayen, A. Jagannathan, G. Foffi, and F. Smallenburg.
Infinite-pressure phase diagram of binary mixtures of (non)additive
hard disks.
The Journal of Chemical Physics, 152(20):204901, 2020.

[12]
S. P. Fekete, P Keldenich, and C. Scheffer.
Packing disks into disks with optimal worst-case density.
Discrete and Computational Geometry, 2022.

[13]
A. Florian.
Ausfüllung der Ebene durch Kreise.
Rendiconti del Circolo Matematico di Palermo, 9:300–312, 1960.

[14]
S. P. Fekete, S. Morr, and C. Scheffer.
Split packing: Packing circles into triangles with optimal worst-case
In Algorithms and Data Structures, pages 373–384. Springer
International Publishing, 2017.

[15]
T. Fernique and D. Pchelina.
Compact packings are not always the densest.
<https://arxiv.org/abs/2104.12458>, 2021.

[16]
L. Fejes Tóth.
Über die dichteste Kugellagerung.
Math. Z., 48:676–684, 1943.

[17]
L. Fejes Tóth.
Compact Packing of Circles.
Studia Sci. Math. Hungar., 19:103–107, 1984.

[18]
L. Fejes Tóth and J. Molnár.
Unterdeckung und Überdeckung der Ebene durch Kreise.
Mathematische Nachrichten, 18:235–243, 1958.

[19]
C. F. Gauss.
Untersuchungen über die Eigenschaften der positiven
ternären quadratischen Formen von Ludwig August Seber.
Göttingische gelehrte Anzeigen, 1831.

[20]
T. C. Hales, M. Adams, G. Bauer, D. T. Dang, J. Harrison, T. L. Hoang,
C. Kaliszyk, V. Magron, S. McLaughlin, T. T. Nguyen, T. Q. Nguyen, T. Nipkow,
S. Obua, J. Pleso, J. Rute, A. Solovyev, A. H. T. Ta, T. N. Tran, D. T.
Trieu, J. Urban, K. K. Vu, and R. Zumkeller.
A formal proof of the Kepler conjecture.
Forum of Mathematics, Pi, 5:e2, 2017.

[21]
T. C. Hales.
A proof of the Kepler conjecture.
Annals of Mathematics, 162(3):1065–1185, 2005.

[22]
A. Heppes.
On the densest packing of discs of radius $1$ and $\sqrt{2}-1$.
Studia Scientiarum Mathematicarum Hungarica, 36:433–454, 2000.

[23]
A. Heppes.
Some densest two-size disc packings in the plane.
Discrete and Computational Geometry, 30:241–262, 2003.

[24]
T. C. Hales and S. P. Ferguson.
The Kepler conjecture.
Discrete Comput. Geom., 36(1):1–269, 2006.

[25]
T. C. Hales and S. P. Ferguson.
A Formulation of the Kepler Conjecture, pages 83–133.
Springer New York, New York, NY, 2011.

[26]
A. B. Hopkins, F. H. Stillinger, and S. Torquato.
Densest binary sphere packings.
Phys. Rev. E, 85:021130, 2012.

[27]
T Kennedy.
A densest compact planar packing with two sizes of discs.
<https://arxiv.org/abs/math/0412418>, 2005.

[28]
T. Kennedy.
Compact packings of the plane with two sizes of discs.
Discret. Comput. Geom., 35(2):255–267, 2006.

[29]
J. Lagarias.
Bounds for local density of sphere packings and the Kepler
Discrete and Computational Geometry, 27:165–193, 2002.

[30]
Miek Messerschmidt.
The number of configurations of radii that can occur in compact
packings of the plane with discs of $n$ sizes is finite.
<https://arxiv.org/abs/2110.15831>, 2021.

[31]
Patrick I. O’Toole and Toby S. Hudson.
New high-density packings of similarly sized binary spheres.
The Journal of Physical Chemistry C, 115(39):19037–19040,

[32]
T. Paik, B. T. Diroll, C. R. Kagan, and C. B. Murray.
Binary and Ternary Superlattices Self-Assembled from Colloidal
Nanodisks and Nanorods.
Journal of the American Chemical Society, 137(20):6662–6669,

[33]
A. Thue.
Über die dichteste Zusammenstellung von kongruenten
Kreisen in einer Ebene.
J. Dybwad, 1910.

§ ALL COUNTER EXAMPLES

In this section, we give all the found counter examples grouped by a
type of binary triangulated packing we used to construct it (the list
of binary triangulated packings is given in
Fig. <ref>). Each counter example is presented as
in Fig. <ref>: to the left, we give the deformed binary packing
using a pair of discs whose radii ratio is close to the one of the
binary triangulated packing; to the right, the triangulated ternary

§.§ Counter examples derived from b1

§.§ Counter examples derived from b3

§.§ Counter examples derived from b4

§.§ Counter examples derived from b5

§.§ Counter examples derived from b6

b619, 63

§.§ Counter examples derived from b7

b760, 64, 70, 73, 80, 95, 98, 99, 104, 133, 137, 139, 142, 152, 154, 159, 161, 163

§.§ Counter examples derived from b8

b892, 97, 100, 126, 132, 136, 138, 162, 164

§.§ Counter examples derived from b9

b920, 25, 135, 141