Origin of Symmetry Breaking in the Grasshopper Model

The planar grasshopper problem, originally introduced in (Goulko&Kent 2017 Proc. R. Soc. A 473, 20170494), is a striking example of a model with long-range isotropic interactions whose ground states break rotational symmetry. In this work we analyze and explain the nature of this symmetry breaking with emphasis on the importance of dimensionality. Interestingly, rotational symmetry is recovered in three dimensions for small jumps, which correspond to the non-isotropic cogwheel regime of the two-dimensional problem. We discuss simplified models that reproduce the symmetry properties of the original system in N dimensions. For the full grasshopper model in two dimensions we obtain quantitative predictions for optimal perturbations of the disk. Our analytical results are confirmed by numerical simulations.


I. INTRODUCTION
Ref. [1] introduced the following problem in geometric combinatorics: A grasshopper lands at a random point on a planar lawn of area one.It then jumps once, a fixed distance d, in a random direction.What shape should the lawn be to maximise the chance that the grasshopper remains on the lawn after jumping?
The original motivation [2] to study the grasshopper problem was to formulate and analyze new Bell inequalities for the setup where two parties carry out spin measurements about randomly chosen axes separated by a fixed angle and obtain the spin correlations.For this purpose Ref. [2] formulated the problem on the surface of a sphere, where the grasshopper's jump corresponds to the separation of the measurement axes on the Bloch sphere, and presented some partial results.More detailed analytical results for the grasshopper problem on the surface of the sphere were presented in Ref. [3].Ref. [3] also analysed the problem on the one dimensional circle, which characterises a physically significant subset of qubit measurements often considered in quantum cryptography and quantum information theory.
Here we discuss the planar grasshopper problem and its analogues in higher dimensions.One surprising discovery from Ref. [1] is that, while the jump setup is rotationally symmetric, the resulting optimal lawn shapes are not rotationally symmetric in two dimensions.This is true for arbitrary small jumps d > 0. The specific shapes of the optimal lawns were computed numerically for a range of d.For small jumps, d ≲ R 0,2 , where R 0,2 = 1/ √ π is the radius of the two-dimensional ball with unit volume, i.e. a unit area disk, the optimal lawns shapes resemble cogwheels.The distance between the cogs is approximately the jump length d.The periodic structure of the lawn boundary implies that the rotational symmetry of the jump setup reduces to a dihedral symmetry in two dimensions.While numerical results for the optimal lawn shapes and an analytical proof that the disk shaped lawn is never optimal have already been obtained, so far a fundamental understanding of this symmetry breaking was lacking.The goal of this paper is to systematically study the nature of the emergent symmetry breaking and its implications on the optimal grasshopper lawn shapes.As we will see, the symmetries of the solutions depend strongly on the dimensionality of the problem.

II. PROBLEM AND METHODS
For a formal statement of the grasshopper problem in N dimensions, we denote the lawn domain as Ω ⊂ R N and define the lawn function µ(r) : R N → {0, 1} as the indicator function of Ω, i.e., such that µ(r) = 1 if and only if (1) Note that in the most general version of the problem, lawns may have continuous density in the range (0, 1].The known numerical and theoretical results are consistent with there always being an optimal lawn with density 1.The version of the problem restricting to density 1 lawns is also natural, and we consider it here.The grasshopper probability functional, p µ (d), is then given by where S(N, d) is the surface area of the N -dimensional sphere with radius d, i.e., the area of the region accessible to one jump.Specifically, in two dimensions, the surface area (circumference) of the two-dimensional sphere (disk) equals S(2, d) = 2πd, and in three dimensions, the corresponding area equals S(3, d) = 4πd 2 .The grasshopper probability (2) can be evaluated exactly when the lawn shape is an N -ball with radius R 0,N , by performing the integral (2) explicitly in generalized spherical coordinates.For the special case of a 2-ball, i.e. disk, the probability equals The success probability approaches one as d → 0 and decreases monotonically with increasing d.
In three dimensions, the corresponding probability for jumps smaller than the 3-ball radius, d ≤ R 0,3 = (3/(4π))1/3 , can be computed as Working with the double area integral of a δ-distribution is difficult in more complicated geometries.To overcome this, we exploit the insight that the general probability functional from Eq. ( 2) can be recast as a boundary rather than a volume integral, following a procedure previously used in the study of magnetic fluid monolayers [4][5][6].Specifically, where ∂Ω denotes the boundary of the lawn Ω.Here we abbreviate the interaction as and define Φ to be a radially symmetric solution to ∇ 2 Φ(r) = f (r) (in spherical coordinates).The vectors n are the normal unit vectors orthogonal to the lawn boundary and dS denotes the surface area element corresponding to the volume element dV .To prove Eq. ( 5) we use the divergence theorem (Gauss's theorem). 1We also use the fact that In the last step we used the gradient version of the divergence theorem.Note that dS = n dS, matching (5).
We can work out explicitly the form of the function Φ in the surface integrals.The Laplacian equation for radially symmetric functions in N dimensions is This can be solved by explicit integration.For 2D, the result is while for all N > 2 the result is where Θ is the Heaviside function.Φ is not unique, since one can add to it any function with zero Laplacian and still satisfy its definition.For N > 2 we have chosen a convention where Φ(r) approaches zero at large r, which will be useful in the next section when we want to neglect the contribution of distant boundaries.

III. A MINIMAL MODEL: THE INFINITE HALF-SPACE
The goal of this section is to analyze the stability of the isotropic lawn in the limit of small d in a way that generalizes to arbitrary dimension N .We will see that in two dimensions there is the potential for instability to small perturbations, as found for the disk in [1], while three and higher dimensions such instabilities do not exist.
The starting point is the isotropic solution.Intuitively, one could think that in the limit d ≪ 1 this can be approximated as having zero curvature, i.e. as a half-plane (2d) or half-volume (3d).One needs to be careful in basing arguments on this intuition, since the problem is defined globally rather than locally, and since optimal lawns for small d need not necessarily have small curvature.Also, since the half-spaces have infinite volume, the lawn volume is not normalizable, and (2) is not well-defined.Nonetheless, as we will explain, a perturbative analysis of the half-space problem in N dimensions gives very useful insights into the standard problem framed above.
One could define a version of the problem on the half-space by assuming that the perturbed lawn boundary is a periodic function of a lattice, then normalizing by restricting to one lattice cell and to jumps from initial points within some given distance D from the unperturbed boundary, then taking the limit of large D.Here we simplify the discussion by considering the contribution of small plane wave perturbations, which are periodic in 1 direction rather than N − 1 independent directions (i.e. are the general periodic perturbations only for N = 2).We calculate their differential contribution to an expression of the form (2) for the half-space, which is a meaningful expression assuming a normalization as described.
Let x ∈ R N −1 parameterize the boundary of the lawn and r ∈ R N .We consider a boundary defined by r(x) = (x, ϵh(x)), where the height function h gives a smooth perturbation from the isotropic solution and its magnitude ϵ is such that ϵ ≪ d.The unnormalized normal vector to this surface is n(x) = ∇ r (r N − ϵh(x)) = (−ϵ∇ x h(x), 1), where r N is the N th component of r.Then, the unnormalized probability of remaining on the lawn is given by the double surface integral The stability of the flat plane to small perturbations is governed by the coefficient for p at second order in ϵ.Expanding p − p 0 in ϵ, we find To examine the effect of plane-wave perturbations, we take h(x) = cos(kx 1 ).We expect the effect of such a perturbation to be proportional to the volume of the surface, so we define the stability coefficient δp k (d) to be where L is the linear dimension of the problem, which we will eventually take to infinity.Details of the calculation of this coefficient can be found in Appendix A. The result is The stability of flat half-space in N dimensions to plane-wave perturbations of wavenumber k.All values are negative except at k = 0 (which corresponds to translation of the interface and is not a probability-conserving perturbation) and for kd that are multiples of 2π in the case N = 2.These values are zero, meaning that there is a marginal stability to perturbations with wavenumber commensurate with d in two dimensions but not an instability.For larger dimension N of space, the result becomes increasingly insensitive to the commensurability of k and d.
where J m is the mth Bessel function.For N = 2 and N = 3, we have δp respectively.Some of the resulting functions can be seen plotted in Fig. 1.First, notice that for N ≥ 3 and k ̸ = 0, the result is always negative.This means that for three and higher dimensions the lawn is stable to plane wave instabilities.For N = 2, there is a point of marginal stability for each integer multiple of kd/2π, which corresponds to wavenumbers that are commensurate with the jump size d.As N is increased, the dependence on commensurability with d becomes steadily weaker.There is a simple intuition for this.The average absolute value of the overlap of two random unit vectors in N dimensions is 1 √ N .Hence the expected projection of the jump onto the wavevector of a plane wave instability shrinks like dN −1/2 .The diminishing probability of alignment between the jump and the wavevector pushes any oscillations to higher wavevectors and the variation of this alignment causes the oscillations to decay.Indeed, the curves for various N collapse like N −1/2 d δp k (N 1/2 x/k) with x = N −1/2 kd, approaching a limit curve of the form lim with oscillations totally suppressed.Before delving into the instabilities inherent to the two-dimensional disk, it is worth commenting on the potential role of higher-order perturbations, i.e., ϵ 4 and greater, in determining the fate of marginal stability for N = 2.One might argue that a complete picture would require incorporating these higher-order terms.However, our study focuses on a finite jump size d, and as we will show, this alone is sufficient to drive a compact system from a marginally stable state to true instability.Thus, while high-order perturbations may offer nuanced insights into stability in the half-space version of the problem, they are not relevant to our present investigation of the stability of compact regions in the plane, which centers on the transition driven by finite-size effects.

IV. INSTABILITIES OF THE TWO-DIMENSIONAL DISK
Even with the simplification of moving from the area integrals to line integrals, there are few two dimensional geometries where the problem is analytically tractable.Besides the half-space analyzed above, the disk can also be treated explicitly, as well as small perturbations around the disk.Here, we will compute the second variation of the grasshopper probability for the disk with respect to small perturbations of various discrete symmetry breaking orders, explicitly confirming the findings of [1].
The information provided by this calculation will not tell when a cog of a given number of teeth will be the most stable configuration, since the cog is a large perturbation to the circle.Instead, this variation of the probability describes a kind of spinodal: when does an infinitesimal deformation of the disk towards a given cog raise the probability instead of lowering it?These spinodal points therefore provide some information about when symmetry breaking of a given kind becomes advantageous to isotropy.
Consider a two-dimensional lawn whose boundary is given parametrically in standard polar coordinates by The lawn function, also in polar coordinates, is Here, the radius R ϵ is picked to preserve the area of the lawn.We are interested in small perturbations, and in particular in the stability of the probability of the disk.Specifically, we will evaluate When the coefficient δp n (d) is negative, the perturbation lowers the probability and the circular lawn is stable to sinusoidal perturbations of n-fold symmetry.When it is positive, the probability is increased by the perturbation and the circular lawn is unstable.
Because we only evaluated the quadratic coefficient of the probability in ϵ, we only need the area R ϵ to quadratic order.Since normalization requires that we must set where R 0,2 = 1/ √ π.To treat the line integral, we change the variable of integration from the line element dS to the angular one dθ.This produces a correction to the measures given by the norm of the tangent to the curve, or where the prime denotes the derivative of a function with respect to its parenthetical argument.Noticing that n1 • n2 = t1 • t2 (the scalar product of unit normals is the same as the scalar product of unit tangents for a 1D curve), and that t = r ′ (θ)/∥r ′ (θ)∥, we find where the scalar product is now between the unnormalized tangent vectors.Finally, we change integration variables to θ = θ 1 , ϕ = θ 2 − θ 1 , and have a convenient expression because the integrand in the limit ϵ = 0 is everywhere symmetric in ϕ → −ϕ.
The detailed expansion of (24) to second order in ϵ and the evaluation of the resulting integrals can be found in Appendix C. The resulting stability coefficient is  14).As n is increased, the disk result asymptotically approaches that of the semi-infinite plane.Since the finite-n curves tend to the zeroes of the limit curve from above, the marginally stable points of the half-plane are destabilized at finite n.
where ϕ 0 is given by . This corresponds to the angle at which ∥r n,0 (θ) − r n,0 (θ + ϕ 0 )∥ = d, i.e., the difference in angle along the circle at which that distance between two points is equal to the jump size.This function is plotted for several n in Fig. 2. There, one can see that the local maxima of each curve are positive, indicating instabilities in their vicinity.
The limit of (25) as n → ∞ asymptotically approaches the result (15) found in the previous section for the infinite half-space, as expected.The convergence can be seen in Fig. 3.The mechanism by which the disk is unstable at small d for a large number of oscillations n is revealed: the approach of the finite-n curves to the zeroes of the limit curve is from above.Therefore, the points at which the wavevector k is commensurate with the jump size d, which are marginally stable perturbations of the flat interface, become unstable in this case.The curves of other finite geometries might approach the limit from the other direction, in which case they would be stable against perturbations of this type.
The first zero of δp n (d) for fixed n is significant because it indicates the smallest jump size d for which the disk The location of the first peak is extremely similar to (3.4) of [1] but not identical.If the cutoff in n is increased, more bands appear at higher n.Black markers denote the corresponding numbers of cogs for the optimal solutions found through numerical simulation in Ref. [1] (the same data are shown in the left panel of Fig. 5 in Ref. [1]).The numerical results from [1] are in very close agreement with the current prediction.
is unstable to infinitesimal perturbations of n-fold symmetry of the form (17).The first zero of δp n (d) for given n corresponds to which corresponds to zeros at whose large-n behavior is d 0 ≃ 2 √ πn −1 , in agreement with the prediction from [1].
The stability is examined by seeing where the function first becomes positive.For instance, δp 2 (d) is negative for small d, and crosses zero at d = 3/π, where the circular lawn becomes unstable to elongation.This expression, along with all even-order perturbations, diverges at d = 2R 0,2 , the diameter of the disk, while odd-order perturbations approach zero.Intuitively, the even-order perturbations allow longer jumps across the disk's diameter, while the odd ones do not.
In addition to the first moment of instability, we can analyze which modes are maximally unstable for a given d, i.e., which modes maximize δp n (d).For almost all d, any purported maximum at some value n will be improved by another much higher value.However, we expect that these high frequency modes will be severely penalized by higher-order terms in the expansion of p n,ϵ (d), and so to find better candidates for the relevant symmetry-breaking modes it is useful to truncate the range of n considered.In Fig. 4, we show the values of the two most unstable modes n at a variety of fixed d for n truncated at 18.As evidenced in that figure, the first or second most unstable mode corresponds fairly closely with the actual ground-state symmetry found numerically in [1].Both roughly follow a trend that is described by the location of the first peak in δp n (d).
The stability coefficient (25) is somewhat predictive of what shapes are favored at nonzero ϵ.Fig. 5 shows the numerical probabilities for lawns with some fixed number of cogs for two different values of d.The peaks of these curves correspond well with those predicted by the peaks in the stability coefficient.Moreover, when the probability is rescaled by the amplitude of the cogs, the entire dependence of its variation with n is roughly predicted by our formula.Of course, they do not correspond exactly, since with nonzero ϵ there are higher-order corrections that must enter.In addition, we cannot predict the optimal amplitude ϵ without going to higher order.25), which corresponds roughly with the finite-ϵ numerical data.

V. NUMERICAL RESULTS IN 3D
In previous work [1] numerical results were presented for the grasshopper problem on the 2d plane.In this section we present a similar discussion for three dimensions.The numerical setup in 2d involves defining a discrete version of the grasshopper problem by dividing the plane into grid cells and assigning a spin variable s i to the center of each cell i, where s i = 0 and s i = 1 correspond to the cell being unoccupied or occupied by the lawn, respectively.This setup easily generalizes to higher dimensions.The discrete version of the functional p µ (d) from Eq. ( 2) becomes where M denotes the number of spins and h the edge length of a unit cell.The normalization is M h N = 1.The discretization involves a smoothed approximation δ(r) → δ h (r) = ϕ(r/h)/h to the (one-dimensional) δ-function, which appears in (2).As in [1], we choose ϕ(r/h) = 1 + cos( πr 2h ) /4 if ∥r∥ ≤ 2h and 0 otherwise, following [7].In the continuum limit, P {s} (d) → p µ (d).To resolve the jump distance d we require h ≪ d.We look for spin configurations that maximize P {s} (d), which represent solutions to the discrete grasshopper problem considered.This problem is equivalent to finding the ground state of a spin system with Hamiltonian H = −P {s} (d).This spin system is a conserved spin Ising model with fixed-range interactions, where the interaction range d is large.The optimal configurations are found using simulated annealing or parallel tempering algorithms [1].
To quantify discretization effects in three dimensions, in Fig. 6 we compare the probability for the discretized solid 3-ball to the exact continuous model probability ( 4).An analogous analysis for two dimensions was presented in [1].We can see that as in [1] the discrete model is an excellent approximation of the continuous grasshopper problem and the discretization errors are of the order of 0.3% or less at the highest resolution considered.
In contrast to 2d, in 3d the isotropic solid ball configuration was found to be optimal for jump lengths d ≤ R 0,3 , in agreement with our analytical results presented in the preceding sections.For jumps that slightly exceed R 0,3 the configuration remains isotropic by developing a spherical hole in the center of the ball.We refer to these configurations as shells, specifically 3-shells for N = 3.It is easy to see that the ball cannot be the optimal solution for d > R 0,3 , since in this case a grasshopper starting out near the center of the ball would land outside of the lawn with probability one.Thus it seems intuitive that an optimal lawn may be obtained by removing the central part with radius of the order of R i,3 = d − R 0,3 from the 3-ball and redistributing it elsewhere.This intuition is confirmed by numerics.
In Fig. 7 we analyze discretization effects for the 3-shell configuration.As for the solid 3-ball, the errors due to discretization are very small.Figure 8 shows cross-sections of numerically found optimal grasshopper lawn shapes for two representative jump lengths below and above R 0,3 .
The precise optimal inner radius of the 3-shell configurations can be computed directly, by evaluating the grasshopper probability integral (2) for the radially symmetric 3-shell of a given inner radius in spherical coordinates, and then finding the maximum of the resulting expression as function of inner radius for a given jump.The 3-shell probability at fixed jump length is a simple function of the inner radius with exactly one maximum whose location can be found using standard maximization methods.Fig. 9 shows the analytical results together with the numerically found values for different jumps.Numerical simulations show that the isotropic 3-shell configuration ceases to be optimal for jumps exceeding a critical value of approximately 1.4R 0,3 , when the shell starts exhibiting periodic perturbations: 8-fold perturbations (resembling a rounded cube) were observed for 1.4R 0,3 ≲ d ≲ 1.58R 0,3 and 6-fold perturbations (resembling a bulbous octahedron) were observed for 1.58R 0,3 ≲ d ≲ 1.68R 0,3 .This breaking of full rotational symmetry into a discrete subgroup parallels the symmetry breaking to dihedral symmetry that occurs in the two dimensional model.However, in three dimensions symmetry is only broken for sufficiently large jumps.Both 3d shapes in this category that were observed numerically appear to have full octahedral symmetry (i.e.appear to be invariant under the 48 element group O h ) and appear to correspond to spherical harmonic modulations of the 3-shell.It is conceivable that additional shapes with higher symmetry orders emerge near the transition.These could potentially be observed Right panel: the corresponding optimal grasshopper probabilities for the optimal inner radius (red solid line for analytical value and black dots for the numerical value) and for d − R0,3 (blue dashed lines).The isotropic 3-shell ceases to be optimal for jumps exceeding a critical value of approximately 1.4R0,3.
numerically at higher resolution.
A quantitative procedure to asses whether a given numerical configuration is rotationally symmetric is depicted in Fig. 10.We histogram the radial coordinates of all boundary points of a configuration, and compare these histograms with the known values of the outer and inner radii of the corresponding optimal isotropic configurations.As can be seen in the left panel of Fig. 10, in the 3-shell regime the boundary histograms indeed overlap with the expected values, and their width is of the order of the lattice cell width h.In the symmetry-broken regimes on the other hand, the distribution widths substantially exceed h.
At even larger jumps the numerically found optimal configurations become disconnected.For 1.7R 0,3 ≲ d ≲ 2.45R 0,3 the optimal shapes consist of a ring with two caps, and for larger jumps of two nested crescents.We show representative shapes obtained numerically and the corresponding phase diagram with optimal probabilities in Fig. 11.For better visualization from different perspectives, animations of the 3d shapes are included in the Supplemental Material.

VI. SUMMARY AND DISCUSSION
The setup of the grasshopper problem is rotationally symmetric: the starting point and the direction of the jump are both chosen uniformly.For small jumps, in order for the grasshopper to be able to land outside the lawn, the starting point needs to be close enough to the boundary.Thus one might intuitively assume that the optimal grasshopper lawn shape for small jumps would minimize the boundary.The shape that minimizes the boundary at fixed volume is the N -ball, which in two dimensions corresponds to the disk.However, numerical simulations in 2d revealed that optimal lawn shapes were not rotationally symmetric for any jump length, and this result was also proven analytically [1].For small jumps, full rotational symmetry was found to be reduced to dihedral symmetry D n , where the number n of remaining rotational symmetries is determined by the number of edges of a polygon inscribed in the unit area disk, such that the polygon edge length is close to the length of the grasshopper jump d.
In this work we answered the open question regarding the origin of the symmetry breaking in the grasshopper model.In the limit of small jumps, d → 0, a disk-shaped lawn can be approximated as having zero curvature.We showed that in this minimal model of an infinite half-plane, periodic perturbations of the boundary that are commensurate to the jump lengths d are only marginally stable.For the same model in higher dimensions we showed that all perturbations are stable and the result becomes increasingly insensitive to the commensurability between the perturbation and jump.Numerical simulations in three dimensions confirm that optimal solutions are indeed isotropic: for small jumps the lawns are 3-balls, which subsequently develop spherical holes as the jump size increases.Rotational symmetry is only broken for large values of the jump (d ≳ 1.4R 0,3 ).
To show that marginally stable perturbations of the minimal model indeed become unstable in the full model, we analyzed the instabilities of the 2d disk with respect to small perturbations of different discrete symmetry breaking orders n.We confirmed the dihedral symmetry breaking described in Ref. [1] and found additional discrete rotationally symmetric solutions that are advantageous compared to the disk.While this perturbative calculation does not directly describe the optimal grasshopper shapes found numerically, there is remarkable agreement between the perturbative predictions and the numerical results.Moreover, the n → ∞ limit asymptotically approaches the result of the minimal model.
The analytical and numerical techniques we have developed hold promise for broader applications, particularly in elucidating long-range interaction effects in ultracold-atom based quantum emulators and prototype quantum computers.Several optimal grasshopper shapes, including cogwheel patterns, suggest intriguing parallels with dipolar Bose-Einstein condensates, where similar patterns at droplet boundaries have been reported [8,9].Another promising application of our techniques is the domain dynamics of quantum Ising spins in Rydberg atom arrays [10].Therefore, our work serves not only as an illuminating analysis of symmetry breaking in the grasshopper model but also as a stepping stone for the study of complex quantum many-body systems with long-range interactions.In addition, recent advances in the creation and control of Rydberg atom arrays with distance-selective interactions [11,12] present an exciting opportunity to directly build and study two and three dimensional grasshopper-type spin systems in ultracold atoms experiments.
It follows that the stability coefficient resulting from the sum of two plane waves is given in this case by the sum of the coefficients due to the individual pure plane waves.Therefore, the effect of more complicated perturbations than those considered here can be reduced to that of the simple ones, with the caveat that these results only strictly hold for perturbations build from sums of finitely many oscillatory terms.12.The points give the difference in probability between the perturbed disk and the disk as a function of perturbation size ϵ for n = 2 and d = 0.98 > d0.They were computed using numeric integration on the full expression.The solid line gives δp2(d)ϵ 2 .The results agree at small ϵ, as expected.
The first piece is proportional to the integral over a Heaviside theta, and this simply changes the limits of integration.We have Next, we have Finally, we treat the last term.First, a helpful identity for evaluating an integral of the derivative of a δ-function convolved with a function.Note that ∂ ∂x δ(g(x)) = δ ′ (g(x))g ′ (x).Using manipulations and integration by parts, we have FIG.1.The stability of flat half-space in N dimensions to plane-wave perturbations of wavenumber k.All values are negative except at k = 0 (which corresponds to translation of the interface and is not a probability-conserving perturbation) and for kd that are multiples of 2π in the case N = 2.These values are zero, meaning that there is a marginal stability to perturbations with wavenumber commensurate with d in two dimensions but not an instability.For larger dimension N of space, the result becomes increasingly insensitive to the commensurability of k and d.

11 FIG. 2 .FIG. 3 .
FIG.2.The coefficient of stability for the disk to small perturbations of n-fold symmetry.A large dot is drawn at the smallest value of d where each curve first becomes positive.The black bar shows the point in d at which the transition to disconnected shapes occurs as measured in[1].

FIG. 4 .
FIG.4.This plot shows the 1st and 2nd most unstable modes n at a given d for n ≤ 18.The blue and yellow lines show the location of the first and second peak in δpn(d), respectively.The location of the first peak is extremely similar to (3.4) of[1] but not identical.If the cutoff in n is increased, more bands appear at higher n.Black markers denote the corresponding numbers of cogs for the optimal solutions found through numerical simulation in Ref.[1] (the same data are shown in the left panel of Fig.5in Ref.[1]).The numerical results from[1] are in very close agreement with the current prediction.

FIG. 5 .
FIG. 5. Left: Grasshopper probability for lawns with fixed cog number for two values of d, obtained numerically for discrete lawns with 10, 000 cells.The maxima corresponding to the two leading unstable modes are clearly visible (n = 9 and n = 17 for d = 0.4; n = 7 and n = 15 for d = 0.46) and are marked with thin vertical lines.These are indeed (local) maxima of the grasshopper problem.Horizontal dashed lines denote the corresponding disk probabilities given by Eq. (3).Right: The same probability rescaled by the squared amplitude ϵ of the cogs (symbols with error bars).The solid lines show δpn(d) from (25), which corresponds roughly with the finite-ϵ numerical data.

5 FIG. 6 .FIG. 7 .
FIG.6.Study of discretization effects.Left panel: The exact continuous probability functional pµ(d) for the solid 3-ball of unit volume given by Eq. (4) (red solid line) compared with the corresponding discrete P {s} (d) (black dots) as function of the grasshopper jump distance d.For d ≤ R0,3 the 3-ball configuration is the optimal lawn shape.Right panel: Relative deviation of P {s} (d) for the solid 3-ball configuration from the corresponding pµ(d) as function of the lattice spacing h for two representative values of the grasshopper jump: d = 0.2 ≈ 0.3R0,3 (blue dots and line) and d = 0.5 ≈ 0.8R0,3 (red dots and line).For the highest resolutions considered (M ≈ 113, 000) the discretization error is well below 0.3%.Lines are to guide the eye.

FIG. 8 . 3 FIG. 9 .
FIG.8.Left panel: cross section of the optimal configuration for d = 0.64R0,3 found numerically for a system with M = 160, 000 spins.The configuration has the shape of a solid 3-ball.This configuration was found to be optimal for all d ≤ R0,3.Right panel: cross section of the optimal configuration for d = 1.32R0,3 found numerically for a system with M = 160, 000 spins.If the jump length exceeds R0,3 the configurations remain isotropic (for d ≲ 1.4R0,3) but develop a spherical hole in the center; the radius of the hole grows with increasing d.Note that the outer radius of the configuration is slightly larger than R0,3 to ensure that it has unit volume.

3 FIG. 10 .
FIG. 10.Histograms of the radial coordinates of the configuration boundary points for, from left to right, d = 1.39R0,3 (isotropic 3-shell), d = 1.42R0,3 (8-fold perturbation), and d = 1.58R0,3 (6-fold perturbation).Results shown were obtained for systems with M = 160, 000 spins.The vertical red lines denote the theoretical values for the optimal (for the respective value of d) inner and outer radii of the corresponding 3-shell configurations.The shaded regions mark the ±h interval around the optimal radii.
FIG.12.The points give the difference in probability between the perturbed disk and the disk as a function of perturbation size ϵ for n = 2 and d = 0.98 > d0.They were computed using numeric integration on the full expression.The solid line gives δp2(d)ϵ 2 .The results agree at small ϵ, as expected.