{'fraction_non_alphanumeric': 0.03579040852575489, 'fraction_numerical': 0.008348134991119005, 'mean_word_length': 4.669939577039275, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Pinning of Diffusional Instabilities by Non-Uniform Curvature

Turing patterns emerge from a spatially uniform state following a linear instability driven by diffusion. Features of the eventual pattern (stabilized by non-linearities) are already present in the initial unstable modes. On a uniform flat surface or perfect sphere, the unstable modes and final patterns are degenerate, reflecting translational/rotational symmetry. This symmetry can be broken, e.g. by a bump on a flat substrate or by deforming a sphere. As the diffusion operator on a two dimensional manifold depends on the underlying curvature, the degeneracy of the initial unstable mode is similarly reduced. Different shapes can pin different modes. We adapt methods of conformal mapping and perturbation theory to analytically examine how bumps and ripples entrain modes of the diffusion operator on cylinders and spheres. We confirm these results numerically, and provide closed form expressions that describe how non-uniformities in curvature pin diffusion-driven instabilities and the resulting patterns.

Turing patterns emerge from a spatially uniform state following a linear instability driven by diffusion. Features of the eventual pattern (stabilized by non-linearities) are already present in the initial unstable modes. On a uniform flat surface or perfect sphere, the unstable modes and final patterns are degenerate, reflecting translational/rotational symmetry. This symmetry can be broken, e.g. by a bump on a flat substrate or by deforming a sphere. As the diffusion operator on a two dimensional manifold depends on the underlying curvature, the degeneracy of the initial unstable mode is similarly reduced. Different shapes can pin different modes. We adapt methods of conformal mapping and perturbation theory to analytically examine how bumps and ripples entrain modes of the diffusion operator on cylinders and spheres. We confirm these results numerically, and provide closed form expressions that describe how non-uniformities in curvature pin diffusion-driven instabilities and the resulting patterns.
In 1953, Turing predicted that spatial patterns can arise from a combination of diffusion and reactions of two chemical specices [1]. Such patterns are but one example of diffusion-driven instabilities, since studied on scales ranging from tens of nanometers on neurons [2], microns in chemical engineering [3], centimeters in zoology [4], to meters in ecology [5]. Substrate curvature plays a role in many of these systems. For example, axons and synapses of neurons are clearly not flat, and this curvature may influence the resulting patterns. Pattern formation has mostly been studied on flat surfaces, threedimensional reaction vessels, or surfaces of uniform curvature such as spheres and cylinders [6,7]. Studies on ellipsoids [8], where the Laplacian is again known in closed form, is a notable exception. The reaction diffusion equation has also been integrated numerically on generic surfaces [9,10]. However, to our knowledge, the intimate link between pinning of patterns and the spectrum of the Laplacian operator has not been pursued. Patterns, whether spots or stripes, break underlying translational/rotational symmetries, much like a crystal. Inhomogeneities in curvature, such as protrusions or cavities, reduce such symmetries and can pin or modify the patterns. To understand how nonuniform curvature can entrain and modify patterns, we study changes to the Laplacian, and its eigenmodes.
We follow a two prong strategy. First, we identify the onset of instabilities, by linearizing evolution equations expressed in terms of the appropriate 'modes', e.g., Fourier, cylindrical, or spherical harmonics. Modes with the largest (positive real part) growth are harbingers of the final patterns molded by nonlinearities. The modes are eigenfunctions of the diffusion (Laplacian) operator on the relevant manifold. Symmetries of the manifold, reflected in degeneracies of the eigenfunctions, must be broken in the final patterns. There is extensive literature on the behavior of the Laplacian on Riemannian manifolds, mostly focused on aggregate quantities, such as the determinant [11][12][13][14] appropriate to field theory. The short-time movement of particles due to perturbations of a manifold in the neighborhood of a given point has also been studied extensively [15][16][17]. We focus instead on how non-uniform curvature breaks the degeneracy of Laplacian modes, pinning eigenfunctions to inhomogeneities. Methods of mathematical physics are useful in this task, and we focus on conformal mappings and perturbation theory. There is no guarantee that patterns resulting from non-linear evolution are similarly entrained, so the second step in our study is to explore patterns with simulations. For this, we implemented the Thomas-Murray set of equations [4] on COMSOL [18].
As a first example, consider a cylinder with a localized axially symmetric deformation, described in cylindrical polar coordinates by ρ = R(z). This geometry can be simplified by a conformal mapping. The line element on the surface is where R z denotes the derivative with respect to z. Introducing a conformal axial coordinate v, we express the line element in the conformally flat form where R 0 is the asymptotic radius and Ω is the conformal factor. By parameterizing R = R 0 (1 + h(z)), we obtain the functional relationship between v and z, as If the deformation is symmetric and analytical at z = 0, the condition v ∼ z at z = 0 sets C = 0. The Laplacian of any scalar field φ is arXiv:1710.00103v3 [cond-mat.soft] 6 Oct 2017 where g ab is the metric of the underlying geometry and g is the determinant of the metric. In the conformal coordinates, √ g = Ω 2 R 0 , and it follows that the Laplacian on the bump geometry ∆ G , and in the conformally flat coordinates ∆ 0 , are related by Therefore, the eigenfunctions ψ k of ∆ G , with eigenvalue −k 2 , satisfy the equation similar to a Schrödinger equation, albeit with a potential proportional to k 2 . (This differs from da Costa's geometric potential, which comes from confining a particle to a surface [19].) For the undistorted cylinder of length L with periodic boundary conditions, the eigenfunctions and eigenvalues are with integers s and k. Now consider an axiallysymmetric Gaussian bump at the center of the cylinder, parametrized by standard deviation σ and height v 0 . To first order in v 0 R 0 ≡ and R 0 h z , the height function h is the same in both the initial and conformal coordinates. In the conformally-flat coordinates, the allowed values of k are modified to first order in as to satisfy the prescribed boundary conditions upon converting back to physical coordinates. Corrections to these eigenfunctions can be calculated perturbatively in using Eq. 6, which to first order reads with φ sk ≡ Z sk (v)e isϕ √ 2πR 0 . Note the sign of our potential: an attractive (repulsive) potential corresponds to an outward (inward) bump. We now apply Rayleigh-Schrödinger perturbation theory to calculate first order eigenvalue corrections, as The positive/negative sign is determined by the modecosine (symmetric) modes obtain a positive correction, whereas sine (antisymmetric) modes receive a negative correction. The sign of depends on the orientation of the Gaussian deformation-positive for a bulge and negative for a constriction. Thus, to first order, a Gaussian bump breaks the degeneracy in the eigenvalues of a Laplacian on the cylinder, causing the eigenvalues of sine (cosine) modes to become more negative in the case of an outward (inward) bump. This correction is largest for 2π k = 4πσ, or when the wavelength is approximately twice the width of the bump. This feature is a more general trait of deformations, explored later for a rippled cylinder. Of course, the bump can be positioned at any point along the cylinder, with corresponding shifts of cosine and sine modes to be centered on the Gaussian. The modified spectrum of the Laplacian affects any physical system involving diffusion. For example, consider a two-component reaction-diffusion system: After linearization, small deviations around a uniform stable fixed point of the reactions, (Ψ 1 , Ψ 2 ) * , evolve as with ω(k) satisfying the eigenvalue equation While the diffusion and the stable reaction matrices both posses negative eigenvalues, Turing showed that their sum can have positive eigenvalues (ω + (k) > 0) signaling finite-wavelength instabilities [1]. The (possibly degenerate) modes, k * , with the largest eigenvalue evolve to the final pattern. In the case of a cylinder with a Gaussian bump, the degeneracy between sine and cosine modes is broken, with one mode becoming more unstable than the other. This degeneracy breaking is dependent on the sign of h-switching from an outward bump to an inward bump changes which mode becomes more unstable. This linear analysis indicates only the onset of instabilities. To understand the patterns formed after nonlinearities stabilize the dynamics, we conducted finite element simulations of reaction-diffusion on curved substrates using the Thomas-Murray model [4]: We used COMSOL Multiphysics [18], which approximates the Laplacian by finite differences on a mesh and computes the fully non-linear reaction terms. (See supplemental material for animation and .mph file.) As shown in Fig. 1, density maxima of Turing patterns on a deformed cylinder are entrained by a Gaussian bump, switching between sine and cosine modes, depending on whether the deformation is inward or outward. This conformal mapping approach can be generalized to a variety of geometries, including a circularly symmetric deformation h(r) on a planar drum (see supplemental material), and an axially-symmetric bump on a sphere, where the spherical radius is given by R(θ). For the latter, the metric is given as where prime denotes derivatives with respect to θ. In conformally flat coordinates, relating the two coordinate systems as The Laplacian eigenvalue equation becomes where ∆ 0 is the Laplacian on the conformally round sphere. Setting R = R 0 (1 + h(θ)), allows us to expand in powers of . To O( ), θ = Θ, and Eq. (14) reduces to Equation (15) can be used to study how Laplacian eigenfunctions Y ,m (θ, φ) on a sphere are modified by deformations. Specifically, one sees that the second term, being only dependent on θ, breaks the full rotational symmetry of Laplacian eigenfunctions and fixes their polar orientation (although azimuthal rotational symmetry is preserved). This can have the effect of entraining Turing patterns to the deformations, although incommensurability of length scales and nonlinear interactions can compete against this, see Fig. 2. We also note that, similar to the deformed cylinder, eigenfunction modifications are sensitive to the sign of h-an inward bump causes an increase in diffusion, whereas an outward bump results in the opposite. Although non-linear effects on the development of Turing patterns on spheres have a tendency to stabilize spotted patterns instead of the exact spherical harmonic patterns predicted by the linear analysis [6], this feature can be observed in the initial development of Turing patterns, where an inward deformation causes the most diffusionally unstable mode to switch to one corresponding to a lower diffusion rate, see Fig. 2.
Multiple deformations can simply be represented by a superposition of "potential" terms in Eq. (15), each corresponding to the contributions of the individual deformations. To first order, eigenvalue corrections from the multiple potential terms also contribute additively. Identical deformations on each pole of a sphere, such as in Fig. 2, give twice the eigenvalue correction as an isolated bump to first order. In the initial development of Turing patterns, an undeformed sphere presents diffusionally-unstable modes corresponding to the eigenfunction Y 5 5 . Upon introducing inward deformations, the diffusion rate of lower harmonics is increased, thus causing the Y 0 3 pattern to become diffusionally unstable.
The above examples provide evidence as well as a mechanism for the pinning of patterns by curvature in-homegeneities. However, the conformal mappings employed are useful only for configurations of sufficient symmetry. More generally, consider a surface described by the embedding vector X(σ 1 , σ 2 ), and tangent vectors t a = ∂ a X. The induced metric is then given by g ab = t a ⋅ t b , and the Laplacian ∆ g is described by Eq. 4. Take any shape for which the eigenvalues, λ (0) k , and eigenfunctions, φ (0) k , of the Laplacian are known in closed form, e.g. a sphere, plane, or cylinder. Now deform the surface by a small amount X → X( ), where is a small bookkeeping parameter. This new surface will possess its own metric, G ab . If the deformation is sufficiently small, the chart σ a will continue to cover the surface X( ), and the Hilbert space of functions will be identical to that on X; in particular, the functions φ (0) k remain a basis for this Hilbert space. Any quantity defined on X( ) can be written in terms of its counterpart on X, plus small corrections proportional to , 2 , ⋯.
The Laplacian on this new shape is where G ab is the matrix inverse of G ab . This decomposition is exact. The ratio of volume measures, Ω, is equal to one plus corrections at each order in . The inner product correction, T ab , is non-perturbative and starts at 1 . To determine the perturbed eigenvalues, λ k , and eigenfunctions, φ k , of the the Laplacian ∆ G , we start with the ansatz: Using relations imposed by the constraint of orthogonality (see supplemental material), we find at first order: For a specific application, consider again a cylinder of length L with periodic boundary conditions. The embedding vector describing the cylinder is given by where φ ∈ [0, 2π). Similar to the conformal mapping approach, we describe longitudinal variations in the tube radius with a height function, R 0 (1 + h(z)). These modifications change the metric to To develop perturbation theory, we take to be smaller than 1 h and 1 h ′ . We also assume that the length and radius are of similar dimensions, L ∼ R 0 .
A rippled cylinder, such as in Fig. 3, is described by h(z) ∝ cos(2πpz L) ≡ cos(pz), where p is the number of ripples that fit along L. At first order, with the basis of eigenfunctions defined in Eq. (7), there is no first-order eigenvalue correction except near the boundaries of Brillouin zones,k =p 2 + δ, where non-degenerate perturbation theory breaks down. A pair of such boundaries exists for each value ofp. At each of these edges, degenerate perturbation theory yields first-order corrections analo-gous to those in a one-dimensional crystal [20]. Ifk ≥p 2, the eigenvalues are given by with momentum-dependent potential, h Fourier coefficients, and F ± : analoguous to Bloch's analysis. As in electron wavefunctions in a one-dimensional crystal, the perturbed eigenfunctions split cleanly into sine and cosine modes at the BZ edge, see Fig. 3. See animation of changing spectrum in supplemental materials rippled-cylinder-spectrum-sweep.mov.
Consider chemicals diffusing on the surface. Near the BZ edge, our analysis predicts an effective increase in the diffusive flux of the sine modes (troughs) and decrease in the cosine modes (ridges). Hence, diffusion is enhanced (diminished) where Gaussian curvature is negative (positive). This agrees with the short-time analysis of diffusion on a Riemannian manifold, in which the leading order correction to diffusion is proportional to Gaussian curvature [15][16][17]. In the context of Turing patterns, the steady-state regions of high concentration switch sharply from ridges to troughs as the most unstable wavelength is slowly dialed past twice the ripple wavelength, see Fig. 3.
The effects of curvature on diffusion in curved spaces may find application in heat flow [21], chemical diffusion [22], the local behavior of a soap film [23], fluid thin films where Marangoni flows result from surface tension gradients [24], in diffusion of electric and magnetic fields within the skin depth of conductor, and other systems. Recently, flocking systems have been studied on curved surfaces, suggesting a link between information transport of the flock and the underlying curvature [25]. Even in a system as simple as a cylinder with a bump, geometry can dramatically affect physical processes that depend on the Laplacian. While we discussed here patterns forming on a deformed manifold, it is interesting to consider membranes where the manifold itself deforms in response to patterns. JRF acknowledges support from the Hertz Foundation. JG acknowledges support from CONACyT grant 180901. MK acknowledges support from DMR-1708280. The technique of using conformal mapping to understand how deformations modify the Laplacian can be applied to a variety of surfaces. In a similar way to how we developed this approach in the case of deformed spheres and cylinders, we study the effects of a circularly symmetric deformation on a planar drum.
A planar drum with a circularly symmetric deformation can be described by a height function, h(r). The line element on this surface is given by Where h r denotes a derivative with respect to r. We assume a reparameterization of the surface geometry (v, ϕ) → (u, ϕ) so that the line element assumes the manifestly conformally flat form Setting these two line elements equal to each other, we obtain the two relations From these equations, we are able to calculate the relationship between u and r u = Cr exp Provided our bump possesses a tangent plane at r = 0, the relationship r ∼ u at r = 0 sets the constant of integration C = 1, and we thus obtain an exact relationship which can be expanded in powers of h r . To O h 2 r , where we have inverted the equation to obtain r(u). As u → ∞, this relationship becomes linear Eq. 3, we additionally obtain Ω, and hence the relationship between the Laplacian defined on the curved coordinates, ∆ G , and conformally flat coordinates, ∆ 0 .
Using Eq. 5 and the definition of Γ, we can express this in a more telling form where Similar to the cases of the deformed sphere and cylinder, we can express the effect of a deformation as the addition of an eigenvalue-dependent potential. As an example, if one considers a planar drum with a circularly symmetric Gaussian bump, Eq. 8 shows that this leads to a Gaussian potential in conformal coordinates.

OPERATOR PERTURBATION ON DEFORMED SURFACES
We wish to obtain a general relation between Laplacian eigenfunctions/values on some surface, given by φ 0 and λ 0 , and their counterparts φ and λ on a slightly deformed version of the same surface. As described in the main paper, the Laplacian on the original surface, parameterized by an embedding vector X(σ 1 , σ 2 ) and tangent vectors t a = ∂ a X, is given by where g ab = t a ⋅ t b and g = det g ab . On the new surface X( ), where characterizes the size of the deformation from the original surface, the Laplacian is where G ab is the matrix inverse of G ab , which is not the same as g ac G cd g db ≠ G ab , where g ab ≡ g (0) ab . This decomposition is exact. The ratio of volume measures, Ω, is equal to one plus corrections at each order in . The inner product correction, T ab , is non-perturbative and starts at 1 .
To determine the perturbed eigenvalues, λ k , and eigenfunctions, φ k , of the the Laplacian ∆ G , we start with the eigenvalue equation, ∆ G φ k = λ k φ k , multiply by an arbitrary wavefunction φ * j , and integrate over the entire manifold: The two terms in Eq. 10 are manifestly self-adjoint. We now adapt Raleigh-Schrödinger perturbation theory. Whereas this is usually applied to perturbations of operators on a fixed geometry involving a potential, here the perturbation is associated with changes in the geometry itself. We start with the ansatz: This structure itself is familiar in quantum mechanics on flat-space with c (n) jk ≡ ⟨j (0) k (n) ⟩. On the curved manifold, however, the inner product changes. This has many consequences. For example, the orthonormality is: The matrix elements of the dimensionless density ratio: In Eq. 12, each order in n (n > 0) must equal zero independently. Thus, shape changes rescale the eigenfunctions even at first order: Re c This differs from non-spatial perturbations, where one can always arrange perturbative corrections such that c (1) kk = 0. In either setting, one can set the imaginary part to zero -it is just a phase. These rescalings affect the eigenbasis of any normal operator in a perturbed geometry.
Inserting the perturbative expansions from Eq. 11 and similarly expanding T ab and Ω, we can group terms order by order in . After integrating by parts and using identities from Eq. 12, we obtain expressions superficially similar to perturbation theory on a Euclidean geometry. As presented in the paper: The difference of course is in the structure ofH (1) jk , which involves the self-adjoint operator −∇ a T ab ∂ b , as well as the addition of a term originating in the non-trivial measure, λ k Ω (1) jk . For completeness, the rippled cylinder's Laplacian: Matrix elements of the perturbed Laplacian density in the unperturbed eigenbasis are:

USING THE COMSOL FILES
To help others repeat our results, we have included several .mph files that you can run in COMSOL Multi-physics. Tuning the meshing and other simulation pa-rameters can affect convergence and interpretation; in particular, mesh irregularities at the ends of the cylinder models can act like defects and pin patterns. We cope with this by using smooth periodic boundary conditions and a rectangular mesh.
To reproduce the rippled cylinder results, you can: 1. Open periodic-cylinder-length-sweep.mph 2. Click the green equals sign to launch the solver, see Fig. 1.
3. The solver should finish after ten to twenty minutes of processing. You can watch progress in the log window or convergence plots. If it fails, try sim-ply relaunching it, because some versions of COM-SOL do not successfully initialize all of the needed components on the first run.
4. Click 3D Plot Group 5 to see solutions. It will default to showing the longest cylinder at the end of the sweep, see Fig. 2.
5. Select the length value, scroll up to a shorter cylinder, like 25, and click Plot, see Fig. 3.
6. The full sequence of length sweep solutions has exported as an animation, see rippled-cylinder-ridge-trough-switch.mov.