Can exotic disordered "stealthy" particle configurations tolerate arbitrarily large holes?

TL;DR

This study investigates the size limitations of large holes in disordered stealthy particle configurations, revealing they cannot tolerate arbitrarily large cavities, which influences their physical and thermodynamic properties.

Contribution

The paper provides numerical evidence that disordered stealthy configurations in up to three dimensions have a maximum hole size, indicating a structural rigidity not present in typical disordered systems.

Findings

Disordered stealthy configurations cannot have arbitrarily large holes.

Hole probability has compact support in these systems.

Maximum hole size varies across the first three dimensions.

Abstract

The probability of finding a spherical cavity or "hole" of arbitrarily large size in typical disordered many-particle systems in the infinite-size limit (e.g., equilibrium liquid states) is non-zero. Such "hole" statistics are intimately linked to the physical properties of the system. Disordered "stealthy' many-particle configurations in $d$-dimensional Euclidean space $\mathbb{R}^d$ are exotic amorphous states of matter that lie between a liquid and crystal that prohibit single-scattering events for a range of wave vectors and possess no Bragg peaks [Torquato et al., Phys. Rev. X, 2015, 5, 021020]. In this paper, we provide strong numerical evidence that disordered stealthy configurations across the first three space dimensions cannot tolerate arbitrarily large holes in the infinite-system-size limit, i.e., the hole probability has compact support. This structural "rigidity" property apparently endows disordered stealthy systems with novel thermodynamic and physical properties, including desirable band-gap, optical and transport characteristics. We also determine the maximum hole size that any stealthy system can possess across the first three space dimensions.

Tables

Table 1. Table 1: Maximum dimensionless hole size, R c ​ K subscript 𝑅 𝑐 𝐾 R_{c}K , for various 2D crystalline structures.

Crystal	$R_{c}K$
Square lattice	4.44
Honeycomb crystal	4.19
Triangular lattice	4.19
Kagome crystal	3.63

Table 2. Table 2: Maximum dimensionless hole size, R c ​ K subscript 𝑅 𝑐 𝐾 R_{c}K , for various 3D crystalline structures.

Crystal	$R_{c}K$
Face-centered cubic	5.44
Simple cubic	5.44
Hexagonal close packed	5.13
Mean centered-cuboidal lattice 38	5.03
Body-centered cubic	4.97
Simple Hexagonal	4.80
Diamond	4.71
Pyrochlore crystal 39	4.51

Equations

$$E_V(r)=\exp [-\rho v_1(r)],$$

$$\underset{|\mathbf{k}|\to 0}{\lim }S(\mathbf{k})=0.$$

$$S(\mathbf{k})=0\mbox{ for all $0<|\mathbf{k}|<K$.}$$

$${\Phi }_s(K;{\mathbf{r}}^N)=\frac{N}{2V}\sum _{0<|\mathbf{k}|<K}\tilde{v}(\mathbf{k})S(\mathbf{k}),$$

$$\chi =\frac{M}{d(N-1)}$$

$$\rho \chi =\frac{v_1(K)}{2d(2\pi {)}^d}.$$

$$\mathcal{G}=\frac{2}{d}{\int }_0^{\infty }r{E}_V(r)dr.$$

$$\Phi ({\mathbf{r}}^N)={\Phi }_s(K;{\mathbf{r}}^N)+{\Phi }_{ex}(R_f;{\mathbf{r}}^N),$$

$${\Phi }_{ex}(R_f;{\mathbf{r}}^N)=\sum _{i}F(R_f;r_{ic}),$$

$$F(R_f;r_{ic})=\{\begin{array}{ll}(R_f/r_{ic}-1),& \text{if ric<Rf,}\\ 0,& \text{otherwise,}\end{array}$$

$${\Phi }_h({\mathbf{r}}^N)=\sum _{i<j}{v}_2(r_{ij}),$$

$$v_2(r_{ij})=\{\begin{array}{ll}(1-r_{ij}{)}^2,& \text{if if rij<1}\\ 0,& \text{otherwise,}\end{array}$$

$$g_1(\mathbf{r})=\sum _{j=1}^{\infty }{c}_j\delta (|\mathbf{r}|-r_j),$$

$$g_1(\mathbf{r})\propto \sum _{j=1}^{\infty }\delta (|\mathbf{r}|-j+b),$$

$${\tilde{g}}_1(\mathbf{k})={\int }_{\mathbf{r}}\exp (i\mathbf{k}\cdot \mathbf{r})g_1(\mathbf{r})\propto \sum _{j=1}^{\infty }\frac{(j-b{)}^{d/2}}{k^{d/2-1}}{J}_{d/2-1}[k(j-b)],$$

$${\tilde{g}}_{1}(k)\propto\sum_{j=1}^{\infty}\cos[k(j-b)]\qquad\mbox{$(d=1)$},$$

$${\tilde{g}}_{1}(k)\propto\sum_{j=1}^{\infty}(j-b)J_{0}[k(j-b)]\qquad\mbox{$(d=2)$},$$

$${\tilde{g}}_{1}(k)\propto\sum_{j=1}^{\infty}(j-b)\cos[k(j-b)-\pi/2]\qquad\mbox{$(d=3)$}.$$

$${\tilde{g}}_{1}(k)\propto\sum_{j=1}^{\infty}\sqrt{\frac{j-b}{k}}\cos\left[k(j-b)-\frac{\pi}{4}\right]\quad\mbox{ $(d=2)$}.$$

$$p(t)=\sum _{n=1}^{\infty }{I}_n\exp (-t/T_n),$$

$$E_V(R_c)=\exp [-\rho v_1(R_c)]=\frac{1}{N_{c}N}.$$

$$\frac{v_1(R_c)v_1(K)}{2d\chi (2\pi {)}^d}=\ln (N_{c}N).$$

Full text

Can exotic disordered “stealthy” particle configurations tolerate arbitrarily large holes?

G. Zhang,a F. H. Stillinger,a and S. Torquato∗a,b

**Received Xth XXXXXXXXXX 20XX, Accepted Xth XXXXXXXXX 20XX

First published on the web Xth XXXXXXXXXX 200X**

DOI: 10.1039/b000000x

The probability of finding a spherical cavity or “hole” of arbitrarily large size in typical disordered many-particle systems in the infinite-size limit (e.g., equilibrium liquid states) is non-zero. Such “hole” statistics are intimately linked to the physical properties of the system. Disordered “stealthy’ many-particle configurations in $d$-dimensional Euclidean space $\mathbb{R}^{d}$ are exotic amorphous states of matter that lie between a liquid and crystal that prohibit single-scattering events for a range of wave vectors and possess no Bragg peaks [Torquato et al., Phys. Rev. X, 2015, 5, 021020]. In this paper, we provide strong numerical evidence that disordered stealthy configurations across the first three space dimensions cannot tolerate arbitrarily large holes in the infinite-system-size limit, i.e., the hole probability has compact support. This structural “rigidity” property apparently endows disordered stealthy systems with novel thermodynamic and physical properties, including desirable band-gap, optical and transport characteristics. We also determine the maximum hole size that any stealthy system can possess across the first three space dimensions.

1 Introduction

Statistical-mechanical studies of disordered many-particle systems often focus on quantifying various statistics of particle locations. This includes n-body correlation functions,1, 2, 3, 4 the structure factor,1, 3, 2 nearest-neighbor probability distributions,5, 6 and various statistics of the Voronoi cells.7, 8, 9, 10, 11, 12 However, rather than considering the particles themselves, it has been suggested that the space between the particles may be even more fundamental and contain greater statistical-geometrical information.13, 14 A major focus of this paper is the study of a particular property of the void space between point particles in disordered “stealthy” systems,15, 16, 17, 18, 19 which are disordered many-particle configurations that anomalously suppress large-scale density fluctuations, endowing them with unique physical properties. 20, 21, 22, 23, 24, 25 The specific question that we investigate is whether disordered stealthy systems can contain arbitrarily large holes. Here we define a “hole” as a spherical region of a certain radius that is empty of particle centers. It is noteworthy that this hole statistic plays a central role in the “quantizer” and “covering” problems that arise in discrete geometry.26, 14

Given a general many-particle system in $d$-dimensional Euclidean space $\mathbb{R}^{d}$, can one find arbitrarily large holes? For disordered systems, the answer to this question is often “yes.” Consider the void-exclusion probability function, $E_{V}(r)$, which gives the probability of finding a randomly located spherical cavity of radius $r$ empty of particles.13 If $E_{V}(r)$ is non-zero for an arbitrarily large $r$, then one can find arbitrarily large holes in the infinite system, even if these are very rare events. For example, as explained in Fig. 1, the void-exclusion probability for a Poisson point process (i.e., an ideal gas) at number density $\rho$ is given by 13

[TABLE]

where $v_{1}(r)=\pi^{d/2}r^{d}/\Gamma(1+d/2)$ is the volume of a $d$-dimensional sphere of radius $r$,14 and $\Gamma(x)$ is the gamma function. Although $E_{V}$ decays exponentially as $v_{1}(r)$ increases, it is always positive for any finite $r$. Thus, no matter how large a hole is desired, the rare event of forming such a hole can always be observed in the infinite system. Similarly, $E_{V}(r)$ is found to be positive for arbitrarily large $r$’s for equilibrium hard-sphere fluid systems across dimensions.5 Therefore, they also allow arbitrarily large holes. It is noteworthy that $E_{V}(r)$ can be expanded as a series involving n-body correlation functions.5 Therefore, $E_{V}(r)$ requires many-body correlation information to quantify the probability of hole formation.

Even for many-particle systems in which $E_{V}(r)$ is not exactly known in the large-$r$ limit, there are often strong arguments indicating that holes of arbitrary sizes can occur. For equilibrium systems of particles interacting with some potentials (e.g., Lennard-Jones potential) at some positive temperature $T$, the free energy cost of creating a hole, $\Delta F$, often scales as the hole volume and/or hole surface area, and is therefore finite. Thus, the probability of finding a large hole [roughly $\exp(-\Delta F/T)$] is also nonzero. Moreover, hard-sphere systems in a glassy or crystalline state away from jamming points possess collective motions that can produce arbitrarily large holes in the infinite-system limit.27

Besides the aforementioned many-particle systems with unbounded hole sizes, we also know of several systems in which the hole radii are bounded from above. A simple class of systems whose hole probability must have compact support are perfect crystalline (periodic) many-particle systems. Spheres large enough to encompass entire unit cells always contain particles. Thus, holes of arbitrarily-large radii cannot exist. A simple disordered class is saturated random sequential addition (RSA) sphere packings across dimensions. RSA is a time-dependent packing process, in which congruent hard spheres are randomly and sequentially placed into a system without overlap. In the infinite-time limit, the system becomes saturated, i.e., spheres can no longer be added to the packing, and hence holes must be finite in size. By contrast, RSA packings below the saturation density were found to have positive $E_{V}(r)$ for arbitrarily large $r$,28 and therefore allow for the presence of very large holes.

So far we have seen that although all perfect crystalline many-particle systems prohibit arbitrarily large holes, many disordered many-particle systems allow them. A promising class of amorphous structures that may not tolerate arbitrarily large holes is disordered hyperuniform systems. Such systems have received considerable attention because they anomalously suppress density fluctuations.20, 21, 22, 23, 24, 25 Specifically, if one places a spherical window of radius $R$ into a $d$-dimensional many-particle system and counts the number of particles in the window, then the number variance, $\sigma^{2}(R)$, scales as $R^{d}$ for large $R$ in typical disordered systems. Any system in which $\sigma^{2}(R)$ grows slower than $R^{d}$ is said to be hyperuniform.29 Equivalently, a hyperuniform many-particle system is one which the structure factor $S({\bf k})$ tends to zero as the wavenumber $|{\bf k}|$ tends to zero,29 i.e.,

[TABLE]

Disordered hyperuniform systems are a good starting point to search for more examples of disordered systems with bounded hole sizes because the formation of large holes might be inconsistent with hyperuniformity, which suppresses large-scale density fluctuations.

However, we know that not all disordered hyperuniform systems prohibit arbitrarily large holes. For example, in a hyperuniform fermionic-point process in $d$ spatial dimensions, $E_{V}(r)$ scales as $\exp(-cr^{d+1})$ (where $c$ is a constant) for large $r$.30 Also, the hyperuniform two-dimensional one-component plasma possesses an $E_{V}(r)$ that scales as $\exp(-cr^{4})$ for large $r$.31, 32 Both of these systems thus allow arbitrarily large holes. Therefore, hyperuniformity alone is not a sufficient condition to guarantee boundedness of the hole size. Nevertheless, different hyperuniform systems have different levels of suppression for large-scale density fluctuations. While any system in which $\lim_{|\mathbf{k}|\to 0}S(\mathbf{k})=0$ is considered hyperuniform, the “stealthy” variants of hyperuniform systems have $S(\mathbf{k})=0$ in the entire interval $|\mathbf{k}|\in(0,K]$ for a certain value of $K$. Stealthy hyperuniform systems are known to possess many unique physical properties, including negative thermal expansion behavior,20 complete isotropic photonic band gaps comparable in size to those of a photonic crystal,21, 22, 23 transparency even at high densities,24 and nearly optimal transport properties.25 The behavior of $S(\mathbf{k})$ near $\mathbf{k}=\mathbf{0}$ in stealthy systems is identical to that in perfect crystals. Since perfect crystals prohibit large holes, could stealthy hyperuniform systems also prohibit large holes?

In this paper, we present strong numerical evidence that disordered stealthy systems indeed prohibit arbitrarily large holes. It is nontrivial to study the existence of large holes not only because formation of large holes is extremely rare, but also because numerical simulations are limited to finite-sized systems and one wants to infer the infinite-volume-limit behaviors. With periodic boundary conditions, such systems are always perfect crystals, even if the repeating units may be very large. As we have mentioned, perfect crystals always have bounded hole sizes. We developed two numerical techniques to overcome these issues to distinguish whether a system can tolerate arbitrarily large holes or not that can be applied to infer the maximum hole size in general disordered systems (whether they are stealthy or not) in the infinite-volume limit. Specifically, we first attempt to determine the maximum size of the holes that naturally emerges in stealthy hyperuniform systems across the first three space dimensions by studying the tail behavior of $E_{V}(r)$. We find that the tail of $E_{V}(r)$ for stealthy systems is qualitatively similar to that for crystals and saturated RSA sphere packings, which have finite holes, and is qualitatively different from that for Poisson point processes with unbounded hole sizes. We then determine the maximum hole size that any stealthy system can possess across the first three space dimensions. To do this, we generate large stealthy systems with largest possible holes by imposing repulsion fields with sizes equal to the desired hole sizes in stealthy systems. We discover that this method can only create holes of certain finite sizes without breaking stealthiness. In stealthy configurations with largest possible holes, particles concentrate in concentric shells around the hole. Analytical studies on this pattern allows us to derive a conjectured upper bound of the hole radius for all stealthy systems.

The rest of the paper is organized as follows: Section 2 defines stealthy point patterns and two associated parameters, $\chi$ and $K$. Section 3 studies maximum hole sizes and the tail behavior of $E_{V}(r)$ in such systems. Section 4 defines the repulsion field we used to create holes, study the pattern of stealthy systems with such holes, and conjecture an upper bound for the hole radius, in one to three dimensions. Section 5 provides concluding remarks and discussions.

2 Mathematical definitions

For a single-component system with $N$ particles, located at $\mathbf{r}^{N}=\mathbf{r}_{1},\mathbf{r}_{2},\cdots,\mathbf{r}_{N}$, in a simulation box of volume $V$ with periodic boundary conditions in a $d$-dimensional Euclidean space $\mathbb{R}^{d}$, the static structure factor is defined as $S(\mathbf{k})=|\sum_{j=1}^{N}\exp(-i\mathbf{k}\cdot\mathbf{r}_{j})|^{2}/N$, where $i$ is the imaginary unit and $\mathbf{k}$ is a $d$-dimensional wavevector (which must be integer multiples of the reciprocal lattice vectors of the simulation box).3, 33

As we have explained earlier, a hyperuniform system is defined as one in which the number variance $\sigma^{2}(R)$ grows more slowly than $R^{d}$ for large window radius $R$, or a system in which $\lim_{|\mathbf{k}|\to 0}S(\mathbf{k})=0$. Stealthiness is a stronger condition than hyperuniformity. For some positive $K$, we call a system “stealthy up to $K$” if

[TABLE]

In this paper, we define the following potential energy function to be a “stealthy potential” ***This definition of $\Phi_{s}$ actually differs from previous definitions of “stealthy potentials” 18 by a constant, which has no effect on the configurational behavior of the system.

[TABLE]

where ${\tilde{v}}({\bf k})$ is a positive function of $\mathbf{k}$. For present purposes, we choose ${\tilde{v}}({\bf k})=1$ for simplicity. Because $S(\mathbf{k})$ is by definition always non-negative, the ground-state energy of this potential is zero, and the set of the ground states is equal to the set of configurations stealthy up to $K$.

Only half of the constraints in Eq. (3) are independent. This is because by definition, $S(\mathbf{k})=S(-\mathbf{k})$. Let the number of independent constraints be $M$, so the parameter

[TABLE]

quantifies the fraction of degrees of freedom that is constrained. Because $\chi$ is proportional to $M$, it is also proportional to $v_{1}(K)$, the volume of a $d$-dimensional sphere of radius $K$. Indeed, we have previously found 17

[TABLE]

It was found that for $\chi<0.5$, the ground states of stealthy potentials are uncountably infinitely degenerate, and possess no long-range order.18 As $\chi$ increases beyond 0.5, the ground states are still uncountably infinitely degenerate, but develop long-range translational and orientational order.19 As $\chi$ increases further, these ground states eventually undergo phase transitions into the integer lattice, the triangular lattice, and the BCC lattice in one, two, and three dimensions, respectively.17 In this paper, we want to study hole sizes of disordered stealthy systems, and will therefore focus on the $\chi<0.5$ range. Because ground states of the stealthy potentials are uncountably infinitely degenerate, one can have different ways to sample the ground states, which assign different weights to different parts of the ground state manifold. We have previously focused on the zero-temperature limit of the canonical ensemble (i.e., define the probability measure $P({\mathbf{r}}^{N})\propto\exp[-\Phi_{s}(K;{\mathbf{r}}^{N})/k_{B}T]$, where $k_{B}$ is the Boltzmann constant and $T$ is the temperature, and then take the $T\to 0$ limit). However, in this paper, we will also assign different weights to bias toward configurations with large holes.

3 Hole Probability and Maximum Hole Size in Unbiased Stealthy Systems

If an upper bound on the hole sizes exists, how should it depend on $K$ and $\chi$? The $K$ dependence can be easily ascertained from a scaling argument: If there exists a configuration with hole size $R$ that is stealthy up to $K$, then by rescaling the real-space configuration by a factor $\alpha$, one can create another configuration with hole size $R\alpha$, stealthy up to $K/\alpha$. Therefore, the maximum hole radius, $R_{c}$, must be inversely proportional to $K$. Therefore, we henceforth study the dimensionless hole size, $R_{c}K$, rather than $R_{c}$ itself.

A different argument can shed light on the dependence of the hole size on $\chi$. A superposition of multiple configurations, each stealthy up to a certain $K$, is also stealthy up to the same $K$.17 Therefore, if there exist $n$ configurations, each with a hole of radius $R$ that is stealthy up to $K$, then one could superpose them with hole centers aligned to create another configuration with the same hole radius $R$ and $K$. However, since the number of particles increases by a factor of $n$, $\chi$ decreases by a factor of $n$. Therefore, if there exists a configuration of a certain hole size and $K$ at some $\chi$ value, then there exists a configuration of the same hole size and $K$ at arbitrarily small $\chi$ values. In other words, $R_{c}K$ as a function of $\chi$ must achieve the global maximum in the $\chi\to 0^{+}$ limit.

With these preliminary analytical results in mind, let us examine the numerical results from unbiased ground states of stealthy potentials (i.e., $T=0$ limit of the canonical ensemble). We have previously generated such ground states in two and three dimensions for various $\chi$ values by performing low-temperature ($k_{B}T=2\times 10^{-6}$ in 2D and $k_{B}T=10^{-6}$ in 3D) molecular dynamics simulations, periodically taking snapshots, and then minimizing the energy starting from each snapshot; see Ref. 25 for more details. For each $\chi$, we generated 20,000 configurations. The number of particles, $N$, is always between 421 and 751 and is detailed in Ref. 25. For each configuration, we rescaled it to unity $K$ and performed a Voronoi tessellation and found out the largest distance between each Voronoi vertex and its neighbor particles. This distance is the maximum hole size for any particular configuration. We then determined the maximum hole size among all 20,000 configurations and plotted them as a function of $\chi$ in Fig. 2. For a comparison, we also present the same quantity for Poisson point processes at the same conditions, derived in Appendix A. As Eq. (6) shows, with $K$ fixed to unity, $\rho$ is inversely proportional to $\chi$. Thus, it is not surprising that $R_{c}$ for Poisson processes increases as $\chi$ increases. In unbiased stealthy ground states, however, $R_{c}$ weakly increases with increasing $\chi$ and saturates at some constant value, suggesting that $R_{c}$ is bounded for stealthy ground states with fixed $K$. The critical radius $R_{c}$ decreases slightly as $\chi$ tends to zero because unbiased stealthy ground states become less ordered. Therefore, although large hole formation is still possible, its probability decreases. When this probability is too low, it becomes computationally more difficult to find such a large hole with only 20,000 configurations.

Examining the large-$r$ tail behavior of $E_{V}(r)$ suggests strongly that $R_{c}$ is finite in stealthy systems. As we have explained in Sec. I, if the hole size is bounded, $E_{V}(r)$ for some value of $r$ must be identically zero, instead of being exponentially small. In Fig. 3, we closely examine the tails of $E_{V}(r)$ of stealthy systems in the first three space dimensions in a semi-log scale. As we showed earlier, numerically found $R_{c}$ suffer from greater sampling errors if $\chi$ is too small. Thus, to study the tail behavior of $E_{V}(r)$, we choose sufficiently large $\chi$ values (0.45-0.46) in Fig. 3. Nevertheless, we will show in the next section that smaller $\chi$ values do not result in any qualitative difference. For purposes of comparison, we compare our results for stealthy systems to $E_{V}(r)$ for systems in which we know that the holes must be finite in size, namely, lattices in which $E_{V}(r)$ is given exactly14 and saturated RSA sphere packings; and contrast our results to Poisson point processes, in which hole sizes are unbounded. As Fig. 3 shows, the tail behavior of stealthy systems resembles that of crystalline structures and saturated RSA packings. For each of these systems, the logarithm of $E_{V}(r)$ must decay to its bounded cut-off value of $R_{c}$ with an infinite slope at which $E_{V}(R_{c})=0$, which may be regarded to be singularity. However, these figures necessarily present $E_{V}(r)$ above certain positive lower limits and hence only nearly-infinite slopes are apparent. By contrast, Poisson point processes and equilibrium hard-sphere fluids (not shown in the figure), which have unbounded $R_{c}$’s, possess $\log[E_{V}(r)]$’s that comparatively have very small slopes on the scale of the figures, without any singularity. Note that although $E_{V}(r)$ of RSA packings have been studied before,28, 34 this is the first study that focuses on its tail behavior.

It is noteworthy that the three lattice structures we chose (integer, triangular, and BCC lattice) are the optimal solutions of the covering and quantizer problems 26 in their respective dimensions. In a specific dimension and density, the covering problem asks for the configuration with the smallest cutoff in $E_{V}(r)$ (i.e., the smallest $R_{c}$), while the quantizer problem asks for the configuration that minimizes the so-called “quantizer error,” defined as14

[TABLE]

As Fig. 3 shows, in two and three dimensions, $E_{V}(r)$ of stealthy systems at $\chi=0.45-0.46$ is quite close to $E_{V}(r)$ of the triangular and BCC lattices. Therefore, stealthy ground states at high $\chi$ values should provide nearly optimal solutions to these two problems.

4 Stealthy configurations with largest possible holes

In the previous section we studied the largest holes naturally occurring in unbiased disordered ground states of stealthy potentials. In this section, we study the maximum hole sizes consistent with stealthiness. To do so, we impose a radial exclusion field at the center of the simulation box to bias the configuration toward ones with largest holes. We combine the stealthy potential with such an exclusion field, and try to find the ground state of the system. We then study the patterns of the resulting ground states.

4.1 Simulation details

To bias toward configurations with large holes, we let the total potential energy be a sum of the stealthy potential contribution and the exclusion field contribution:

[TABLE]

where $\Phi_{s}(K;\mathbf{r}^{N})$ is the stealthy potential given in Eq. (4), and $\Phi_{ex}(R_{f};\mathbf{r}^{N})$ is the exclusion-field contribution, given by

[TABLE]

where $r_{ic}$ is the radial distance from particle $i$ to the center of the simulation box,

[TABLE]

and $R_{f}$ is the radius of the exclusion field. By varying $R_{f}$, we can probe the largest possible hole size in a particular system. Before $R_{f}$ reaches $R_{c}$ (the upper bound of the hole radius), $\Phi_{ex}$ can be zero. However, once $R_{f}$ surpasses $R_{c}$ for a particular system, $\Phi_{ex}$ must be positive.

If we can find a configuration for which $\Phi(\mathbf{r}^{N})=0$, then both $\Phi_{s}(K;\mathbf{r}^{N})$ and $\Phi_{ex}(R_{f};\mathbf{r}^{N})$ must be zero, and therefore this configuration is stealthy up to $K$ while simultaneously having a hole radius $R_{f}$. To test if there are such configurations, we perform energy minimizations using the L-BFGS algorithm,35, 36, 37 starting from many random initial configurations, and finding if the ending $\Phi(\mathbf{r}^{N})$ in any configuration dropped below a strong tolerance of $10^{-10}$. We consider a certain number, $R_{c}$, to be the numerically found maximum hole size if a zero-energy configuration is found within $N_{trial}$ energy minimization trials for $R_{f}=R_{c}$, but not found for $R_{f}=R_{c}+\delta_{R}$. Here we choose $N_{trial}=100$ and $\delta_{R}=0.01$. For a two-dimensional system at $\chi=0.10$, and $N=400$, with this choice of $N_{trial}$ and $\delta_{R}$ we find $R_{c}K=4.58$; while using $N_{trial}=1000$ and $\delta_{R}=0.0001$, we find $R_{c}K=4.5903$. Therefore, our choice of $N_{trial}$ and $\delta_{R}$ produces $R_{c}K$ values with approximately $10^{-2}$ precision. As explained in our previous work,18 to minimize boundary effects for the stealthy potential, we use a rhombic simulation box with a $60^{\circ}$ interior angle in 2D and a simulation box in the shape of a fundamental cell of a body-centered cubic lattice in 3D with periodic boundary conditions.

As a test for this methodology, we combined the exclusion field [Eq. (9)] with following pair potential

[TABLE]

where

[TABLE]

and performed energy minimizations in two dimensions. For this potential to be zero, any pair of particles cannot be closer than distance 1. Therefore, the ground state of this potential corresponds to an equilibrium hard disk system of diameter 1. As we have mentioned in Sec. I, any such system in the infinite-volume limit must possess an unbounded hole size. Nevertheless, the formation of very large holes is still very rare and may be difficult to observe if one simply samples unbiased configurations. We performed our simulation on an $N=400$ system with volume fraction $\eta=0.5$. As shown in Fig. 4, the energy minimization algorithm is capable of creating a hole of of radius $R=9.2$, although the probability of finding such a hole in an unbiased system is extremely small. According to Eq. (4.21) of Ref. 5, $E_{V}(9.2)=4\times 10^{-279}$. This demonstrates that if the hole size is unbounded in the infinite-system-size limit for some system, this numerical protocol can indeed create very large holes in a finite-size simulation. Figure 4 also shows that in creating such a large hole, the particles are pushed to each other as closely as possible (i.e., up to interparticle contacts). Therefore, even larger holes should be possible if we simulated larger systems at the same volume fraction.

4.2 One-dimensional study

We first examine $R_{c}K$ values found by the above-mentioned algorithm in 1D, since this is computationally the easiest dimension to study and will shed light on corresponding results in higher dimensions. Our result for several different $\chi$’s and system sizes are summarized in Fig. 5. It appears that $R_{c}K$ as a function of $\chi$ is chaotic and displays no systematic trend. Nevertheless, Fig. 5 does show that $R_{c}K$ is always close to $\pi$ but never exceeds it. As we will see later, $\pi$ is the upper bound of $R_{c}K$ in 1D.

Examining stealthy configurations with hole sizes $R_{f}\approx\pi/K$ reveals a more interesting behavior. Such a configuration is shown in Fig. 6. At exclusion-field size $R_{f}=3.1/K$, 100 particles self-assemble into 10 clusters, each containing 10 particles. These clusters then form a one-dimensional integer lattice.

As we have explained in Sec. 3, a superposition of multiple integer lattices, with hole centers aligned, have the same $R_{c}K$ as a single integer lattice. It is straightforward to calculate $R_{c}K$ of an integer lattice: If the distance between neighboring lattice sites is $L$, then the maximum hole radius is $L/2$, and the stealthy range $K$ is equal to the location of the first Bragg peak, $2\pi/L$. Therefore, $R_{c}K$ of any integer lattice is simply $\pi$. To summarize, the numerically found hole radius is never above $\pi/K$; and superposed integer lattices can indeed achieve hole radius $\pi/K$. Therefore, we expect that $\pi/K$ is an upper bound of the hole size for stealthy 1D structure at any $\chi$.

4.3 Two- and Three-dimensional studies

We now move on to study maximum hole sizes in two and three dimensions. As we will see, these higher dimensions are computationally more challenging than 1D because the structures that maximize the hole size is not periodic. The $R_{c}K$ values found by the algorithm mentioned in Sec. 4.1 is presented in Fig. 7. Similar to the 1D case, the dependence of $R_{c}K$ on $\chi$ or $N$ is weak and non-systematic. However, the 2D configurations, one of which is shown in Fig. 8, exhibit a more complicated pattern, in which particles concentrate in a lower-dimensional manifold. Although this pattern is non-crystalline, it is still much more ordered than unbiased stealthy ground states at this $\chi$ value.18 To better reveal this pattern, we computed the one-body correlation function, $g_{1}(\mathbf{r})$, of a 2D system of $\chi=0.1$ and $N=400$, shown in Fig. 9A. The plot shows high-intensity concentric shells around the exclusion field (located at the center of the simulation box) and honeycomb network structures away from the exclusion field. Figure 9B also shows $g_{1}(\mathbf{r})$ of a larger 2D system, which exhibits the same pattern. Figure 9C shows $g_{1}(\mathbf{r})$ of a 3D system, which again has concentric shells around the exclusion field, but the structure away from the center is not obvious.

By pushing $R_{c}K$ to its numerical limit, we obtain periodic structures in 1D but non-periodic structures in 2D and 3D. Is it possible that this transition from periodic structures to non-periodic structures arises from increased numerical difficulties in higher dimensions? To eliminate this possibility, we analytically calculated $R_{c}K$ values for various 2D and 3D periodic structures for comparisons. In 2D, crystal structures achieve $R_{c}K=4.44$ but the system shown in Fig. 9 achieved $R_{c}K=4.6$; while in 3D crystal structures achieve $R_{c}K=5.44$ but the system shown in Fig. 9 achieved $R_{c}K=5.85$. Therefore, these non-periodic structures indeed have the largest known value of $R_{c}K$.

It would be useful to analytically model these $g_{1}(\mathbf{r})$ functions to find the maximum dimensionless hole size in the infinite-system-size limit. We will focus on the rings before considering the honeycomb-like structure away from the hole center. Comparing Fig. 9A with Fig. 9B, we see that increasing $N$ increases the number of rings. Therefore, we expect infinitely many rings in the infinite-system-size limit.

It is instructive to model an isotropic collection of concentric shells, for which we can write

[TABLE]

where $c_{j}$ is the intensity of the shells, $\delta$ is the Dirac delta function, and $r_{j}$ is the location of the shells. To determine $c_{j}$ and $r_{j}$, we computed the angular average of $g_{1}(\mathbf{r})$ shown in Fig. 9B, and identified five peaks from it. As Fig. 10 shows, $r_{j}$ appears linear with $j$, for which linear regression produces $r_{j}=0.0612j-0.01478$. By rescaling the configuration, we can eliminate one fitting parameter and get $r_{j}=j-b$, where $b=0.242$.

To find $c_{j}$, we have computed the fraction of particles located on each ring, $p_{j}$. We find again $p_{j}$ is linear with $j$, with linear regression result $p_{j}=0.0275(j-0.242)\propto r_{j}$. Because $p_{j}$ is proportional to $r_{j}$, and is therefore proportional to the circumference of the rings, each ring has the same intensity. Neglecting a constant factor, we can then set $c_{j}=1$.

To summarize, numerical results suggest that in the infinite-system-size limit,

[TABLE]

where constant $b$ is numerically measured as 0.242 in 2D. Note that this equation also applies to the 1D numerical result (an integer lattice of particle clusters) if we let $b=1/2$. The hole radius of this system is simply $R_{c}=1-b$, the radius of the first ring. After determining $R_{c}$, we should then ascertain $K$. Since $S(\mathbf{k})$ is zero for all $\mathbf{k}$ such that $0<|\mathbf{k}|<K$, the collective coordinates ${\tilde{\rho}}(\mathbf{k})=\sum_{j=1}^{N}\exp(-i\mathbf{k}\cdot\mathbf{r}_{j})$ should also be zero. Thus, the Fourier transform of $g_{1}(\mathbf{r})$, which we denote by ${\tilde{g}}_{1}(\mathbf{k})$, should also be zero in this range. Fourier transforming Eq. (14) gives

[TABLE]

where $J_{\nu}$ is the Bessel function of order $\nu$. In Eq. (15), letting $d=1$, 2, and 3 respectively yields

[TABLE]

[TABLE]

and

[TABLE]

For large $x$, $J_{0}(x)$ is asymptotically $x^{-1/2}\cos(x-\pi/4)$. Substituting this into Eq. (17) gives

[TABLE]

We have already seen in the previous section that, the solution to maximizing $R_{c}K=(1-b)K$ in 1D is $b=1/2$ and $K=2\pi$. Comparing Eq. (19) with Eq. (16), in light of the numerical result b\approx 0.242\mbox{ (d=2)}, suggests that $b=1/4$ in 2D. Somehow the $\pi/4$ phase factor in Eq. (19) changes $b$ to $1/4$. If $K$ is still $2\pi$, then in 2D we have $R_{c}K=(1-b)K=3\pi/2\approx 4.71$, which is slightly above the numerically observed maximum dimensionless hole size $R_{c}K=4.65$. Similarly, in 3D, the $\pi/2$ phase factor in Eq. (18) probably changes $b$ to [math]. If so, the maximum dimensionless hole size in 3D would be $R_{c}K=2\pi$. The difference between $2\pi$ and the numerically observed maximum $R_{c}K=5.86$ is nontrivial, but this can be explained by the increased numerical difficulty in 3D; for example, fewer concentric shells can be formed with the same number of particles in higher dimensions.

5 Conclusions

In this paper, we have investigated the possibility of creating large holes in stealthy hyperuniform many-particle systems using numerical and analytical techniques. We demonstrated that hole sizes in such systems are bounded, first by examining the tail of $E_{V}(r)$ in unbiased ground states of stealthy potentials, and then by imposing radial exclusion fields to bias stealthy configurations toward ones with the largest possible holes. These results suggest that holes larger than a certain upper bound cannot exist in such systems. We then found that $R_{c}K$ is bounded from above by $\pi$, $3\pi/2$, and $2\pi$ in one, two, and three dimensions. A conjectured formula for the upper bound on the dimensionless hole size in $d$ dimensions is $(d+1)\pi/2$. An outstanding problem is a rigorous proof that stealthy infinite systems cannot tolerate holes of arbitrarily large sizes.

Our methods should be applicable to study the existence of arbitrarily large holes in other disordered many-particle systems. This is useful because maximum hole sizes and hole probabilities are related to several other important quantities, including the principal relaxation time $T_{1}$ associated with diffusion-controlled reactions among traps. Specifically, consider a reactive chemical species that can diffuse in the void space between particles, and can be absorbed when it is within a certain distance to any particle. The fraction of such species, released at time $t=0$, that is not absorbed at time $t$ (in other words, the survival probability of the molecules of such species), can be expanded as a series of exponential functions 40

[TABLE]

where $I_{n}$ are coefficients and $T_{n}$ are relaxation times. The largest relaxation time $T_{1}$ is called the “principal relaxation time.” The relaxation times can be measured directly by NMR experiments, in which proton magnetization decays at the phase boundary.41, 42, 43 It has been demonstrated that $T_{1}$ is determined by the largest holes in the configurations, and is therefore divergent if arbitrarily large holes can occur.40 Indeed, for a reactive species in equilibrium hard-sphere systems, the large-$t$ behavior of its survival probability is actually $p(t)\sim\exp[-t^{3/5}]$ in three dimensions.40 It is noteworthy that stealthy trap model systems that prohibit arbitrarily large holes would have finite $T_{1}$ values.

It is noteworthy that the maximum hole size of a solvent is also related to the largest solute particles that it can dissolve. In a solvent with a finite value of $R_{c}$, particles with exclusion radius larger than $R_{c}$ would create intolerably large holes, and would therefore not dissolve. Solute particles smaller than $R_{c}$ would dissolve in a strictly stealthy solvent, but the effective interactions between them deserve future research. Would particles larger than $R_{c}/2$ refuse to touch each other in order to avoid combining the holes they create? Also, if the solute particles are only slightly smaller than $R_{c}$, solvent particles should be concentrated in concentric-shell regions around the solute particles. Could the interference between these concentric shells induce very complicated effective interactions?

Appendix A Expected $R_{c}$ for a finite number of finite-sized Poisson configurations

Although there is no theoretical limit on the hole radii in Poisson configurations (ideal gas), one still expects to find a finite $R_{c}$ if one only studies a finite number of finite sized configurations. If one studies a total of $N_{c}$ configurations of $N$ particles, one expects to see roughly $N_{c}N$ uncorrelated holes. Of these $N_{c}N$ holes, one expects to find the largest hole once. Therefore

[TABLE]

This equation predicts the largest hole size, $R_{c}$, as a function of $\rho$, $N_{c}$, and $N$. To find $R_{c}$ presented in Fig. 2, notice that for stealthy systems of a given $\chi$ and $K$, $\rho$ is given in Eq. (6). Substituting Eq. (6) into Eq. (21) yields

[TABLE]

Here we use $K=1$, $N_{c}=20000$, and $N=500$ to be consistent with stealthy results.