Does the Round Sphere Maximize the Free Energy of (2+1)-Dimensional QFTs?

We examine the renormalized free energy of the free Dirac fermion and the free scalar on a (2+1)-dimensional geometry $\mathbb{R} \times \Sigma$, with $\Sigma$ having spherical topology and prescribed area. Using heat kernel methods, we perturbatively compute this energy when $\Sigma$ is a small deformation of the round sphere, finding that at any temperature the round sphere is a local maximum. At low temperature the free energy difference is due to the Casimir effect. We then numerically compute this free energy for a class of large axisymmetric deformations, providing evidence that the round sphere globally maximizes it, and we show that the free energy difference relative to the round sphere is unbounded below as the geometry on $\Sigma$ becomes singular. Both our perturbative and numerical results in fact stem from the stronger finding that the difference between the heat kernels of the round sphere and a deformed sphere always appears to have definite sign. We investigate the relevance of our results to physical systems like monolayer graphene consisting of a membrane supporting relativistic QFT degrees of freedom.


Introduction
The equilibrium configuration of a physical membrane is often determined by a competition between several physical effects. For instance, the round shape of a soap bubble arises from a competition between the bubble's intrinsic surface tension, which energetically prefers to collapse it, and a pressure differential between the air inside and outside of the bubble. Likewise, the bending modulus of simple lipid bilayers tends to flatten them, with thermodynamic effects potentially causing deformations.
Here we are specifically interested in membranes supporting relativistic quantum degrees of freedom living on them. A fiducial example is that of a graphene monolayer, whose energetics in a Born-Oppenheimer-like approximation can be split into a sum of two contributions: one from the atomic background lattice and another from relativistic excitations that propagate on this background. At scales well above the lattice spacing, the lattice can be treated as a continuous membrane, and the free energy will depend on its geometry. The contribution to this free energy from the background, interpreted as a classical contribution F c , is then captured by a Landau free energy constructed from its embedding into an ambient flat space [1][2][3] (see [4] for a review). The effective relativistic excitations are two free massless Dirac fermions (with effective speed of light c eff given by the Fermi velocity ∼ 10 6 m/s), and their free energy -interpreted as a quantum contribution F q that depends on the membrane's instrinsic geometry -is computed via an appropriate path integral. Since the classical contribution to the free energy is relatively well-understood, our goal is to understand the contribution from the Dirac fermions.
In fact, while graphene is our motivating physical system (and one to which we will often turn for physical interpretation), there exist other more exotic examples of relativistic quantum fields supported on membranes -for instance, domain walls in cosmology [5] or braneworld models of our universe [6,7]. Moreover, one could also imagine engineering other graphene-like two-dimensional crystalline materials that exhibit relativistic excitations living on them. Consequently, in this paper our main objective is to study the free energy of more general classes of relativistic QFTs living on (2 + 1)-dimensional geometries, with the massless Dirac fermion corresponding to graphene as a special case.  Table 1. A summary of results for the free energy F q (or just energy E in the fourth row) for various types QFTs at various temperatures. "Holographic CFTs" refer to CFTs dual to smooth geometries obeying Einstein's equations; the scalar refers to a scalar field of any mass and curvature coupling; and the fermion refers to a Dirac fermion of any mass.ḡ is always taken to be the maximally symmetric geometry on the manifold Σ (soḡ is the round sphere metric when Σ = S 2 and the flat metric when Σ = R 2 or T 2 ), and when we writeḡ + h it is understood that the result holds to leading nontrivial order in the perturbative expansion parameter . "Long wavelength" refers to metrics whose curvature is small compared to the temperature and/or mass of the field. All inequalities on ∆F q and ∆E are saturated if and only if g =ḡ.
Interestingly, previous work has found that such relativistic quantum fields tend to energetically prefer deformed geometries. To briefly summarize, let us assume that this field theory lives on R × Σ, where Σ is a two-dimensional spatial manifold with metric g. Per the discussion above, F q will depend on the intrinsic geometry of Σ, i.e. on g. Because F q is a free energy, it is extensive, and therefore in order to sensibly discuss its dependence on the shape of Σ we should imagine keeping the volume of Σ (computed with respect to g) fixed as we vary g 1 . We therefore consider the "background-subtracted" free energy ∆F q ≡ F q [g] − F q [ḡ], withḡ taken to be some fiducial reference metric which endows Σ with the same volume as g. With this understanding, the relevant extant results are summarized in Table 1. The key takeaway is that ∆F q is negative for many different theories, a result most well-established when the geometry on Σ is perturbatively close to the round sphere or the flat plane. Note that ∆F q remains nonzero even at zero temperature, when it can be interpreted as a Casimir energy.
These observations naturally lead to the following question: if the free energy F q governs the equilibrium configuration of a membrane, does the fact that ∆F q is always negative lead to an instability of the round sphere or flat space? If so, will a membrane settle down to some less-symmetric equilibrium configuration, or does this instability ultimately lead to a runaway process (which presumably breaks down once a UV scale is reached)? Answering this question will of course depend on how ∆F q competes with other contributions to the free energy; returning to the case of graphene, in [14] we performed a parametric comparison of the competition between ∆F q and the classical bending free energy ∆F c , finding that the typical curvature scale l crit at which the negative contribution of ∆F q becomes dominant over the (positive) contribution of ∆F c agrees with the "rippling" length scale l rip of graphene measured in experiments [15]. However, the order of magnitude of this scale is only slightly above that of the lattice spacing, where the effective Dirac fermion description breaks down, so the validity of our estimates for ∆F q and ∆F c in this regime is suspect. Moreover, the order-of-magnitude analysis of [14] was made subtle for two reasons. First, since Σ is a plane, its volume is infinite; hence one must be careful in defining precisely what is meant by the condition that the volumes computed from g andḡ match. Second, the leading perturbation to ∆F c is linear in the deformation amplitude , while ∆F q is quadratic in ; hence balancing these two contributions requires a careful accounting of various orders-of-limits.
Our first purpose here is therefore to repeat the perturbative analysis of [14] -that is, the perturbative computation of ∆F q for the nonminimally coupled free scalar and the free Dirac fermion -in the case where Σ is a topological sphere, rather than a plane. This modification alleviates both of the issues just mentioned, since when Σ is a sphere it has finite volume and we also find that both the contributions ∆F c and ∆F q are quadratic in . Again we find that the round sphere locally mazimizes F q .
Our second purpose is then to investigate the questions posed above: namely, is the round sphere a global maximum of ∆F q ? Does ∆F q eventually find some new equilibrium configuration after a sufficiently large deformation to the sphere, or can it decrease indefinitely? To address these questions, we numerically compute ∆F q for large (axisymmetric) deformations of the round sphere, finding that it is always negative, and in fact that it can be made arbitrarily negative as the geometry becomes singular. Our conclusion, therefore, is that the energetics of ∆F q favor geometries that are not smooth. Surprisingly, we also find that the behavior of ∆F q for such large deformations is remarkably similar for the scalar and the fermion when normalized by its perturbative expression.
In fact, our main result is stronger: not only is ∆F q always negative for the deformations we study, but the heat kernel ∆K(t) which computes ∆F q has definite sign for all t. This heat kernel will be introduced in more detail in Section 2, but in short it is related to the eigenvalues λ I of the operator L that defines the equations of motion of the free fields: ∆K L (t) = I e −tλ I − e −tλ I , (1.1) where λ I andλ I are the eigenvalues of L on the deformed sphere and the round sphere, respectively. The fact that ∆K L (t) apparently has fixed sign for all t is therefore a nontrivial statement about the behavior of the eigenvalues of L. The universality of this result leads us to conjecture that ∆K L (t) has fixed sign for any free field theory and area-preserving deformation of the sphere. The order-of-magnitude analysis of the competition between ∆F c and ∆F q is provided immediately below, in Section 1.1, for the sake of illustrating more clearly some of the concepts discussed so far. We then establish our setup, conventions, and formalism in Section 2, focusing specifically on the computation of ∆F q using heat kernels. We then present the perturbative calculation of ∆F q for the free nonminimally coupled scalar and the free Dirac fermion in Section 3, along with some checks showing that our results reproduce the CFT result of [13] and the flat-space result of [14] in appropriate limits. We then present the numerical calculation of ∆F q for large deformations of the sphere in Section 4, focusing for simplicity on axisymmetric perturbations. Finding that F q seems to decrease monotonically as the amplitude of the deformation is increased, in Section 5 we analyze its behavior on extremely deformed geometries, showing that approaching a conical singularity allows ∆F q to become arbitrarily negative. Section 6 concludes with a summary of our main conclusions and unexpected results.

A Perturbative Example: Graphene
For illustrative purposes, let us now study in some more detail the competition between ∆F c and ∆F q for two-dimensional crystalline materials such as monolayer graphene, focusing on the case where Σ is a small deformation of a round sphere (for large deformations, this competition will be discussed in Section 5.4). We note there is considerable technological interest in producing spherical monolayer graphene (see for example [16]). As mentioned above, at scales much larger than the lattice spacing this crystal can be described as a smooth membrane, and the free energy will depend on the geometry of this membrane. For simplicity, we will further assume that this effective description is diffeomorphism-invariant (although of course this is not expected to be the case for a crystalline material like graphene [4]).
For an order-of-magnitude estimate of perturbations to the round sphere, we also assume that the sphere minimizes the classical bending contribution ∆F c to the free energy. Keeping only up to second derivatives, we may then write the Landau free energy as where K is the mean curvature of Σ, κ is a bending rigidity, r 0 is the radius of the sphere that minimizes ∆F c , and no term containing the scalar curvature of g appears because such a term is topological. Working perturbatively around the sphere of radius r 0 , we write the ambient flat space in the usual spherical coordinates and take {θ, φ} as coordinates on Σ and embed Σ as r = r 0 (1 + f (θ, φ)), where is a dimensionless expansion parameter. To linear order in , the induced metric on Σ is in a gauge conformal to the round sphere, where ∇ 2 denotes the Laplacian on the round sphere of unit radius. We then decompose f in spherical harmonics as with the condition that the volume of Σ remain unchanged imposing that f 0,0 = 0. We thus obtain The general contribution of quantum scalar or Dirac fermionic fields to ∆F q is obtained in Section 3 below. To streamline the present analysis, let us take the relativistic quantum fields living on Σ to be a CFT; this is the case for graphene when it is slightly perturbed from a flat plane (additional gauge fields associated to the underlying lattice structure vanish when the the metric is in a conformally flat form) [17][18][19][20]. Here we assume the effective CFT description remains valid even for small perturbations of the 2 Breaking diffeomorphism invariance would allow for more general coefficients in front of the f 2 , f ∇ 2 f , and (∇ 2 f ) 2 terms. round sphere. At zero temperature, the contribution of these degrees of freedom to ∆F q is [13] where c eff is the effective speed of light for these relativistic degrees of freedom and c T is the central charge (defined as the coefficient in the two-point function of the stress tensor); in our conventions, the central charges of a conformally coupled massless scalar field and of a massless Dirac fermion are c T = (3/2)/(4π) 2 and c T = 3/(4π) 2 , respectively [21,22]. Importantly, at large the coefficients in the sum grow like 3 , indicating that this contribution is non-local: that is, unlike ∆F c it does not arise from some local geometric functional. Also note that while technically (1.8) is only valid at zero temperature, the leading corrections to it go like e −l 2 is a thermal length scale and l is the typical length scale of the perturbation f ; hence (1.8) holds for l T 2l. (The corrections to the zero-temperature result will be discussed in Section 4.3.) The combined contribution to the free energy from the classical and quantum contributions therefore goes like 3 is some characteristic length scale and (1.11) We would like to investigate whether this combined expression can ever be negative in its regime of validity. This question can be investigated as follows: first note that at large , A (c) goes like 4 while A (q) only grows like 3 , so the positive classical free energy will always dominate at sufficiently high angular momentum quantum number. We must therefore investigate the behavior of the lowest modes: = 0 does not contribute since f 0,0 = 0, while the contribution of = 1 modes to both ∆F c and ∆F q vanishes due to the fact that such deformations correspond to infinitesimal diffeomorphisms. However, since A (c) vanishes quadratically around = 1 while A (q) only vanishes linearly, it is clear that for sufficiently small − 1 > 0, A (c) − (γ/r 0 )A (q) < 0. Since is an integer, making ∆F negative therefore requires this to be true all the way to = 2, and hence Since our analysis is only valid at scales well above the lattice spacing a, we also require r 0 a, which implies γ a. For the particular case of graphene, typically the bending rigidity is taken as κ ∼ 1 eV, a ∼ 2.5 Å, c eff ∼ 10 6 m/s [23], and c T = 2 × 3/(4π) 2 (the factor of two coming from the two Dirac points in graphene's band structure), from which one finds γ/a ∼ 0.1 (note that the numerical prefactors matter: a purely parametric estimate would give γ/a ∼ c eff /aκ ∼ 10). Hence for graphene it does not seem likely that the quantum effect we have identified can ever compete with the classical bending energy to render the round sphere unstable, even if one were to keep more careful track of the precise form of the Landau free energy. In the absense of fine-tuning, this result could have been expected: with no fine-tuning, the energy scale κ should be set by the lattice spacing and hence κ ∼ c eff /a, from which it would follow that γ/a is order unity 4 . We therefore interpret the condition γ a as the required fine-tuning of the membrane parameters (e.g. κ, c eff ) that makes it possible for ∆F q to dominate over ∆F c . Given the great current interest in monolayer graphene-like materials, conceivably such fine-tuned crystalline membranes could be engineered in a lab.

Setup
We consider thermal states of (2+1)-dimensional (unitary, relativistic) QFTs on the geometry R × Σ, where Σ is a two-dimensional manifold with sphere topology. The Euclidean continuation of this geometry is with the period of Euclidean time τ given by the inverse temperature β = 1/T , and we have made explicit the fact that the spatial metric g ij (x) on Σ is independent of τ . The free energy is thus a functional of g ij and of β; to simplify notation, we will denote this free energy simply as F [β, g] (i.e. without the subscript q as was used above).

Free Energy
The desired free energy F is determined by the Euclidean partition function Z, which will depend on both β and the spatial geometry g ij : where S E is the Euclidean action and Φ schematically stands for the QFT fields in the system. Of course, as written Z (and thus F ) is UV-divergent, so we must regulate it.
Since we are only considering relativistic QFTs, any UV regulator (like, say, a lattice) cannot break diffeomorphism invariance in the IR, and hence for simplicity we may use a covariant UV regulator to ultimately compute UV-finite quantities. To that end, note that for a UV cutoff Λ, the most general covariant counterterms that can be added to the Euclidean action are where µ schematically stands for any parameter in the QFT with dimensions of energy, if one exists (for instance, a mass), R is the Ricci scalar of g, and the theory-dependent coefficients c i are dimensionless and independent of Λ and of the geometry. Hence the most general divergence structure of the free energy takes the form where χ Σ is the Euler characteristic of Σ and F fin [β, g] is finite as Λ → ∞. Note in particular that the divergence structure depends on g only through the volume Vol[g] of Σ; in the context of two-dimensional crystalline lattices discussed in Section 1, one can think of these terms as contributing to some (UV cutoff-dependent) tension in the classical membrane action. In other words, we may interpret the volume preservation condition as merely a convenient way of grouping the leading-order divergences in (2.4) with the couplings in the classical membrane action. Physical information about the free energy is contained in the finite part F fin , but this object is not uniquely defined by the expansion (2.4) (since a general change in the UV cutoff can induce a change in F fin ). However, the differenced free energy ∆F ≡ F [β, g] − F [β,ḡ] discussed above (in whichḡ is a reference metric such that Vol[ḡ] = Vol[g]) is scheme-independent. As shown in [14], this differenced free energy can be defined via where ∆S E is the difference of the Euclidean actions constructed from g andḡ and the expectation value on the right-hand side is taken in the thermal vacuum state (of inverse temperature β) associated to the geometryḡ.

Heat Kernels
Let us now restrict to the case where the QFT fields Φ are free; in such a case, the Euclidean action is quadratic, and the partition function reduces to a functional determinant. The free energy is then conveniently evaluated via heat kernel methods, which we now review. The massive free scalar fields and Dirac fermions on which we focus have actions where ξ is the curvature coupling of the scalar, M is a mass, and the spinor conventions are as in [14]. Performing the path integral on the geometry (2.1), one obtains [14] where σ = −1/2 (+1) for the scalar (fermion) and L is a differential operator on Σ.
For the non-minimally coupled scalar we have simply L = −∇ 2 + ξR, which acts on functions with spin weight zero. The case of the Dirac fermion is slightly more complicated and we will give the full expression for L in (3.22) below, but the key idea is that the square of the Dirac operator on the ultrastatic geometry (2.1) is diagonal in the spinor indices and L is one of these two diagonal components, which acts on functions of spin weight 1/2. We now define the heat kernel K L (t) ≡ Tr(e −tL ), in terms of which the free energy is This form of the free energy makes manifest its UV divergence structure, as UV divergences are associated with small t in the above integral. More explicitly, by the heat kernel expansion [24] the small-t behavior of K L (t) goes like where the coefficients b 2n can be expressed as integrals of local geometric invariants on (Σ, g ij ). UV divergences are controlled by the leading and subleading coefficients b 0 and b 2 , which depend only on the volume and topology of (Σ, g ij ) (though they are otherwise theory-dependent): Thus the differenced free energy can be obtained directly from the difference ∆K L (t) between heat kernels corresponding to the spatial geometries (Σ, g ij ) and (Σ,ḡ ij ): It is clear from (2.9) and (2.10) that as long as (Σ, g ij ) and (Σ,ḡ ij ) have the same volume and topology, ∆K L is O(t 1/2 ) at small t, and thus that ∆F is UV-finite, as expected from the arguments above. Now we may specify to our case of interest: using the decomposition (2.7) for L, is the difference of heat kernels of the operators L on the two-dimensional geometries (Σ, g ij ) and (Σ,ḡ ij ) and we have defined which arises from a sum over Matsubara frequencies on the thermal circle. Moreover, we will be concerned with the case where Σ is a (topological) sphere, in which case it is natural to take the reference metricḡ ij to be that of a round sphere. Finally, we note that the heat kernel expansion (2.9) for ∆K L (t) takes the form where the ∆b 2n are the differences of the heat kernel coefficients between the geometries (Σ, g ij ) and (Σ,ḡ ij ).

Perturbative Results
The expression (2.12) for the differenced free energy in terms of the heat kernel of L is convenient because it simply requires computing the variation in the spectrum of L as the spatial geometry g is varied: where I indexes the eigenvalues of L and L. An explicit computation of this perturbed heat kernel was performed for deformations of flat space in [14], with the key result that for both the fermion and the scalar, to leading nontrivial order σ∆K L (t) is negative for all t (and hence ∆F is negative for all perturbations). In order to compare to our later results, we now repeat this calculation on the perturbed round sphere (3.2). We remind the reader that the reasons for working on the sphere are twofold: first, since the sphere is compact we don't have to deal with IR divergences; second, we will find that for small perturbations of the round sphere, the free energy of quantum fields is of the same order as the contribution from the classical membrane free energy, and hence the two can consistently be compared.
In this Section, we will take the metric on Σ to be conformal to the round sphere, as in (1.4): where f is some scalar field on the sphere and where we are using units in which r 0 = 1. We ; the reference metric corresponds to taking = 0. The volume preservation condition thus requires that 5 Similarly, we write the resulting expansion of L and of its eigenvalues λ I and eigenvectors h I as where the eplicit expressions for L (1) and L (2) in terms of f (1) and f (2) are provided in Appendix A. Hence from (3.1), the perturbed heat kernel is Homogeneity of the round sphere implies that the leading variation of ∆K L is quadratic, so ∆K (1) = 0 which we indeed find shortly [11]. Now, defining the matrix elements a standard consistency condition in degenerate perturbation theory requires that L (1) IJ be diagonal on any degenerate subspaces of L (that is, we must have L (1) IJ = 0 for any I, J with I = J butλ I =λ J ) 6 . Then standard perturbation theory yields the perturbations of the eigenvalues: It is important to note that while consistency of the perturbation theory requires an appropriate choice of the unperturbed eigenfunctionsh I , in fact the final expression for the heat kernel is insensitive to this choice. To see this, let us write the index I as the pair ( , m), with labeling each degenerate subspace of degeneracy d and m indexing its elements 7 . Then we may relate the eigenfunctionsh ,m to any other basish ,m by a unitary transformation on each degenerate subspace: where c m m are the components of a unitary matrix chosen to ensure that L ,m, ,m = 0 for m = m . We then have where bold characters denote matrices on the degenerate subspaces, so that e.g. c is the d × d -dimensional matrix with elements c m m , L , is the d × d -dimensional matrix with elements L (1) ,m, ,m , etc. Hence with the final expression following from the basis-independence of the trace. Likewise, we have , is required to be diagonal, the first sum in the square brackets can be written simply as Tr((L (1) , ) 2 ). Then again using (3.9) and cyclicity of the trace, we find that All dependence on c has vanished due to the traces, and hence for the purposes of computing the heat kernel we may compute the matrix elements L (n) ,m, ,m in any desired basish ,m .

Scalar
For the scalar, the operator L for general f is with ∇ a the covariant derivative on the round sphere (f = 0). The unperturbed operator L is −∇ 2 + 2ξ and has eigenvaluesλ = ( + 1) + 2ξ, with ∈ {0, 1, 2, . . .} a non-negative integer. For the computation of the matrix elements L (n) ,m, ,m , we may take the eigenfunctionsh ,m to just be the usual spherical harmonics Y ,m . Since the calculation is rather cumbersome and unilluminating, we relegate it to Appendix A; in short, expanding f (1) in spherical harmonics as for the non-minimally coupled scalar one ultimately obtains ∆K (1) = 0 and with the general expressions for α , and β , given in (A.15b) and (A.18) in the Appendix. For the special case of odd , the expressions simplify substantially to 8 where (x) n ≡ Γ(x + n)/Γ(x) are Pochhammer symbols. We will comment further on this expression in Section 3.5 below.

Dirac Fermion
For the benefit of the reader, let us briefly summarize how to obtain the operator L for the fermion; more details can be found in [14]. We first evaluate Levi-Civita connection on the full (three-dimensional) Euclidean geometry, ω aµν the spin connection, and S µν the generators of the Lorentz group. Evaluating this object in the ultrastatic geometry (2.1), one finds that it is diagonal in its spinor indices: where L is the operator introduced in (2.7) and P L,R are projectors onto left-and righthelicity Weyl spinors on the two-dimensional geometry Σ; is it this decomposition that allows us to compute the fermion partition function from just the spectrum of the (non-spinorial) operator L. The explicit form of the operator L defining L can be given most easily by working in conformally flat coordinates on Σ, (3.20) The expression adapted to the spherical coordinates of (3.2) can be obtained easily by transforming from the conformally flat coordinates {x 1 , x 2 } to the spherical coordinates {θ, φ} via sin θ = sech x 1 , φ = x 2 ; then sincef = f + ln sin θ, in terms of the conformal factor f one ultimately obtains 9 where as before ∇ a is the covariant derivative on the round sphere. L acts on functions with spin weight 1/2, and hence the unperturbed eigenfunctionsh ,m can be taken to be the spin-weighted spherical harmonics 1/2 Y ,m of spin weight 1/2, where ∈ {1/2, 3/2, 5/2, . . .} is a positive half odd integer and as usual m ∈ {− , − + 1, . . . , }.
The corresponding unperturbed eigenvalues areλ = ( + 1/2) 2 . Again we relegate the details of the computation of the heat kernel to Appendix A; ultimately we obtain ∆K (1) = 0 and with the expressions for α , and β , given in (A.29). As for the scalar, taking to be odd substantially simplifies them:

Check: Conformal Field Theories
As a simple check of our results, let us compare to the results of [13], which computed the zero-temperature perturbative energy difference ∆E (2) for any unitary conformal 9 A more covariant expression can be given by introducing the spin weight raising and lowering operators ð,ð, in terms of which more details are presented in Appendix A.
field theory. There it was found that in any CFT, this leading-order energy difference is with c T the central charge defined as the coefficient in the two-point function of the stress tensor. We now show that our expressions (3.15) and (3.23) reproduce (3.25) with the correct central charges when the fields are conformal; we note that in this case,λ = ( + 1/2) 2 for both the scalar and the fermion (though the allowed values of of course still differ). To do so, first note that the free energy difference is given by inserting (3.15) and (3.23) into (2.12); for simplicity we will restrict to perturbations f ,m with odd , so that we may use the more compact expressions (3.16) and (3.24). In the zerotemperature limit, the integral over t can be performed explicitly by noting that Poisson resummation gives (for both the scalar and the fermion) for any t > 0; hence the zero-temperature perturbative energy difference for odd is where we used the fact that the sum over is finite to integrate term-by-term (it is understood that the sum over runs over integers or half-integers depending on whether we are considering the scalar or the fermion, with α , the corresponding expression; we also remind the reader that σ = −1/2 for the scalar and σ = 1 for the fermion). We therefore have where in the expression for A (ferm) we shifted the index of summation by 1/2. While we are unable to analytically show that these expressions reproduce the form (3.25) predicted by CFT perturbation theory, by computing these sums exactly we find that they do, with the correct central charges c T = (3/2)/(4π) 2 and c T = 3/(4π) 2 (we have checked up to = 1001). Interestingly, if is even then evaluating A by integrating term-by-term produces a divergent sum, presumably due to the fact that the (now infinite) sum over in a (t) doesn't commute with the integration over t. Nevertheless, the behavior of a (t) for even makes clear that the integral is indeed finite when performed after the summation, and we have confirmed numerically that it reproduces (3.25) for a range of even .

Check: The Flat Space Limit
As a final check of our results, let us consider the limit in which the radius of the sphere is taken to be very large, and only modes with large , m are excited. In this limit, we expect the theory to be insensitive to the curvature of the sphere, and thus the heat kernel should reproduce its flat space behavior. This behavior was computed for both the scalar field and the Dirac fermion in [14], in which it was found that when the perturbed metric is in the conformally flat form ds 2 = e 2f δ ij dx i dx j , the perturbation to the heat kernel is where k is a wave vector defined by the Fourier decomposition of f (1) as

31)
k = | k| is its magnitude, and the functions I(ζ) are given for the scalar and fermion as We now introduce an appropriate flat-space scaling limit in which our expressions for ∆K (2) reproduce (3.30). To do so, let us explicitly reintroduce the radius r 0 of the sphere, so that the deformed sphere metric (3.2) becomes ds 2 = r 2 0 e 2f dθ 2 + sin 2 θ dφ 2 . (3.33) The scaling limit is defined by "zooming in" on a point on the equator of the sphere by introducing new coordinates x = r 0 (θ−π/2), y = r 0 φ and then taking the limit r 0 → ∞ with x, y held fixed. The resulting metric is in the desired conformally flat form, with x and y having infinite range. Restoring r 0 to the expressions (3.15) and (3.23), we obtain with α , and β , unchanged. Now let us again focus on the case where f (1) only contains modes with odd , so β , vanishes and the sum over runs to /2. In order to consider modes with large , we define k = /r 0 and k = /r 0 and keep k and k fixed as we take r 0 → ∞. As we show in Appendix B, in this limit we find that with I(ζ) precisely the functions given in (3.32); assuming a r 0 k (t) is continuous in k in this scaling limit, we may now remove the restriction to modes with odd r 0 k. We also find that as long as , the upper (lower) signs correspond to even (odd) (k + k y )r 0 , and we are neglecting an overall phase that will cancel out. Inserting these expressions into (3.35) and decomposing the sum over = kr 0 into sums over even and odd kr 0 , we finally obtain precisely the flat-space expression (3.30) given in [14]: It is perhaps worth emphasizing that computing the perturbation to the free energy of a perturbation of flat space is rather subtle due to the requirement that the perturbed and unperturbed geometries have the same volume: since the volume of flat space is infinite, an IR divergence is introduced, and the volume preservation condition is interpreted as controlling this IR divergence to yield a finite differenced free energy. In [14], this problem was addressed by computing the heat kernel on a torus and then taking the limit in which the cycles of the torus go to infinity; this is analogous to the procedure performed here, where we computed the heat kernel on the sphere and then took a flat space scaling limit. In these regularization schemes, the "extra" bits of the torus or the sphere that get sent to infinity in the flat space limit essentially deform in such a way as to ensure that the leading-order UV divergences in (2.4) cancel out between the deformed and undeformed geometries. It is, however, possible to ensure that the UV divergent terms in (2.4) cancel out even without such a compactification: as shown in [11], one can introduce a one-parameter family of large diffeomorphisms on flat space (that is, diffeomorphisms that don't vanish in the asymptotic region) in order to ensure that the differenced free energy is UV and IR finite.
Though the final result obtained in the flat-space limit of both the torus and the sphere is the same, we note the interesting difference that on the finite-size torus, the Dirac fermion has a negative mode for which ∆K L does not have fixed sign, and in fact even renders ∆F positive 10 ; as we now discuss, this is not the case for the sphere.

Negativity of ∆K
Our results imply that for any nontrivial deformation of the round sphere, σ∆K is strictly negative for all t to leading nontrivial order in , and hence so is a (t), and thus ∆F . This is easiest to see when f (1) contains only modes with odd : in this case, it is clear from the expressions (3.16) and (3.24) that σα , is negative when is odd and greater than one, and hence so too is σ∆K (2) (t) for all t (recall that σ = −1/2 for the scalar and σ = 1 for the fermion). When = 1, a (t) = 0 for both the scalar and the fermion, and hence ∆K (2) vanishes; this is due to the fact that = 1 deformations generate infinitesimal diffeomorphisms of the sphere and therefore do not change its intrinsic geometry to leading order in . We will show this explicitly in Section 4.2.
The case of even is more subtle. For the scalar, it follows from the full expression (A.18) that σα , can be positive for ≥ /2; likewise, for the fermion it follows from (A.29b) that σβ , is positive. Hence in both cases the sign of a (t) is not immediately clear. However, note that the large-behavior of a (t) can be obtained in the flat-space scaling limit discussed above, and is given in (3.36); assuming a r 0 k (t) is 10 Explicitly, consider the deformed torus , where x 1 and x 2 both have periodicity ∆x, and we take f = cos(2πx 1 /∆x) + O( 2 ). The perturbative heat kernel for the Dirac fermion on a deformed torus is computed in [14], and in this case comes out to be This expression is positive for t > t * for some t * , and with T = 0, M = 0 the differenced free energy (2.12) comes out positive as well.  The function a (t) for even for the minimally coupled scalar (left) and the Dirac fermion (right); from dark to light gray, the curves correspond to = 2 to = 40. The dashed red curve is the flat-space limit given by (3.36). The convergence of a (t) to the flat-space limit for the nonminimally coupled scalar is analogous.
continuous in k as r 0 → ∞, we therefore conclude that for all large (whether even or odd), σa (t) is negative. We therefore need only investigate the sign of σa (t) for small even (i.e. before the transition to the flat-space behavior). The result is shown in Figure 1, which verifies that σa (t) < 0 for all t.
Thus σ∆K (2) L (t) is indeed negative for all t. This implies, of course, that small, nontrivial deformations of the round sphere all lower the free energy of the scalar and of the fermion (for any mass, temperature, and curvature coupling), but it is in fact a much stronger result: negativity of ∆F does not require that the heat kernel ∆K L (t) itself be everywhere negative. We now investigate whether this stronger result continues to hold even for large area-preserving deformations of the sphere.

Nonperturbative Results
We have thus found that small perturbations of the round sphere always yield a negative σ∆K L (t) (and hence also free energy ∆F ) for both the scalar and fermion at any temperature, mass, or curvature coupling. This observation naturally prompts a question: is the free energy maximized globally by the round sphere, as it is for holographic CFTs at zero temperature [8]? Or do there exist sufficiently large deformations of the round sphere at which the free energy eventually increases above its value on the round sphere? If the free energy is globally maximized, does the stronger result that the differenced heat kernel has fixed sign continue to hold? Our purpose now is to examine these questions. To do so, we will ultimately need to resort to numerics in order to evaluate the heat kernel (3.1) for large deformations of the round sphere. However, we will first examine the behavior of ∆K L (t) at large and small t, which is tractable analytically even for large deformations.

Heat Kernel Asymptotics
Recall that the small-t behavior of the differenced heat kernel is given by the heat kernel expansion (2.14): The leading-order coefficient ∆b 4 is given by [24] ∆b (scal) 4 where R is the Ricci scalar of the round sphere. But it follows from volume preservation, the Gauss-Bonnet theorem, and the fact that R is constant that with equality if and only if R = R, i.e. if g ij is the metric of the round sphere. Hence for both the scalar and fermion, σ∆K L (t) is strictly negative at sufficiently small t for any nontrivial deformation of the sphere, regardless of the size of the deformation.
To inspect the large-t behavior, we instead recall that the differenced heat kernel can be expressed in terms of the eigenvalues λ I ,λ I of the operators L and L: The large-t behavior of this expression -particularly its sign -is clearly dominated by the smallest eigenvalue of either L or L, so we must compare the low-lying spectra of these two operators. For the scalar, this comparison can be performed by using a Rayleigh-Ritz formula for the lowest eigenvalue of L: with the infimum taken over all (square-integrable) test functions φ. Thus λ min can be bounded from above by taking φ to be a constant function; then again using the Gauss-Bonnet theorem and volume preservation, we have with equality if and only if a constant function is an eigenfunction of L, which for ξ = 0 is only the case if R is a constant and thus g ij is the metric of the round sphere. Hence for ξ = 0 and a nontrivial perturbation of the sphere, the lowest eigenvalue of L is always strictly less than any eigenvalue of L, and ∆K L (t) is positive at sufficiently large t. On the other hand, when ξ = 0 constant functions are always eigenfunctions of L = −∇ 2 , and hence the lowest eigenvalues of L and L are identical. The large-t behavior of ∆K L (t) is then controlled by the next-lowest eigenvalue λ next of −∇ 2 , which is known to be bounded by [25] λ next ≤ 8π Hence for the case ξ = 0, we again find that ∆K L (t) is positive at sufficiently large t.
We come to a similar conclusion for the fermion by invoking a theorem from [26]: namely, given that Σ is a two-dimensional manifold of genus zero, all eigenvalues of the squared Dirac operator are bounded below by 4π/Vol[g] =λ min , with equality only holding if g ij is the metric of the round sphere. Hence again we conclude that at sufficiently large t, σ∆K L (t) is negative.
We have therefore established that σ∆K L (t) is always negative at sufficiently small or large t, regardless the size of the perturbation to the sphere. To analyze the intermediate-t regime, and in particular to determine whether ∆F decreases arbitrarily as the size of the perturbation grows, we turn to numerics.

Numerical Results
The advantage of using heat kernels to evaluate the differenced free energy is that computing the (differenced) heat kernel (3.1) amounts to computing the spectrum of L. Moreover, only the smallest few eigenvalues of L are needed to obtain a good approximation for ∆K L (t) everywhere except near t = 0 -but small t is precisely the region in which the heat kernel expansion gives a good approximation. The heat kernel expansion therefore provides both a check of the numerics as well as a tractable way of computing the differenced free energy (which requires the behavior of ∆K L (t) to be known for all t). In short, we compute ∆K L (t) numerically to a sufficient accuracy that at sufficiently small t it agrees with the leading linear behavior ∆b 4 t of the heat kernel expansion, and we then sew these two behaviors together to perform the integration over all t that gives ∆F . We present more information on the numerical method used, as well as details of these checks, numerical errors, and computation of ∆F , in Appendix C. Here we instead describe the setup and the results.
First, note that on sufficiently deformed backgrounds the Ricci scalar will become negative somewhere, and hence the spectrum of L for the non-minimally coupled scalar may become negative. If these eigenvalues are sufficiently larger in magnitude than M 2 (as will always occur if M is fixed and the sphere is deformed more and more extremely), their presence introduces tachyonic instabilities, implying that the theory becomes illdefined. Consequently, we will restrict to numerical analysis of only the minimally coupled scalar ξ = 0, which as we showed in the previous section always has a nonnegative spectrum. No restriction is required on the fermion, since as mentioned above the spectrum of L for the fermion is always positive.
A numerical analysis can only by used to study a specific subset of deformations of the round sphere. Here we will consider certain classes of axisymmetric deformations. We begin by considering deformed spheres embedded in R 3 via r = R(θ), corresponding to the induced metric 11 Specifically, we will take where c , is a (positive) constant that ensures the volume of the sphere remains unchanged as is varied. It is straightforward to see that to linear order in , the metric (4.8) obtained from these embedding functions is in the form (3.2) conformal to the round sphere, and hence the behavior of ∆K L (t) to leading nontrivial order in should be the same as that obtained in Section 3 (with f (1) = Y ,0 ). However, higherorder effects in break the conformal form of the metric. We will consider deformations (4.9) with = 1, 2, . . . , 6, while the range of ∈ ( min , max ) is fixed by the condition that R , > 0 everywhere; we show cross-sections of the embeddings of these surfaces into R 3 in Figure 2. Note that for odd it suffices to consider only > 0, since positive and negative are related by a parity transformation: the transformation (θ, φ) → (π − θ, φ + π) sends Y ,0 → (−1) Y ,0 = −Y ,0 , and thus since Y ,0 always appears with a factor of , we find that R , → R − , for odd . It is also worth noting that the = 1 embedding doesn't appear to change the shape of the sphere much at all until is relatively large; as mentioned above, this is because the = 1 deformation is an infinitesimal diffeomorphism, and thus the deformation of the intrinsic geometry is trivial to linear order in . This can be seen explicitly by noting that since Y 1,0 (θ) = p cos θ (with p = 3/4π), the induced metric (4.8) with R = R 1, becomes  Figure 2. Cross-sections of the geometries we consider; these should be rotated around the dotted axis to generate the corresponding deformed sphere. The dotted blue circle is the unperturbed sphere; from light to dark gray, each curve corresponds to ranging in steps of 0.1 max (or 0.1 min ) from 0.1 max (0.1 min ) to 0.9 max (0.9 min ). For odd , negative is related to positive by a parity transformation which turns the cross-section "upside-down", leaving the geometry unchanged.
converting to a new coordinate ϑ defined by θ = ϑ − p sin ϑ + p 2 2 sin ϑ cos ϑ + O( 3 ), we get so the induced metric to linear order in is diffeomorphic to the round sphere, as claimed. The nontrivial perturbation comes in at order 2 and takes the form of those considered in Section 3 with f (1) = −Y 2,0 / √ 20π; the differenced heat kernel should thus be O( 4 ).
In Figures 3 and 4, we show the differenced heat kernels ∆K L (t) for the minimallycoupled scalar and for the Dirac fermion normalized by 2 (or by 4 in the case of = 1) along with the perturbative results derived in Section 3. Note that we only plot ∆K L down to t = 0.0005; this is because in the small-t regime more and more eigenvalues of L contribute to ∆K L leading to difficulty in controlling the numerics. But as discussed above, the small-t regime is controlled by the heat kernel expansion, which guarantees the sign of σ∆K L to be negative there. We therefore see that σ∆K L (t) is negative for all t even for large deformations of the sphere. Interestingly, ∆K L (t) appears to grow with at sufficiently small fixed values of t; this is due to the fact that as the geometry becomes more singular, its Ricci curvature grows, causing the heat kernel coefficient ∆b 4 defined in (4.2) to grow as well. This growth is especially pronounced in the deformations with odd and those with even and < 0; comparing to Figure 2, these deformations all limit towards a connected geometry with a cusp-like defect (the geometries with even and > 0, on the other hand, pinch off into separate disconnected components as → max ). This growth of ∆K L at small t should lead to a corresponding growth in the free energy ∆F ; we now investigate this free energy, and then more carefully investigate the divergence structure associated to the limiting singular geometries.

Behavior of the Free Energy
At zero mass and temperature, the differenced free energy ∆F T =0 = ∆E may be computed by using (3.26) and then integrating the heat kernel with (2.12); more details on the computation can be found in Appendix C. In Figure 5 we show ∆E for the deformations described above, normalized by the perturbative result ∆E pert . As expected, |∆E| grows monotonically with increasing | |, though this may not be apparent from Figure 5 as the curves are normalized by a factor of 2 (or 4 for = 1) contained in ∆E pert . Due to the growth of ∆K at small t as approaches min or max , we only show ∆E for a range of within which the error in ∆E is no greater than a few percent (this corresponds to up to 0.8 max for = 1 and up to 0.5 max for = 6). Nevertheless,  the growth in the small-t behavior of the heat kernel makes clear that ∆E should continue to grow as the geometry is successively deformed; we will investigate this growth in more detail in the following Section. For now, let us note the remarkable feature that ∆E/∆E pert looks extremely similar for both the scalar and fermion, despite the fact that the corresponding heat kernels in Figures 3 and 4 are more substantially different. It therefore appears that the theory-dependence of ∆E is contained almost completely in the perturbative contribution ∆E pert : the ratio ∆E/∆E pert is almost entirely theory-independent (we highlight almost: the difference between the curves is larger than numerical error, so they are genuinely different). This is feature is interestingly reminiscent of the results of [22], which study the free energy of the massless  Dirac fermion, the conformally coupled scalar, and holographic CFTs on a squashed Euclidean three-sphere; they found that for small and modest squashings, the free energies of all of these theories agree more closely than should be expected from CFT considerations alone.
In fact, this theory-independence becomes exact in a long-wavelength limit. Specifically, let l be the typical curvature scale of the deformed geometry; then the heat kernel coefficient ∆b 2n scales like l −2n for n ≥ 2. The heat kernel expansion (2.14) then implies that for l T −1 , M −1 , the free energy (2.12) can be expressed as an expansion where J (n) is the nth derivative of the function given by , fermion.
which is an expansion in 1/(M l) 2 . On the other hand, for M T we have which is an expansion in 1/(M l) 2 for the scalar and 1/(T l) 2 for the fermion. The point is that as long as l T −1 , M −1 , the leading-order behavior of the differenced free energy is governed by the lowest heat kernel coefficient ∆b 4 : where · · · denotes subleading terms. The theory-dependence of ∆F can be seen by expanding ∆b 4 = ∆b is the same for the fermion and the scalar, so to leading order ∆F/∆F pert is independent of the theory (as well as of the mass and temperature).
At intermediate masses and temperatures, ∆F interpolates between the massless zero-temperature behavior shown in Figure 5 and the behavior given by (4.16). As a representative example, we show this interpolation in Figure 6 for the case of the fermion and the deformed spheres (4.9) with = 2 (results for the scalar and higher are analogous). The takeaway is that for any mass and temperature, large deformations of the sphere appear to decrease ∆F arbitrarily. The deformations considered here tend to "pinch off" the sphere somewhere, and hence to better understand the behavior of ∆F under such extreme deformations, we now examine more closely the behavior of the heat kernel near these transitions.

Towards Singular Geometries
As remarked above, the deformations shown in Figure 2 fall into roughly two classes: the left two columns (corresponding to odd and even with < 0) limit to a connected geometry that "pinches" somehwere, while the geometries shown in the right column (corresponding to even with > 0) tend to disconnect as → max , with the individual connected pieces each potentially having a defect near the transition. In both classes, we expect the gradient of ∆K L to diverge at t = 0 as → min,max because the heat kernel coefficient ∆b 4 diverges as the geometry becomes singular due to the Ricci scalar becoming unbounded near the pinchoff 13 . However, the behavior of ∆K L at small nonzero t differs between these two classes. The class with even and > 0 is perhaps most intuitive: the case = 2 looks like a change in topology from one sphere to two, while for ≥ 4 the singular geometry also exhibits conical defects near the transition (in addition to the divergence of the Ricci scalar there). An isolated conical defect (with no curvature singularity) can be studied analytically, so we begin with a discussion of the associated divergences.

Conical Defects
In the vicinity of a conical defect on some manifold Σ, the geometry takes the form where φ has periodicity α (with α = 2π corresponding to a smooth geometry). Recall that a conical deficit (corresponding to α < 2π) can be embedded in R 3 , while an excess (corresponding to α > 2π) cannot. The differenced free energy, of course, depends only on the intrinsic geometry, so we may still analyze its behavior regardless of the in this expansion cannot be given by integrals of successively higher-derivative curvature invariants, because these diverge.
existence of any embedding. In the presence of such a defect, the corresponding heat kernel expansion exhibits an additional constant term associated to it [27]: The differenced heat kernel thus satisfies which clearly leads to a UV divergence in ∆F . Importantly, note that this divergence has fixed sign: it always contributes negatively to σ∆K L (t), and hence to ∆F . Interestingly, for a cone (that is, the geometry (5.1) with vanishing subleading corrections), the sign of the divergence of the energy depends on whether the defect corresponds to a conical excess or deficit. For example, in the case of a conformally coupled scalar at zero temperature, the energy density of a cone is [28] ρ ≡ T 00 = G(α) r 3 , where G(α) < 0 for α < 2π and G(α) > 0 for α > 2π. Hence the differenced free energy between a cone and a planar geometry with no conical defect is negatively UVdivergent 14 when α < 2π, and positively divergent when α > 2π. One might have naïvely expected the behavior (5.4) to have been universal near conical defects (at least for QFTs with UV fixed points, which are CFTs in the UV), but the heat kernel expansion (5.2) shows that the behavior of the stress tensor near such defects must be sensitive to the global properties of (Σ, g) (and in particular, if Σ is compact, it follows from (5.3) that the difference ∆F is always negatively UV-divergent, whether the defect is an excess or a deficit).
To manifestly illustrate such deformations, as well as to connect to the deformations considered in Section 4, consider a one-parameter family of spatial geometries that interpolates from a smooth geometry to one with a conical defect at a pole. An explicit axisymmetric example of such a family is given by the embedding where c is a volume-preserving constant that fixes the volume to 4π, i.e. the volume of the round unit sphere. For any < 1, this geometry is everywhere smooth, and for = 1, it exhibits a conical defect at θ = 0 with angle α = 2π/ √ 3 and a Ricci scalar 14 We are ignoring potential IR divergences associated with the fact that a cone is not compact.
which is bounded everywhere excluding the defect; see Figure 7. For < 1 we may therefore numerically compute the heat kernel as described in the previous section; we show these in Figure 8. As expected, the differenced heat kernel vanishes linearly at small t for any < 1, but its gradient there diverges as → 1. In the limit → 1, the heat kernel clearly approaches a function that goes to a nonzero value at t = 0 consistent with the expectation from (5.3): with the right-hand side just the α = 2π/ √ 3 case of (5.3). We can investigate a conical excess analogously by specifying the deformed sphere geometry directly rather than considering an embedding. To that end, consider the family of deformed spheres given by where c is again a volume-preserving constant. For < 1, these geometries are smooth, while the = 1 geometry exhibits a conical defect of angle α at the pole θ = 0 and a Ricci scalar which is bounded everywhere excluding this pole. The small-t behavior of the differenced heat kernels for α = 3π is shown in Figure 9; note that in the limit → 1, these too approach the t = 0 value expected from (5.3). Morevoer, one again finds σ∆K L appears positive for all t. In particular, these results confirm that on a topological sphere with a conical defect, both a deficit and an excess contribute negatively to the free energy, in constract with the expectation from (5.4) for planar geometries.

Even , > 0
Let us now return to the case of the deformed spheres shown in the right-hand column of Figure 2. As a representative example, in Figure 10 we show the small-t behavior of ∆K for the = 2 deformation (4.9) near = max . As expected, the heat kernel always vanishes linearly at t = 0 for any < max but its gradient there diverges as → max . More interestingly, ∆K appears to stabilize to a function that approaches a finite nonzero value at t = 0. This behavior is quite evident in the case of the fermion, though it is a bit less obvious for the scalar as the successive change in ∆K L appears to grow with each successive step in . We might hope to understand this behavior by using the heat kernel expansion, but as mentioned above, for = max this expansion breaks down due to the unbounded Ricci scalar near the pinchoff point. We therefore should not expect the expansion (5.2)  to capture quantitative details of the small-t behavior near the transition. However, it is interesting to note that it does capture some qualitative features; for instance, the differenced heat kernel appears to be approaching a function that limits to a nonzero constant at t = 0, similarly to the hear kernels shown in Figures 8 and 9. Moreover, note that the = 2, = max geometry does not have a conical defect and can be thought of as a transition from one to two topological spheres. This transition doubles the Euler characteristic from χ = 2 to χ = 4, which from the heat kernel expansion  conclusively because the scalar heat kernel does not appear to be growing linearly in near = max ).

Odd and even , < 0
For odd and even with < 0, the limiting geometry instead has a cusp. The corresponding behavior of the differenced heat kernels (for = 2) is shown in Figure 11; note that now the heat kernel itself, rather than just its gradient, appears to grow at small t as the geometry becomes singular. As shown in Figure 12, at intermediate values of t this growth appears to go roughly like t −1/2 , but does not appear to be maintained to arbitrarily small t. Indeed, the difficulty in inferring the limiting small-t behavior is presumably due to the breakdown of the heat kernel expansion in the singular limit -that is, it is unclear whether or not √ t ∆K L vanishes at t = 0 in the singular limit, and therefore whether ∆K L actually approaches a finite nonzero constant at t = 0 like it does for the cone or whether ∆K L genuinely diverges there. Nevertheless, we note that the scaling as t −1/2 is interesting as such a scaling of the heat kernel is expected on manifolds with boundary [24]. This behavior suggests that perhaps the cusp can be interpreted as a sort of boundary.

Implications for Graphene-Like Materials
The fact that ∆K L (t) approaches a nonzero constant at small t could have interesting consequences for materials like graphene, which as discussed in Section 1 may exhibit a competition between a classical membrane free energy ∆F c and the contribution ∆F q from effective QFT degrees of freedom. Indeed, note that (2.12) implies that on a  Figure 12. The same as Figure 11, normalized by a factor of √ t. At small-but-not-toosmall t, ∆K L appears to go like t −1/2 . It is unclear what happens at much smaller t due to lack of numerical precision there.
geometry with a conical defect, ∆F q has a linear UV divergence: where ρ is a short-distance cutoff that resolves the conical singularity (imposed by restricting to t > ρ 2 ) and B is a positive constant. On the other hand, the Landau free energy ∆F c given in (1.2) merely has a logarithmic divergence, which is due to its scale-invariance: near θ = 0 the mean curvature of the embedding (5.5) with = 1 diverges as θ −1 , and hence with B a positive constant. Interpreting (1 − ) ∼ ρ as the resolution parameter of the cone, it therefore follows that ∆F c must grow more slowly with than ∆F q , and hence a deformation with sufficiently close to 1 will have ∆F c + ∆F q < 0. This argument will fail, of course, once is so close to 1 that UV effects from the "tip" of the cone -which presumably are how the divergences in (5.9) and (5.10) are resolved -change the relative growth of ∆F c and ∆F q with . On the other hand, per the analysis of Section 1.1, for graphene we have that at small , ∆F c + ∆F q > 0. So although competition between ∆F c and ∆F q does not render the round sphere locally unstable, it appears that sufficiently large deformations of the round sphere may be preferred to the round sphere itself, even after accounting for the Landau free energy of the membrane. Whether this is actually the case will depend on the details of when our analysis breaks down.

Conclusion
We have provided evidence that for a (minimally or nonminimally coupled) free scalar field and for the Dirac fermion living on R × Σ, with Σ a two-dimensional manifold with sphere topology and endowed with metric g, the free energy is maximized when g is the metric of the round sphere. This observation applies to any mass and at any temperature. We demonstrated this result perturbatively around the round sphere for any nontrivial perturbation to the geometry, while for nonperturbative deformations we focused on a class of axisymmetric deformations. We found, in fact, not just that the free energy difference ∆F is negative, but that the (differenced) heat kernel ∆K L (t) itself has fixed sign -a much stronger result than merely negativity of ∆F . We have also shown that the free energy difference ∆F between an arbitrary g and the round sphere metric is unbounded below, diverging as g develops a conical defect. This property implies that any dynamics of the membrane driven by this free energy will tend to drive the membrane to a singular geometry (which presumably gets regulated by UV effects).
As an application of our results, we have also briefly investigated their relevance to (2 + 1)-dimensional crystalline systems like graphene, in which there are several contributions to the free energy. Specifically, in a Born-Oppenheimer approximation, the free energy we have calculated is that of the low-energy effective field theory of quantum excitations propagating on a fixed background determined by the atomic lattice; it simply corresponds to studying QFT on a curved background. The contribution to the free energy from this background is governed by a classical Landau free energy, and in this approximation it is the sum of these two that gives the total free energy of the membrane configuration. Consequently, understanding whether the negative free energy of QFTs on the membrane is sufficient to render the geometry singular depends on how well this quantum effect can compete with the Landau free energy of the underlying lattice. In the case of graphene at sufficiently low temperatures, we performed an analysis of this competition and found that with the simplifying assumption of a diffeomorphism invariant Landau free energy, the membrane energy dominates for small perturbations of the round sphere, rendering the round sphere stable. However, in Section 1.1 we provided a more general diagnostic (1.12) for when the QFT free energy can make the round sphere unstable. This constraint depends on parameters like the bending rigidity, lattice spacing, and details of the effective QFT fields living on the membrane. Since these parameters are presumably experimentally tunable from one type of crystalline membrane to another, it is plausible that a system can be engineered in which the QFT free energy does dominate that of the Landau free energy, and hence one should be able to experimentally observe the preference of such a membrane to deform to a singular effective geometry. Even without such engineering, the fact that integrating out the QFT fields gives a non-local contribution to the free energy (as opposed to the classical contribution from the geometric membrane action, which is local) suggests that perhaps the QFT contribution could be experimentally teased out from the classical piece, even if it never actually dominates the free energy.
More importantly, for large deformation of the sphere, we have shown that both ∆F q and ∆F c can be made arbitrarily large as the geometry becomes singular, with ∆F q negative and growing faster than ∆F c . Hence it is conceivable that even if the round sphere is locally stable, it is not globally stable, and large deformations are preferred. Verifying whether this is indeed the case requires understanding, for instance, how large must be in the one-parameter family (5.5) before our analysis breaks down due to UV effects.
Regardless of any competition with a Landau free energy, the results we have presented here prompt several questions, a few of which we would like to highlight. First, does the differenced heat kernel of any free field theory on a deformed sphere always have fixed sign? This is essentially a purely geometric inquiry, as one can define a heat kernel associated to any elliptic differential operator L. Presumably arguments like the ones we used in Section 4.1 to show that σ∆K L is negative at sufficiently small and large t could be used to gain control over the asymptotics of ∆K L for general L, but we do not know how to extend that analysis to intermediate values of t even in the case of the Dirac fermion and scalar studied here. Nevertheless, we conjecture that (under the condition that the spectrum of L is positive, to ensure that the free field theory defined by L is stable) the differenced heat kernel does indeed always have a fixed sign. Second, is the differenced free energy ∆F of any unitary, relativistic QFT on a deformed sphere negative? Note that here we include interacting field theories, for which a heat kernel cannot be defined. As shown in Table 1, this differenced free energy has been shown to be negative for holographic CFTs and (perturbatively) for general CFTs, which is tantalizing evidence that perhaps it is some universal feature of general QFTs. The generality of this picture suggests that there must be a universal underlying mechanism; it would be extremely interesting to uncover what this mechanism must be. Finally, how feasible would it be to engineer materials that actually exhibit this negative ∆F , either as a genuine instability of the round sphere, or merely as a non-local contribution to an effective description of a membrane's equilibrium dynamics?
the Simons Foundation (385602, AM). LW is supported by the STFC DTP research studentship grant ST/R504816/1. TW is supported by the STFC grant ST/P000762/1.

A Details on the Perturbative Results
In this Appendix we present the details on the perturbative calculation of the heat kernel for the scalar and the fermion.

A.1 Spin-Weighted Spherical Harmonics
In what follows, we will make use of the spin-weighted spherical harmonics s Y ,m . We refer to the original papers [29,30] for more details and explicit formulae; here we merely list the properties of these functions needed to make this Appendix self-contained. Essentially, a function η associated to a tensorial structure on the sphere is said to have spin weight s if under a local rotation of orthonormal frame by angle ψ, η transforms like η → e isψ η. Scalar fields of course are not tensorial and have spin weight zero, while the components of the Dirac spinor have spin weight 1/2. The spin-weighted spherical harmonics s Y ,m constitute an orthonormal basis for the space of spin weights functions on the sphere (hence the usual spherical harmonics are just the special case of spin weight zero: It will be convenient for later to introduce the spin weight raising and lowering operators ð andð, which act on a function η with spin weight s as ðη then has spin weight s + 1 andðη has spin weight s − 1. These operators obey the Leibnitz rule (even on products of functions of different spin weights) and hence are bona fide derivative operators, and are also total derivatives in the sense that when ðη orðη has spin weight zero, its integral over the round sphere vanishes. Moreover, they relate spin-weighted spherical harmonics of different spin weight to each other: (It then follows that the spin-weighted spherical harmonics with integer s can be generated from the ordinary spherical harmonics Y ,m by successive applications of ð.) The triple overlap of spin-weighted spherical harmonics can be expressed in terms of 3j symbols (see e.g. [31,32] 15 ): where here and in what follows we will leave the volume element sin θ dθ dφ in integrals implied. The 3j symbol

A.2 Scalar
For the non-minimally coupled scalar, the operator L was given in (3.13), which we repeat here: hence the operators L (n) introduced in (3.4a) are Decomposing f (1) in spherical harmonics as in (3.14), reality of f (1) requires that f * ,m = (−1) m f ,−m , while the volume preservation condition (3.3) requires that f 0,0 = 0. Then the matrix elements L (1) ,m, ,m are given by , vanishes, and hence so does the linear correction to the heat kernel: ∆K (1) = 0.
To compute the former, we use the addition theorem (A.1), and hence Tr L where in the first line the Laplacian ∇ 2 f (2) vanishes since it's a total divergence, in the second line we used the volume preservation condition (3.3) to replace the remaining f (2) , and in the final line we used (3.14). To compute Tr L , L , , we use (A.10): Tr L where we used (A.5d). Now, because the 3j symbols vanish unless the sum of the m quantum numbers is zero, m and m must be related by m = m + m 1 . We may therefore replace the phase (−1) m+m with (−1) m 1 , which we may combine with f 1 ,m 1 to give f * 1 ,−m 1 . Then using the orthogonality relation (A.5b) to evaluate the final sum, we obtain We may now insert (A.11c), (A.13), and (A.14) into (3.12) to obtain the secondorder correction to the heat kernel; one obtains the expression (3.15) given in the main text with The sum in the expression for α , can be simplified slightly by noting that by (A.5c), and hence one can write The sum over 1 in α , is a finite sum due to the triangle condition | − | ≤ 1 ≤ + ; by calculating this sum exactly for several values of , and using some sequencefinding functions in Mathematica, we are able to infer the closed-form expression where (x) n ≡ Γ(x + n)/Γ(x) are Pochhammer symbols and H n are harmonic numbers. Though we are unable to provide a general derivation, we have verified that this result agrees with (A.17) for all values of and from zero to 100. Finally, note that the 3j symbols in the expressions (A.15b) and (A.18) are only nonvanishing if + 2 is an even integer, implying that whenever is odd, β , vanishes for all and α , vanishes for all > /2. Thus the case of odd reproduces the expressions (3.16) given in the main text.

A.3 Dirac Fermion
For the fermion, we obtain the L (n) by expanding L given in (3.22). In fact, having introduced the spin weight raising and lowering operators ð,ð in (A.2) above, it is now natural to re-express L in terms of them using the fact that f has spin weight zero and L acts on the space of functions with spin weight 1/2. We ultimately obtain (we could also write ∇ 2 f = ððf and (∇ a f ) 2 = (ðf )(ðf ), but this rewriting will not be needed), and hence the unperturbed operator and the corrections L (n) are where to simplify L (2) we took f (2) to be a constant; there is no loss of generality in this simplification, since the purpose of f (2) is only to ensure that the volume preservation condition ,m, ,m by taking the unperturbed eigenfunctions to beh ,m = 1/2 Y ,m . Proceeding in this manner, first we obtain (again by expanding f (1) in spherical harmonics as in (3.14)) 16 where to slightly compactify notation we have defined with the second expression obtained from the first by using (A.6) and (A.5d).
As for the scalar, it follows immediately from (A.5a) and the fact that f 0,0 = 0 that Tr L (1) , = 0, and hence ∆K (1) = 0 as well. To get the second order term, we first compute and hence also Inserting these into the expression for the second-order correction to the heat kernel and using (A.24b), we thus obtain (3.23) with Note that as for the scalar, β , vanishes whenever is odd. The sum in α , is again a finite sum since the 3j symbol vanishes unless | − | ≤ ≤ + ; by evaluating the sum exactly for various values of , we are able to infer the expression

B Flat Space Scaling Limit
Here we provide some more details on the flat-space scaling limit performed in Section 3.4. First, to obtain the limiting behavior (3.36), note that for odd r 0 k we have from (3.35) where we have defined H(k, k ) = lim r 0 →∞ α r 0 k,r 0 k /r 3 0 and re-expressed the sum as an integral in the r 0 → ∞ limit. The function H(k, k ) can be obtained from the closed-form expressions (3.16) and (3.24) by using the fact that for x > 0, from which we find , fermion for even (k + k y )r 0 , and the same expression with a sine (instead of a cosine) for odd (k + k y )r 0 . Here N kr 0 ,kyr 0 is some new constant, which can be obtained up to an overall phase by the normalization condition for any = kr 0 and m = k y r 0 , where the domain of integration in k and k y comes from the restriction that > 0 and − ≤ m ≤ . Performing the integrals finally gives that up to an overall phase, for even (k + k y )r 0 and the same expression with a sine for odd (k + k y )r 0 . Hence in the r 0 → ∞ scaling limit, the coefficients f ,m give the desired expression (3.37): with the upper (lower) sign for even (odd) (k + k y )r 0 (and we are neglecting an overall phase that will cancel out). We remind the reader thatf (k x , k y ) is the Fourier transform of f (1) (x, y), and the assumption that f (1) (x, y) vanishes at large (x, y) (i.e. f (1) (θ, φ) vanishes away from (θ = π/2, φ = 0)) is what allows us to take the scaling limit of the integrand before evaluating the integral. Finally, decomposing ∆K (2) into contributions based on the parity of (k + k y )r 0 as odd , (B.9) we have ∆K (2) even,odd = t r 2 0 k,ky even,odd (k+ky)r 0 with the upper (lower) sign for ∆K (2) even (∆K (2) odd ). Note that the factor of 1/2 comes from the fact that for fixed kr 0 , we are only summing over k y r 0 with a given parity, and thus the spacing in k y is 2/r 0 . We now first switch the integrals around by taking the range of the k y integral to be (−∞, ∞) and the range of the k integral to be (|k y |, ∞), after which we change to a new variable k x = k 2 − k 2 y which has range (0, ∞). Since dk x = (k/k x )dk, we thus have with the latter expression obtained by expanding out the square and redefining k x → −k x as appropriate. Hence adding ∆K (2) even and ∆K (2) odd we obtain the flat-space expression (3.38).

C Numerical Method
In this Appendix, we describe the numerical methods used to compute the heat kernels and free energy in Section 4 and Section 5. We will always restrict to axisymmetric deformations of the sphere; almost all will be metrics of the form (4.8) obtained from an embedding r = R(θ), though we will also consider metrics of the form (5.7) to allow us to consider conical excesses (which cannot be embedded in R 3 ). In what follows, it will be convenient to work with the function f = ln R, so that the deformed metric (4.8) can be written as ds 2 = e 2f 1 + f (θ) 2 dθ 2 + sin 2 θ dφ 2 . (C.1) We will always use spherical coordinates {θ, φ} with ranges θ ∈ [0, π] and φ ∈ [0, 2π).
To compute the differenced heat kernel we use the form given in (3.1), Due to the axisymmetry, the eigenfunctions of L are separable and can thus be written as h(θ, φ) = w(θ)e imφ , where m takes integer values for the scalar and half-integer values for the fermion. We then have that L(we −imφ ) = e −imφ D m w, with D m a second-order ordinary differential operator given explicitly on the geometry (C.1) by where for the scalar while for the fermion The operator D m on the conical geometry (5.7) can be obtained analogously, so we do not explicitly write it here. We index the eigenvalues λ I as follows. Since λ I is an eigenvalue of L if and only if it is also an eigenvalue of D m for some allowed m, it is natural to take m to index the corresponding subspaces of eigenvalues. Within each subspace (that is, for each fixed m), we then introduce an integer l ≥ 0 to index the eigenvalues of the operator D m in ascending order. We therefore label the eigenvalues as λ m,l : for any given m, the λ m,l for l = 0, 1, . . . are all the eigenvalues of D m . With this notation, the eigenvalues of L (on the sphere with unit radius) are given bȳ λ m,l = (|m| + l)(|m| + l + 1) + 2ξ, scalar, in other words, the usual quantum number is replaced by |m| + l, enforcing that for a fixed m, ≥ |m|. The eigenvalues λ m,l are determined numerically by discretizing the operators D m over the interval θ ∈ [0, π] using standard pseudospectral differencing with a Chebyshev grid of N + 2 lattice points including the two boundary points θ = 0, π (see for example [33]). One subtlety is the appropriate treatment of the poles θ = 0 and π. Regularity of the metric (C.1) requires that f (0) = 0 = f (π), and hence an expansion of f around these points has a vanishing linear term. Using (C.3) for the scalar field, it then follows from a Frobenius expansion that near θ = 0 any solution to D m w = λ m,l w admits a regular behavior that goes like w = θ |m| w 0 + θ 2 w 2 + · · · (C.7) with w 0 = 0, in addition to a singular behavior that goes like θ −|m| . At the other pole, we have an analogous behavior: Hence when we difference the operators D m with m = 0, we impose Dirichlet boundary conditions at the poles, while when we difference the operator D 0 we impose Neumann boundary conditions. On the other hand, for the fermion we instead have the allowed behaviors Thus for |m| = 1/2 we discretize with Dirichlet boundary conditions at both poles, while for m = 1/2 we take a Neumann boundary condition at θ = 0 and Dirichlet at θ = π, and likewise for m = −1/2 a Neumann condition at θ = π and Dirichlet at θ = 0. Thus for each m we obtain an N × N matrix representing the discretization of D m ; for large N , the N eigenvalues of this matrix should approximate the eigenvalues λ m,l for sufficiently low l < N . Of course, a finite N will not be able to keep track of eigenvalues with m too large, so some cutoff on m must be imposed. A natural one is suggested by the spherical harmonics: on the round sphere, we might wish to keep all eigenvalues up to a fixed = |m|+l; since l < N , the strongest constraint is obtained by considering the lowest allowed |m|, which fixes a cutoff < N . Implementing this same cutoff procedure on the deformed sphere leads us to keeping all eigenvalues satisfying |m| + l < N : for each allowed (i.e. integer or half-integer) m with |m| < N , we compute the eigenvalues of the discretized operator D m and keep only the ones with l < N − |m|. The actual computation of the eigenvalues of the discretized D m is conveniently done with the Arnoldi algorithm which is implemented in the Mathematica matrix eigenvalue finder. We also note that for the minimally coupled scalar (i.e. ξ = 0), we explicitly drop the lowest eigenvalue m = 0, l = 0 because as discussed in Section 4.1 it is the same for both L and L and thus cancels exactly in the differenced heat kernel.
For a given N , a truncated differenced heat kernel can then be defined: m,l are the eigenvalues of the discretized operators D m , as described above. Increasing N should yield a better approximation to the exact heat kernel. We expect this approximation to be best at large t, in which the sum is dominated by the smallest eigenvalues, while the approximation should fail for sufficiently small t, when many eigenvalues make nontrivial contributions to the sum. Since the heat kernel time t can be thought of as an inverse square of a length scale, we expect that for a fixed N the agreement should fail for t smaller than order ∼ −2 max ∼ N −2 . However, the differenced free energy is sensitive to the small-t behavior of the differenced heat kernel. To accurately compute the free energy, we therefore implement a cutoff time t cut above which we integrate (2.12) with the truncated heat kernel ∆K (N ) L , and below which we integrate (2.12) using the leading-order behavior ∆b 4 t from the heat kernel expansion (2.14). For each N and choice of cutoff t cut , this gives an approximation to the free energy:  The accuracy of this approximation relies on ∆K (N ) L being well-approximated by the linear behavior ∆b 4 t around t = t cut , so that ∆F (N,tcut) is in fact independent of t cut . With our choice of N = 600, the truncated heat kernel ∆K (N ) L gives a good approximation down to t ∼ 2×10 −4 . For moderate deformations of the sphere (up to around = 0.5-0.7 for the deformations (4.9), depending on ), ∆K (N ) L agrees well with the leading-order behavior ∆b 4 t around this lowest value of t, and we may therefore compute the free energy as described. In this case, typically we take t cut = 2.5 × 10 −4 , and then varying t cut gives an estimate of the systematic error in ∆F (N,tcut) (for all plots in the main text, this error is no greater than a few percent). For larger deformations, however, ∆K (N ) L is not well-approximated by the leading-order behavior of the heat kernel expansion around t ∼ 2 × 10 −4 , and therefore we are unable to accurately compute the differenced free energy for such deformations.
We now discuss in more detail the convergence of ∆K (N ) L with N , agreement with the heat kernel expansion at small t, and agreement with the perturbative results for small deformations of the sphere.

C.1 Convergence
Since we are using pseudospectral differencing, we expect the error in a given eigenvalue to fall exponentially with N until a limit from machine precision is reached. Since the truncated differenced heat kernel is constructed directly from the eigenvalues, it too should converge to ∆K L exponentially until hitting machine precision. This convergence is fast at large t, since there ∆K L is sensitive to only the smallest eigenvalues, whereas at small t convergence (which is still exponential) requires larger N to achieve the same accuracy.
To exhibit this convergence, let us order the eigenvalues of L in ascending order; for a given resolution N , we then define the fractional error in the i th eigenvalue as Err (N ) where the maximum resolution we use is N max = 600. We plot this fractional error in Figure 13 for several eigenvalues in the geometry corresponding to the = 3, = 0.5 max embedding (4.9). This corresponds to a non-linearly deformed sphere, although one that still is not very close to being singular. We see that all the eigenvalues converge exponentially with N until reaching machine precision around N ∼ 100. As we would expect, it is the lower eigenvalues that suffer most from machine precision limitations in terms fractional error since they have a smaller absolute value (roughly the magnitude of the eigenvalues goes as λ i ∼ i). The data indicate that at the resolution N max = 600 used in this paper, the eigenvalues have a fractional error less than ∼ 10 −8 compared to their exact values. We may likewise define the fractional error in the differenced heat kernel as  .

(C.13)
This fractional error is shown in Figure 14 (in the same = 3, = 0.5 max geometry (4.9)) for several different values of t. Again, we observe initial exponential convergence before we become machine precision limited by N ∼ 80. Note that smaller t requires a larger N to reach the same accuracy, but the rate of convergence is roughly independent of t. We can estimate that for t > 0.05 the fractional error in the differenced heat kernel at N max = 600 is better than ∼ 10 −7 , which is commensurate with the error in the individual eigenvalues.

C.2 Comparison to Heat Kernel Expansion and to Perturbative Results
In addition to allowing us to compute the differenced free energy via (C.11) as described above, verifying that the heat kernel approaches the behavior predicted from the heat kernel expansion at small t also provides a check of our numerical methods. To that end, in Figure 15 we compare the small-t behavior of ∆K (N ) L for various N to the linear behavior ∆b 4 t expected from the heat kernel expansion; again we are taking the = 3, = 0.5 max embedding (4.9) as a typical example. There are two features to highlight. First, even the lowest value N = 40 recovers the heat kernel well above t ∼ 0.1, but computing the heat kernel accurately at very small t clearly requires using larger values of N . In particular, with the choice of N max = 600 used in this paper, we can reliably compute the heat kernel down to t ∼ 2 × 10 −4 (with some variation  for various N ; as in Figures 13 and 14, here we show the result for the geometry given by the embedding (4.9) with = 3 and / max = 0.5. From light to dark gray, the curves correspond to N = 40, 60, 100, 200, 400, while the dashed red line shows the linear behavior ∆b 4 t expected from the heat kernel expansion. Note that the linear behaviour only approximates the differenced heat kernel for quite small t, and we need to take N 400 to reach this linear regime.
depending on the deformation). Second, while the linear approximation ∆b 4 t does agree with the truncated heat kernel for sufficiently large N , for even moderate deformations of the sphere this agreement is only valid for very small t (for the case shown here, the fractional error between the linear behavior and the heat kernel is less than about two percent for t < 5 × 10 −4 , but grows much larger for larger t). For larger deformations this agreement moves to smaller and smaller t, eventually leaving the domain in which we can reliably approximate the exact heat kernel. It is for this reason that (C.11) cannot be used to approximate the differenced free energy for very large deformations.
As an additional check of our numerical method, we may compare the truncated heat kernel for very small deformations of the sphere to the perturbative heat kernels (3.15) and (3.23). We show this agreement in Figure 16, again for the = 3 embedding (4.9) but now only with a weak deformation of = 0.01. Even for modest N , ∆K (N ) L is very close to the perturbative result for reasonably large t. Increasing N gives agreement with the perturbative results to smaller t, as expected. We also show a comparison with the leading-order heat kernel expansion; unlike the moderate deformation / max = 0.5 shown in Figure 15, here we see good agreement with the expected linear behavior up to almost t ∼ 0.1.  Figure 15 we see that accuracy at small t requires larger N . We also show the leading linear behavior ∆b 4 t from the heat kernel expansion (dashed, red). order, J. Phys. France 48 (1987), no. 7 1085-1092.