Strain and defects in oblique stripe growth
Kelly Chen, Zachary Deiman, Ryan Goh, Sally Jankovic, Arnd Scheel
aa r X i v : . [ m a t h . A P ] F e b Strain and defects in oblique stripe growth
Kelly Chen , Zachary Deiman , Ryan Goh , Sally Jankovic and Arnd Scheel Massachusetts Institute of Technology, Department of Mathematics, 182 Memorial Drive, Cambridge, MA 02139, USA University of Minnesota, School of Mathematics, 206 Church St. S.E., Minneapolis, MN 55455, USA Boston University, Department of Mathematics and Statistics, 111 Cummington Mall, Boston, MA 02215, USA
Abstract
We study stripe formation in two-dimensional systems under directional quenching in a phase-diffusion approx-imation including non-adiabatic boundary effects. We find stripe formation through simple traveling waves forall angles relative to the quenching line using an analytic continuation procedure. We also present comprehensiveanalytical asymptotic formulas in limiting cases of small and large angles as well as small and large quenchingrates. Of particular interest is a regime of small angle and slow quenching rate which is well described by theglide motion of a boundary dislocation along the quenching line. A delocalization bifurcation of this dislocationleads to a sharp decrease of strain created in the growth process at small angles. We complement our results withnumerical continuation reliant on a boundary-integral formulation. We also compare results in the phase-diffusionapproximation numerically to quenched stripe formation in an anisotropic Swift Hohenberg equation.
We investigate the influence of boundary conditions on the formation of striped patterns. Striped patterns occur inmany experimental setups [33, 9, 6, 34, 5, 32, 24, 35, 13, 1] and their existence and stability is quite well studied.In particular, idealized periodic striped patterns in unbounded, planar systems occur in families parameterized bythe wavenumber, the orientation, and a phase encoding translations. Stability depends only on the wavenumberand instability mechanisms include Eckhaus and zigzag instabilities. Away from instabilities, striped phases are welldescribed by a phase diffusion equation for a phase ϕ which encodes the (local) shift of a fixed reference pattern. Localwavenumbers and orientation are encoded in the gradient ∇ ϕ . Rigorous derivations are possible in a slow modulationapproximation [12]. In a homogeneously quenched pattern-forming system, posed with small noisy initial conditions,the observed pattern indeed locally resembles a suitably rotated and stretched periodic pattern, away from isolatedpoints or lines where defects form. More regular patterns emerge when the pattern-forming region expands in time,either through apical growth at the boundary of the domain, or through directional quenching where a parameter inthe system is changed spatio-temporally such that the parameter region where pattern formation is enabled growstemporally. Our interest here is with this growth scencario in an idealized situation.A prototypical model equation for the the formation of striped patterns is the Swift-Hohenberg equation u t = − (∆ x,y + 1) u + µu − u , ( x, y ) ∈ R , u ∈ R , (1.1)which possesses families of stable striped patterns for µ > u per ( kx ; k ) = u per ( kx + 2 π ; k ), close to p µ/ kx ) for small µ and k ∼
1. Directional quenching here refers to the situation where µ = − µ sign ( x − c x t )for some µ &
0. For patterns with trivial y -dependence and c x = 0, there exists a family of “quenched” periodicpatterns u with | u ( x ) − u per ( k x x − ϕ ; k x ) | → , x → + ∞ , | u ( x ) | → , x → −∞ , (1.2)for wavenumbers obeying the strain-displacement relation k x = g ( ϕ ) ∼ µ sin ϕ ; see [25, 31].For positive speeds c x >
0, one observes the formation of stripes with a selected wavenumber. This stripe formationis enabled by time-periodic solutions u ( t, x ) = u ∗ ( x − c x t, k x x ), with u ∗ ( ξ, ζ ) = u ∗ ( ξ, ζ + 2 π ) and u ∗ ( ξ, ζ ) → u per ( ζ ; k x ) , ξ → −∞ , u ∗ ( ξ, ζ ) → , ξ → + ∞ . These solutions represent stripes parallel to the quenching interface x = c x t , with trivial y -dependence. The wavenum-ber k x of stripes selected by this directional quenching process can be computed in terms of the strain-displacement1elation and effective diffusivities d eff as k x ∼ k min + k c / x + O ( c / x ) , k = − ζ (1 / p k min /d eff , where k min denotes the minimum of the strain-displacement relation; see [14].Including possible y -dependence, one would be interested in solutions that create periodic patterns at a given anglerelative to the quenching interface. This problem was analyzed in [2] when stripes are nearly perpendicular to thequenching interface and in [15] when stripes are almost parallel to the boundary for fixed c x >
0. Our focus here is onthe case of stripes almost parallel to the quenching interface and small speeds. Most of our results are concerned with aphase-diffusion approximation but we demonstrate numerically good agreement with Swift-Hohenberg computations.The phase-diffusion approximation for stripes relies on writing solutions u to (1.1) in the form u ( t, x ) = u per ( ϕ ; k ),with ϕ = ϕ ( x, y, t ) slowly varying and |∇ x,y ϕ | ∼
1, and ϕ t = ∆ ϕ, after possibly scaling x and y so that effective diffusivities agree. Of course, this assumes that the patterns consideredhere are away from possible instabilities, where for instance the Cross-Newell equations would be more appropriate.In a context of directional quenching, such an approximation is meaningful only in the pattern forming region x < c x t .The equation therefore needs to be supplemented at the quenching line x = c x t, y ∈ R , with an effective boundarycondition, which in particular should reflect the strain-displacement relation in the parallel case with c x = 0. Wethen arrive at ϕ t = ∆ ϕ + c x ϕ x , x < ϕ x = g ( ϕ ) , x = 0 , (1.3)where g reflects the strain-displacement relation, g ( ϕ ) = g ( ϕ + 2 π ) , g ( ϕ ) > , (1.4)for instance g ( ϕ ) = 1 + κ sin( ϕ ) for some 0 ≤ κ <
1. Clearly, setting ϕ = ϕ ∗ ( x ) and c x = 0, we find simple affineprofiles ϕ ∗ ( x ) = ϕ + g ( ϕ ) x, corresponding to the solutions in (1.2) compatible with the strain-displacement relation. Note that (1.3) possessesa gauge symmetry that maps solutions ϕ ( t, x ) to solutions ϕ ( t, x ) + 2 π , reflecting the periodicity of the underlyingperiodic pattern that is modulated through ϕ . It does not possess a continuous symmetry ϕ ( t, x ) ϕ ( t, x ) + ¯ ϕ ,¯ ϕ ∈ R , which would result in g ≡ const and reflect boundary conditions insensitive to the crystalline microstructure.This latter situation arises at leading order when one derives averaged amplitude or phase equations and one canthen think of the presence of a nontrivial flux g as a non-adiabatic effect, not visible in averaged approximations.The equation (1.3) was analyzed in [14] for y -independent solutions, deriving in particular universal asymptoticsfor solutions in the cases c x ≪ c x ≫
1. For c x ≪
1, excellent agreement with solutions in (1.1) and severalother prototypical examples of pattern-forming systems was found, including reaction-diffusion, Ginzburg-Landau,and Cahn-Hilliard equations. For bounded initial conditions and c x >
0, solutions eventually become time-periodicup to the gauge symmetry, and converge locally uniformly to linear profiles for large negative x , ϕ ( t + T, x ) = ϕ ( t, x ) + 2 π, | ϕ ( t, x ) − ( k x x − ωt ) | → , x → −∞ , ω = c x k x , for some T = πω >
0, for given g >
0. The existence and stability of such solutions with the minimal, 1:1-resonant period T = πω was established generally in [27]. Here, the resonance refers to the frequency of the periodicsolution 2 π/T relative to the frequency of patterns generated in the far field ω . In particular, subharmonic solutions2 πℓ/T = ω , ℓ >
1, are ruled out.In the two-dimensional, oblique case, these simplest resonant solutions correspond to traveling waves; see Figure1.1. In the far field, x → −∞ , we are interested in oblique stripes which are represented by values of the phase2 igure 1.1: Schematic plot of patterns with found through (1.5)–(1.8), with k y = O (1) (left), k y ≫ k y ≪ c y = − k x c x /k y ; seetext for details. Colors chosen to show contours of u = u per ( ϕ ( x, y, t )) with u per ( φ ) = sin( ϕ ); see Figure 1.3 for computed profiles. ϕ ∼ k x ( x + c x t ) + k y y = k x x + k y ( y − c y t ) with c y = − k x c x /k y . Such solutions are in fact traveling waves in the y -direction. We therefore focus on solutions ϕ ( x, k y ( y − c y t )) to (1.3), periodic up to the gauge symmetry in thesecond argument, that is, solutions to0 = ϕ xx + k y ϕ ζζ + c x ϕ x − k x c x ϕ ζ , x < , ζ ∈ R , (1.5)0 = ϕ ( x, ζ + 2 π ) − ϕ ( x, ζ ) − π, x ≤ , ζ ∈ R , (1.6)0 = ϕ x − g ( ϕ ) , x = 0 , ζ ∈ R , (1.7)0 = lim x →−∞ | ϕ ( x, ζ ) − ( k x x + ζ ) | , ζ ∈ R . (1.8)All solutions are in fact classical solutions since we shall assume g to be smooth. We will also see later that theconvergence in (1.8) is in fact uniform.In addition to ϕ , the system (1.5)–(1.8) includes 3 variables: the lateral periodicity k y , which we will assume to bepositive, without loss of generality; the quenching speed c x which we assume to be non-negative; and the strain k x ina direction perpendicular to the quenching line, which we think of as a Lagrange multiplier that compensates for thephase shift induced by ζ -translations. Given k x = k x ( c x , k y ), one can then determine angle and wavenumber fromthe wave vector ( k x , k y ).Our main results are as follows. Existence for all c x ≥ , k y > . Assuming g is smooth and 2 π -periodic, we have existence. Theorem 1 (Existence) . Suppose g > . Then for all c x ≥ , k y > , we have existence of solutions to (1.5) – (1.8) with k x = K x ( k y , c x ) , smooth. Moreover, solutions are strictly monotonically increasing in ζ . Using reflection symmetry, one can also find monotonically decreasing solutions. Solutions are unique within thisclass of solutions up to the trivial translation symmetry in ζ .We computed the function K x ( c x , k y ) numerically and show the resulting graph in Figure 1.2, using an appropriatecompactification of the positive quadrant c x , k y ≥
0. One sees quite distinct limiting behaviors of the surface andmuch of this paper is concerned with exploring these limits. Figure 1.2 includes a guide to the asymptotics and howthey are reflected in this surface.
Asymptotics c x → ∞ . Solutions ϕ and wavenumbers converge as c x → ∞ with limiting wavenumber K x ( c x = ∞ , k y ) independent of k y , given through the harmonic average of g . At finite but large c x , wavenumbers decreasefrom the harmonic average for small k y and increase for large k y , proportional to c − x at leading order. Asymptotics c x → , k y > fixed. Solutions and wavenumbers are smooth at c x = 0 with limit k x given by theaverage of g , and linear asymptotics for c x small. We establish asymptotics for the linear coefficient as k y → Asymptotics k y → , c x > fixed. Solutions are smooth (albeit likely not analytic) near k y = 0, c x >
0, aregime explored also in [15]. We numerically compute a leading-order quadratic coefficient and explore asymptoticsof this coefficient as c x → igure 1.2: Computed values of k x as a function of k y and c x in a compactified scale including the limits k y = ∞ and c x = ∞ . Surfaceplot (left; see § Asymptotics k y → ∞ . In this limit of perpendicular stripes, we find again the average of g as the limit andasymptotics with leading-order term k − y . Asymptotics k y ∼ c x → . In the most striking regime close to the origin, the sharp peak in the surface inFigure 1.2, we use an inner expansion to arrive at a reduced problem which amounts to describing the glide motionof a dislocation-type defect in the y -direction under an externally imposed strain. Most interestingly, we identifya qualitative “phase transition” where this defect changes type, explaining qualitatively the shape of the surface k x ( c x , k y ) close to the origin. Profiles of solutions in this regime on the boundary and in the whole domain are shownin Figure 1.3, demonstrating in particular the phase transition corresponding to the delocalization of a defect near k y /c x ∼ . § Numerical continuation.
We illustrate results and explore the approximation quality of theoretical asymptoticsusing numerical continuation for solutions of (1.5)–(1.8), and also for corresponding solutions of the Swift-Hohenbergequation. We find good agreement with asymptotics in the phase-diffusion equation, and a qualitatively similartransition near c x , k y ∼ Consequences for homogenized descriptions.
Thinking of the gradient of the phase as a macroscopic, ho-mogenized strain variable for a crystalline phase, our results provide corresponding effective boundary conditionsthrough a micropscopic analysis of the boundary layer. The dependence k x = K x ( k y ; c x ) provides mixed boundaryconditions, such that the renormalized strain φ = ϕ − K x ( k y ; c x ) x solves φ t = ∆ φ + c x φ x , x < , φ x = 0 , x = 0 , eliminating variations on the microscopic scale 1 /K x . Such a description is not possible for c x = k y = 0, sincethe derivative ϕ x at the boundary depends on the microscopic phase variable ϕ and, at steady-state, there aremultiple compatible equilibrium strain configurations. The presence of a spatial defect, k y = 0, or a temporal defect, c x = 0, forces selection of a unique normal strain at the boundary and allows this macroscopic description. From thisperspective, our work establishes existence of a unique normal strain and analyzes in detail properties of this normalstrain in various limiting regimes, in particular relying on properties of the spatio-temporal defect at the boundary.4 igure 1.3: Profiles of ϕ on x = 0 for k y = 2 . × − , k x = 0 . k y = 6 . × − , k x = 0 . c x = 10 − . Note the different scales on the horizontal axis, showing that the jump is stronger localized for larger k y . Associatedprofiles of sin( ϕ ) in the x − y -plane (only part of y -region shown), showing a sharply localized defect for larger k y (top right) and adelocalized defect for small k y (bottom right). Outline.
We introduce a boundary integral formulation together with a priori estimates and numerical setup in § k y = 0 , c x ≥
0, in §
3. We derive asymptotics in the limits c x → c x → ∞ , k y →
0, and k y → ∞ in §
4. We present an analysis near the origin k y , c x ∼ § § Acknowledgment.
KC, ZD, SJ, and AS gratefully acknowledge partial support through grant NSF DMS-1907391.RG was partially supported through NSF-DMS 2006887.
To solve (1.5),(1.6), and (1.8), we first set ψ ( x, ζ ) := ϕ ( x, ζ ) − ( k x x + ζ ) , (2.1)which gives 0 = ψ xx + k y ψ ζζ + c x ψ x − k x c x ψ ζ , x < , ζ ∈ R , (2.2)0 = ψ ( x, ζ + 2 π ) − ψ ( x, ζ ) , x ≤ , ζ ∈ R , (2.3)0 = ψ x − g ( ψ + ζ ) + k x , x = 0 , ζ ∈ R , (2.4)0 = lim x →−∞ ψ ( x, ζ ) , ζ ∈ R . (2.5)Next, writing Fourier series ψ ( x, ζ ) = P ℓ ∈ Z ψ ℓ ( x )e i ℓζ transforms (2.2) intod d x ψ ℓ + c x dd x ψ ℓ − k y ℓ ψ ℓ − k x c x i ℓψ ℓ = 0 , ψ ℓ ( x ) = X ± ψ ± ℓ e ν ℓ ± x , ν ℓ ± = − c x ± r c x k y ℓ + c x k x i ℓ, (2.6)where we use the standard cut at R − in the square root and restrict to c x ≥
0. For ℓ = 0, decay (2.5) requires ψ ℓ − = 0. For c x = ℓ = 0, solutions are affine, ψ ( x ) = ψ + ψ x , and we can set ψ = 0 since this part of the solutionis already parameterized by the ansatz (2.1) through the parameter k x . Evaluating ψ x at x = 0 and substitutinginto (2.4) then reduces (2.2)–(2.5) to the boundary-integral equation0 = D + ( ∂ ζ ; c x , k x , k y ) ψ − g ( ψ + ζ ) + k x , ψ ( ζ ) = ψ ( ζ + 2 π ) , D + (i ℓ ; c x , k x , k y ) = ν ℓ + , (2.7)5here the operator D + is understood as a Fourier multiplier acting through multiplication by ν ℓ + on Fourier series.One readily confirms that D + : H ⊂ L → L is a closed, sectorial operator as a relatively compact perturbationof k y | ∂ ζ | , with compact resolvent and spectrum with strictly positive real part except for the simple eigenvalue λ = 0associated with constant functions. The definition of D + extends to c x = 0 in natural agreement with our problem.For later puposes, we also introduce the associated pseudo-differential operator D − through D − (i ℓ ; c x , k x , k y ) = ν ℓ − . Lemma 2.1.
Suppose g is periodic and smooth. Then there exists a constant C ∞ ( g, c x , k y , m ) such that any solutionto (2.7) with ψ (0) ∈ [0 , π ) satisfies k ψ k C m + | k x | ≤ C ∞ . Moreover, C ∞ is uniformly bounded for fixed m and δ > such that | k y | > δ , k g k C m ≤ /δ . Proof.
Since − R D + ψ = 0 and | g | ∞ ≤ C g , we find an a priori bound | k x | ≤ − R | g ( ψ ( ζ ) + ζ ) | . This in turn gives an L ∞ a priori bound on D + ψ and, using the regularizing properties of D + and a bootstrap, the desired a priori boundon ψ . Uniformity of C ∞ follows readily from the fact that the pseudo-inverse of D + is uniformly bounded from L into H / as long as k y is outside a neighborhood of the origin. Numerical setup.
We solve (2.7) numerically for the variables ψ and k x , with parameters c x and k x , and adding aphase condition R ψ ( ζ ) exp( − ζ /δ )d ζ = 0. The resulting nonlinear equation is evaluated using fast Fourier transform.A Newton method, using gmres to solve the linear equation in each Newton step was found to converge robustlyeven for poor initial guesses. Most of the solutions were then computed using secant continuation in k y for fixed c x with adaptive control of the continuation step. During each step, we control for the number of Fourier modes byensuring that amplitudes in high Fourier modes is below a tolerance, which we found to have little effect once below10 − . Step sizes are very small and numbers of Fourier modes grow when c x , k y ∼
0, due to large gradients in theprofile. We address this regime directly using an inner expansion and a slightly different ansatz function in §
5. Thecode was implemented in matlab and Newton iterations for large sizes N ≥ were carried out on a Nvidia GV100GPU. All numerical results use g ( ϕ ) = 1 + κ sin( ϕ ) with κ = 0 . We prove Theorem 1. For this, we perform a homotopy, introducing g τ ( u ) := τ g ( u ) + (1 − τ ) − R g . Clearly, g τ satisfiesall the assumptions of Theorem 1 for τ ∈ [0 , g τ >
0. Let I ⊂ [0 ,
1] be the set of values where theconclusion of Theorem 1 holds. We will show below that(i) 0 ∈ I ; (ii) I is closed; (iii) I is open.Together, this implies that I = [0 ,
1] and establishes Theorem 1. This general strategy of proof was used in [27] for thecase k y = 0, although the proof there was based directly on the parabolic equation rather than the boundary-integralformulation which we shall exploit here.To show (i), we set k x = − R g and ψ = 0, such that ϕ is strictly monotone.To show (ii), take a sequence of solutions ψ n with wavenumbers k nx for converging values τ n → τ ∗ . We may assume,possibly adding multiples of 2 π , that ψ n (0) ∈ [0 , π ). By Lemma 2.1, we can assume that ψ n → ψ ∞ and k nx → k ∞ x ,possibly passing to a subsequence. The limit then solves (2.2)–(2.5). It remains to show that the limit ϕ ∞ = ψ ∞ + ζ is strictly monotone. Clearly, ψ ′∞ ≥ − ψ ′∞ ( ζ ) = −
1. Note that v = ψ ′∞ + 1 solves (2.2), (2.3), and (2.5), together with the linearized boundaryconditions 0 = v x − g ′ τ ∞ ( ψ ∞ + ζ ) v, x = 0 , ζ ∈ R , and has v ( ζ ) = 0, v ζ ( ζ ) = 0, v ζζ ( ζ ) ≥
0. Extending into x < v x ( ζ ) = 0and, using the equation, v xx ≤
0. On the other hand, since − R v > x = 0, v ( ζ, x ) → − R v | x =0 >
0, a constant. Since6nterior minima are excluded by the maximum principle, the minimum of v is necessarily located at the boundary x = 0 , ζ = ζ , which however implies v x ( ζ ) > ψ ∗ , L bi , ∗ v = D + ( k x ) v − g ′ ( ψ ∗ + ζ ) v, is Fredholm of index zero with ψ ′∗ ( ζ ) + 1 belonging to the kernel. We claim that the kernel is indeed one-dimensionaland that the derivative of (2.7) with respect to k x , D ′ + ( k x ) ψ ∗ , does not belong to the range. Together, this thenestablishes (iii) via the Implicit Function Theorem since the linearization with respect to ( ψ, k ) is onto. Suppose firstthat there is a function v in the kernel that is not a multiple of ψ ′ ( ζ ) + 1. Then we can find a linear combination thatis non-negative but not strictly positive, that is, a function w in the kernel with w ( ζ ) = 0, w ( ζ ) ≥
0, and − R w ≥ D + ( k x ) v − g ′ ( ψ ∗ + ζ ) v = −D ′ + ( k x ) ψ ∗ − , (3.1)where we suppressed the dependence of D + on its arguments other than k + , and the dependence of g on τ . Notethat in the case c x = 0, D ′ + = 0, D + is self-adjoint, with cokernel ψ ′∗ + 1, such that the right-hand side of (3.1) hasnonzero scalar product with the cokernel and hence does not belong to the range. We shall therefore assume in thesequel that c x >
0. The boundary integral equation (3.1) is equivalent to the elliptic equation0 = v xx + k y v ζζ + c x k x v ζ − c x v x , x > , (3.2)0 = v x − g ′ ( ψ + ζ ) v + D ′ + ( k x ) v + 1 , x = 0 . (3.3)We claim that the existence of a solution to (3.3) is equivalent to the existence of a generalized eigenvector inan associated elliptic problem, which will then lead to a contradiction. Consider therefore the eigenvalue problemassociated with our linearization 0 = v xx + k y v ζζ + c x k x v ζ − c x v x − λv, x > , (3.4)0 = v x − g ′ ( ψ + ζ ) v, x = 0 , (3.5)with solution v = ψ ′∗ + 1 at λ = 0. Existence of a generalized eigenvector then amounts to a solution v to0 = v xx + k y v ζζ + c x k x v ζ − c x v x + c x ( ψ ′∗ + 1) , x > , (3.6)0 = v x − g ′ ( ψ + ζ ) v, x = 0 , (3.7)or, setting v = w + x , 0 = w xx + k y w ζζ + c x k x w ζ − c x w x + c x ψ ′∗ , x > , (3.8)0 = w x − g ′ ( ψ + ζ ) w + 1 . x = 0 . (3.9)Solving the first equation using Fourier series in ζ and a variation-of-constant formula exploiting boundedness as x → ∞ , we find after a short calculation w x (0) = D + w (0) + ( D + − D − ) − c x ψ ′∗ | x =0 , which is equivalent to (3.3). This however contradicts the simplicity of the first eigenvalue of the elliptic operatordefined in (3.3). { k y > , c x > } We derive asymptotics in the regular and singular limits when either c x or k y tend to 0 or infinity.7 .1 The case c x = 0 In this case, we can multiply (2.7) by ψ ′∗ + 1 and integrate over ζ ∈ [0 , π ] to find0 = Z ζ (( ψ ′∗ + 1) D + ψ ∗ − ( ψ ′∗ + 1) g ( ψ ∗ + ζ ) + ( ψ ′∗ + 1) k x )= 2 π ( k x − − Z ϕ g ( ϕ )) , where we used that D + is a symmetric operator with kernel spanned by the constant functions to see that the firstsummand vanished, and monotonicity of ψ ∗ + ζ to transform the second summand into an integral over ϕ . As aconsequence k x = − R g is a priori known; see also [22, 2, 3], where this wavenumber selection mechanism was derivedfrom Hamiltonian identities. c x → We suppose that k y > c x →
0. Since the operator D + ( c x ) is continuous in the limit c x = 0 asa map from H into L , this limit is a regular perturbation problem. Using in addition that the linearization at a Figure 4.1: Left: Strain as a function of small growth rates comparing numerical continuation with theory (4.1), where the ˚ H / -normwas computed numerically. Right: Asymptotics for the ˚ H / -norm as k y →
0, comparing with theory (4.2), best fit for O (1) terms. profile, including the parameter k x as a variable, is onto, we conclude that we can formally expand the solution in c x , ψ ∗ ( ζ ; c x ) = ψ ∗ ( ζ ; 0) + c x ψ ( ζ ) + O ( c x ) , k x = k x, + c x k x, + O ( c x ) , k x, = − Z g. Inserting this expansion into the equation and taking the scalar product with the kernel of the linearization ψ ′∗ + 1gives at order c x , expanding D + = D + c x D + O ( c x ), D ( ℓ ) = ( − k x, k y sign( ℓ )), ℓ = 0, sign(0) = 0,0 = Z ζ (cid:0) ( ψ ′∗ + 1)( D ψ ∗ + k x, (cid:1) = 2 π (cid:18) k x, + k x, k y k ψ ∗ k H / (cid:19) , where we used R ψ ∗ = 0 and set k ψ ∗ k H / = − R ψ | ∂ ζ | ψ , which gives k x = k x, + (cid:18) − k x, k y k ψ ∗ k H / (cid:19) c x O ( c x ) , k x, = − Z g. (4.1)Formally setting k y = c x = 0, we find that D + = 0 and a solution ψ , ∗ = − ζ mod 2 π which does not belong to˚ H / . Writing the equation as ψ = (1 + D + ) − ψ + g − k x leads to the prediction ψ ∗ ∼ (1 + k y | ∂ ζ | ) − ψ , ∗ with k ψ ∗ k H / ∼ − k y ) + O k y (1) . (4.2)8n particular, we expect a strong initial stretching , that is, a decrease in k x with c x proportional to − c x | log( k y ) | /k y .Computed solutions k x are compared with the asymptotic prediction in Figure 4.1, where we also show agreementbetween the asymptotic prediction for the linear coefficient and the asymptotic formula (4.2). c x → ∞ We suppose that k y > c x → ∞ . We therefore set c x = ε − and formally expand D + ( ℓ ; ε ) = i k x ℓ + ( k x + k y ) ℓ ε + (2i ℓkx ( ℓ k x + ℓ k y )) ε + O ( ε ) . We start by considering the case ε = 0, where D + ( ∂ ζ ; 0) = k x ∂ ζ . As a consequence, at ε = 0, the solution ψ = ψ + ζ solves the ordinary differential equation k x ψ ,ζ = g ( ψ ) , ψ ( ζ + 2 π ) = ψ ( ζ ) + 2 π, (4.3)with implicit solution from separation of variables. In particular, the wavenumber at infinity is the harmonic averageof the nonlinearity, k x, = (cid:18) − Z ( g ( v )) − (cid:19) − . The linearization at ε = 0, ψ is L v = k x, v ζ − g ′ ( ψ ) v, which we consider as a Fredholm operator of index zero from H into L . The derivative of (4.3) with respect to k x is ψ ,ζ which does not belong to the range, so that the linearization is, as in the case of finite c x discussed in § ε , one needs to besomewhat careful. We therefore first expand formally, k x = k x, + k x, ε + k x, ε + O (3) , ψ = ψ + ψ ε + ψ ε + O (3) , where ψ j , j > L ψ + (cid:0) k x, ψ ,ζ − (cid:0) ( k x, ) + k y (cid:1) ψ ,ζζ (cid:1) = 0 . (4.4)Integrating against the adjoint kernel 1 /ψ ,ζ we see that k x, = 0 since, using the chain rule to compute ψ ,ζζ andchanging integration to ψ instead of ζ , Z π ψ ,ζζ ψ ,ζ d ζ = Z π g ′ ( ψ ) g ( ψ ) d ψ = 0 , by periodicity of log( g ( ψ )). We can then solve for ψ as ψ = ( k x, ) + k y k x, log( ψ ,ζ ) ψ ,ζ = ( k x, ) + k y ( k x, ) log (cid:18) g ( ψ ) k x, (cid:19) g ( ψ ) . (4.5)At order ε , we find L ψ + (cid:18) k x, ψ ,ζ − g ′′ ( ψ )( ψ ) + 2 k x, (cid:0) ( k x, ) + k y (cid:1) ψ ,ζζζ + k x, ψ ,ζ − k x, k x, ψ ,ζζ − (cid:0) ( k x, ) + k y (cid:1) ψ ,ζζ (cid:19) . (4.6)Using that k x, = 0, integrating against the kernel of the adjoint 1 /ψ ,ζ , and changing variables of integration gives k x = k x, + k x, c − x + O ( c − x ) , k x, = − Z (cid:26) g ′′ ( ψ )( ψ ) − k x, (cid:0) ( k x, ) + k y (cid:1) ψ ,ζζζ + (cid:0) ( k x, ) + k y (cid:1) ψ ,ζζ (cid:27) ψ ,ζ ) d ψ , (4.7)9here one substitutes ψ ,ζ = 1 k x, g ( ψ ) , ψ ,ζζ = 1( k x, ) g ′ ( ψ ) g ( ψ ) , ψ ,ζζζ = 1( k x, ) (cid:0) g ′′ ( ψ )( g ( ψ )) + ( g ′ ( ψ )) g ( ψ ) (cid:1) , and uses equation (4.5).The resulting integrals can be evaluated numerically for specific choices of g ( v ). We found that for g ( v ) = 1+ κ sin( v ), | κ | < k x, is monotonically increasing as a a function of k y , k x, < k y = 0 and 0 < k x, ∼ k y for k y large.More explicitly, the integrals can be evaluated to order κ for g ( ϕ ) = 1 + κ sin( ϕ ), yielding k x ( c x ) = p − κ + 12 (cid:18) κ ( − k y ) + 14 κ (3 + 5 k y ) + O ( κ ) (cid:19) c x + O ( c x ) . (4.8)This proves in particular that, at least for small κ , the monotonicity of k x as a function of c x changes, that is, k x, changes sign, to leading order at k y = 1.Figure 4.2 shows numerically computed values of k x compared with asymptotics for large c x , for several values of k y ,and demonstrates the sign change of the second-order coefficient k x, in a comparison with (4.8). Figure 4.2: Left: k x for large c x for several k y -values, compared with theory (4.7); inset shows comparison on of k x − k x, ∞ and c x on log-scales. Right: Leading-order coefficient k x, as a function of k y through numerical evaluation of (4.7), (solid), and explicit approximation(4.8), In order to make this expansion rigorous, we rewrite the equation as(1 − D + , ( ε, ζ )) D ψ − g ( ψ ) + k x = 0 . (4.9)The operator (1 −D + , ( ε, ζ )) is bounded invertible on L as a direct inspection of the Fourier symbol shows. Moreover,it is continuous at ε = 0 as an operator from H to L , again via a direct inspection of the Fourier symbol, withlimit the identity. Therefore, (4.9) can be written as F ( ψ, k x ) := D ψ − (1 − D + , ( ε, ζ )) − ( g ( ψ − k x ) = 0 , where F : H × R → L is continuous in ε at ε = 0. The Implicit Function Theorem then guarantees the existenceof solutions for ε >
0, small, with leading-order terms ψ , k x, . Substituting subsequently higher-order expansion,one can proceed in a similar fashion to establish validity of the expansion to any fixed order. k y → We follow the strategy from the previous section and find at O (2),0 = − Z ψ ad (cid:0) c x + 4 c x k x, ∂ ζ (cid:1) − / ( ψ ,ζζ − c x k x, ψ ,ζ ) d ζ, ψ ad is the (unique up to scalar multiples) periodic solution to the adjoint equation D + ( − ∂ ζ ) ψ − g ′ ( ψ ) ψ = 0.Unfortunately, the solution to the adjoint equation does not appear to be readily expressible in terms of ψ so thatwe will rely on numerical methods to evaluate the integral and obtain coefficients k x, and k x, in the expansion k x = k x, + k x, k y + O ( k y ) . (4.10)The numerically computed results shown in Figure 4.3 show good agreement up to a sharp transition value that weshall discuss in § k x, decreases with c x in a monotone fashion, converges to 0 as c x → ∞ and to ∞ for c x →
0, with power law asymptotics k x, ∼ c − βx , β ∼ /
2. Asymptotics are well capturedthrough k x, = c − / x ( c log( c x ) + c ); (4.11)fitting c and c for c x ∈ [5 · − , · − ] provides excellent agreement for a wide range of c x -values; see Figure 4.3.We did not attempt to justify asymptotics but provide a conceptual explanation in § k y → ∞ Expanding in inverse powers ε = 1 /k y , we find formally O ( −
1) : | ∂ ζ | ψ = 0 (4.12) O (0) : | ∂ ζ | ψ + k x, − c x ψ − g ( ψ + ζ ) = 0 , (4.13) O (1) : | ∂ ζ | ψ − c x ψ + 18 | ∂ ζ | − (cid:0) c x + 4 c x k x, ∂ ζ (cid:1) ψ + ( k x, − g ′ ( ψ + ζ ) ψ ) = 0 . (4.14)At order -1, (4.12), we set ψ = 0, which gives at order 0, in (4.13), k x, = − Z g, ψ = | ∂ ζ | − ( g − − Z g ) . Substituting the result into the equation at order 1, (4.14), yields | ∂ ζ | ψ − c x ψ + ( k x, − g ′ ( ζ ) ψ ) = 0 , which upon averaging gives k x, = − Z g ′ ( ζ ) | ∂ ζ | − ( g ( ζ ) − − Z g ) = 0 , which can be readily seen upon expanding g in Fourier series, and ψ = | ∂ ζ | − (cid:18) ( g ′ + 11 c x ) ψ (cid:19) . Assuming that g ′ is even, for instance g = 1 + κ sin( v ), we see that ψ and ψ are both odd. At the next order, wefind k x, = − Z (cid:18) ( − c x − g ′ ( ζ )) ψ − g ′′ ( ζ )( ψ ) (cid:19) , which vanishes when g ′ is even. Continuing further the expansion, we find that the even part of ψ is nonzero, ψ , e = | ∂ ζ | − (cid:18) − c x k x, g ′ ( ζ ) (cid:19) , and therefore k x, = − Z ψ g ′ = 0 .
11n the specific case g ( v ) = 1 + κ sin( v ), we find k x = 1 + k x, k y + O ( k y ) , k x, = − c x κ ; (4.15)see Figure 4.3 for comparison with directly computed solutions. Note in particular that the asymptotics becomesteeper as c x increases, accommodating thus for the mismatch of limiting values, − Z g = lim k y →∞ lim c x →∞ k x = lim c x →∞ lim k y →∞ k x = (cid:18) − Z g − (cid:19) − ;compare also the graphs in Figure 4.4. -3 -10 -5 00500100015002000
20 60 100 140 1800.99940.99950.99960.99970.99980.99991 -10 -10 Figure 4.3: Left: Selected k x vs k y for k y small, compared with numerically computed quadratic approximation (4.10); note the goodfit, albeit on increasingly small k y -ranges as c x decreases. Center: Quadratic coefficient k x, in (4.10) vs c x and comparison with best fit c = − . c = − . k x vs k y for k y large and sample values of c x , compared with (4.15); inset log-logplot of 1 − k x vs k y confirming the good cubic approximation for moderate values of c x . In the specific case of g ( v ) = 1 + κ sin( v ), the asymptotics described above coincide well with numerical computationsand predictions from the asymptotics give a good qualitative overall picture. Behavior for fixed k y . Fixing k y small, we discuss the curve k x ( c x ). By the integral identities above, k x (0) = 1and k ′ x (0) <
0, while k x ( ∞ ) <
1, monotonically increasing for k y < k ∗ y and monotonically decreasing for k y > k ∗ y , c x ≫
1. The asymptotics are therefore compatible with globally monotonically decreasing k x ( c x ) for k y > k ∗ y andwith k x ( c x ) having a unique minimum for some finite c x ( k y ) for k y < k ∗ y . This simple behavior with unique minimumor simple monotonicity is indeed what we observe numerically. Behavior for fixed c x . From the analysis above, we found k x ( ∞ ) = 1 and k x monotonically increasing for large k y (4.15). For k y = 0, the asymptotics and numerical analysis in [14] predict 1 − κ < k x (0) <
1. The asymptoticswith numerical evaluations of the relevant integrals predict that k x is monotonically increasing for k y ∼
0, as well.Curves k x ( k y ) computed for numerically are in fact monotonically increasing on k y ≥
0, albeit with a characteristictransition that we will discuss in the next section.
Behavior as k y → . One notices that the limit of curves k x ( c x ) as k y → k y = 0, theresults in [14] show a monotone curve k x = 1 − κ + O ( √ c x ), and k x ∈ [1 − κ, κ ] for c x = 0. For k y > k x ( c x )are non-monotone and appear to converge to this limiting set ( c x , k y ) ∈ × [1 − κ, κ ] ∪ { ( c x , k x ( c x )) , c x > } .12 ummary. Rephrasing our findings in terms of strain, measured through the deviation of k x from the equilibriumstrain k x = 1, induced on stripes through forced growth at rate c x and imposed angle determined by k y , we cansummarize our findings as follows.(i) for small angles, k y ∼
0, slow growth creates the largest residual strain in the stripes. For zero angles, k y = 0,the strain decreases with increased growth rate, but for small angles the residual strain first increases with c x before faster growth reduces strain;(ii) for fixed growth rate, residual strain decreases with increasing angles;(iii) for larger angles, strain increases with growth rate.The induced strain at k y = 0 can be understood as a non-adiabatic effect, proportional to κ which measures thenon-adiabaticity, that is, the size of terms that do not commute with the phase averaging symmetry ϕ ϕ + const .Stripes are stretched maximally for small speeds, repeated stripe nucleation helps release stress with increased growthrate as described in [14]. For small angles, an effect similar to zero angle can be observed, with the caveat that forvery small speeds, the gliding of a localized boundary defect along the growth interface can mediate the growthprocess with little residual stress. Increasing the rate of growth increases the glide speed of the defect and therebyresidual strain. Yet stronger growth leads to a phase transition in the nature of the boundary defect that leads todelocalization and decreased strain.Increasing the angle through k y reduces the non-adiabaticity, up to the point where stripes perpendicular to theboundary can grow without deformation at the interface, k y = 1, not creating any strain. Figure 4.4 shows the surface k x ( k y , c x ) from different angles, exhibiting the singularities that occur in the compactification at the boundaries c x , k y ∈ { , ∞} . Figure 4.4: Surface k x as a function of k y and c x . Plots use k y / (1 + k y ) and c x / c x ) as coordinates to include the limits c x = ∞ and k y = ∞ at 1; see also mod space all.mp4 in the supplementary materials. The strain in a large region of parameter space is simply monotone and fairly simple asymptotics explain the behavior.The most intriguing, non-monotone dynamics occur in a vicinity of c x = k y = 0. In this regime, profiles ϕ converge tostep-like functions in ζ ; see Figure 1.3. An inner expansion of the layer-type solution reveals an interesting transitionthat sheds light on the asymptotics in this region.We scale in (2.2)–(2.5) for an inner expansion at the heteroclinic k y = ˜ k y ε , c x = ε and ∂ ζ = ε∂ z , and obtain,expanding the Fourier symbol D + , at leading order D ψ = g ( ψ ) − k x , y ∈ R , ψ ( −∞ ) + 2 π = ψ (+ ∞ ) = ψ ∗ , D = q − ˜ k y ∂ zz + k x ∂ z , (5.1)13here D now is defined as a Fourier multiplier for functions on the real line rather than periodic functions. Thisequation does have a local interpretation as a traveling-wave solution ψ = ψ (˜ k y y + k x t, x ) to the heat equation withnonlinear boundary flux, ψ t = ∆ ψ, x < , y ∈ R , ψ x = g ( ψ ) − k x , x = 0 , y ∈ R . Such traveling waves have been studied in [8], establishing in particular existence and monotonicity properties forsolutions ψ ( y − ct ), with c = c ( k x ) for | k x − | < κ when g ( ψ ) = 1 + κ sin( ψ ). Rescaling y = z/k y shows that thesetraveling solutions give solutions to (5.1) whenever k y = k x c ( k x ) . (5.2)Moreover, monotonicity of c in k x from [8] implies that k y is monotonically increasing as a function of k x withminimum k ∗ y , such that we can rewrite (5.2) as k x = k f x ( k y ) , for k y > k ∗ y . (5.3)For 0 < k y < k ∗ y , we conjecture the existence of heteroclinic solutions with k x = min g ( ϕ ), asymptotic to argmin g ( ϕ )and argmin g ( ϕ ) + 2 π . In particular, the selected k x is constant at leading order.Below, we provide numerical evidence for our predictions. Computing heteroclinic orbits in (5.1) . We focus on the specific case g ( ψ ) = 1 + κ sin ψ . In order to solve(5.1), we rely on Fourier transform. We therefore write ψ = ψ s + ˜ ψ with ψ s ( z ) = ψ ∗ + 2 arctan( z ), where g ( ψ ∗ ) = k x , g ′ ( ψ ∗ ) >
0. The choice of the arctan( z ) is motivated by the fact that the action of the integral operator is explicit, R ( z ; k x , ˜ k y ) := D ψ s ( z ) = 2 √ π ˜ k y z Re (1 + i)U − , , k x ( − i + x )˜ k y !! , where U is the confluent hypergeometric Kummer-U function. We then solve D ( k x , ˜ k y ) ˜ ψ + R ( k x , ˜ k y ) − g ( ψ s + ˜ ψ ) + k x = 0 , with periodic boundary conditions on a large domain | z | ≤ L together with a phase condition R ˜ ψ ( z )e − z d z = 0 andwith k x as a Lagrange multiplier using a Newton method and secant continuation in k y . The spectral discretizationgives accuracy of 10 − for moderate effective discretization sizes of 0 .
1. Solutions decay however only weakly with z − / , z → −∞ , and z − / for z → + ∞ . We found accuracy of 10 − for domain sizes L ∼ using N = 2 ∼ Fourier modes. The code was implemented in matlab and ran on an Nvidia GV100 graphics card allowing forfast evaluation of the large discrete Fourier transforms. The Kummer-U function was evaluated and tabulated in mathematica and interpolated in matlab , since direct evaluation in matlab is slow.Results from the computation of heteroclinic orbits are shown in Figure 5.2, left upper panel, showing a characteristictransition from increasing values k x (˜ k y ) for moderate ˜ k y to constant k x for small ˜ k y . At the transition value, theheteroclinic orbit delocalizes, the amplitude of ψ y decreases. In the limit k y → ∞ , we find the “Hamiltonian” picture,with k x = 1.The computed values of k x compare well with the selected values in the selection problem periodic in y , as shownin Figure 5.1. Selected wavenumbers k x as function of the scaled wavenumber k y /c x , computed for fixed values of c x ≪ k y →
0, converge to the limiting curve given by the heteroclinic orbit.The nonlocal problem is related to the Weertman equation that is used to describe the glide motion of dislocations;see[18] and references therein. In fact, the nonlocal Weertman equation can be obtained by replacing our nonlinearfluxes by a dynamic (Wentzel) boundary conditions, ϕ t = ∆ ϕ, x < ϕ t = − ϕ x + g ( ϕ ) , x = 0 . -3 -10 -8 -6-5-4-3 Figure 5.1: Left: selected k x in (2.7) plotted against k y /c x = ˜ k y , for c x decreasing geometrically by factors 8 / k y from the heteroclinic continuation for comparison. Right: section through the diagram on the left, plotting k x as a function of c x for fixed ˜ k y , showing in particular that the values are almost independent of ˜ k y <
2, below the heteroclinic bifurcation. For such smallvalues, k x ∼ . . √ c x in good agreement with [14]. Our numerical methods in fact resemble the approach taken in [19], although pseudo-differential operators are moredifficult in our case and the emphasis in [18] is on the time-dependent initial-value problem. We conclude this analysiswith a heuristic explanation of the transition from a sharply localized defect selecting strains k x to a delocalizedheteroclinic selecting minimal values of k x , through analogy to a local differential equation. Comparison with local heteroclinic bifurcations.
A qualitatively equivalent picture emerges when the non-local pseudo-differential operator D is replaced by a local operator D loc = − ˜ k y ∂ zz + k x ∂ z . In this case, elementaryphase plane analysis establishes the existence of heteroclinic orbits to D loc ψ = 1 + κ sin( ψ ) − k x . Rescaling k y ∂ z = ∂ y and writing c = k x /k y , we find the traveling-wave equation to the (asymmetric) parabolic Sine-Gordon equation, u yy + cu y = 1 + κ sin( u ) − k x . For k x = 1 we have c = 0 and a heteroclinic between u = 0 and u = 2 π . The heteroclinic is transverse in theparameter c and we can in fact continue the heteroclinic with c = c ( k x ) monotonically increasing as k x is decreasing,until k x = 1 − κ . For c ≫
1, we find at leading order, after a a reduction to a slow manifold, cu y = 1 + κ sin( u ) − k x , which possesses heteroclinic orbits for k z = 1 − κ , connecting the saddle-node equilibria u = π/ π . Theseheteroclinics between saddle-node equilibria are robust up to a heteroclinic codimension-two bifurcation [11, 4]. Theassociated phase-portraits in the u − u x –plane are shown in Figure 5.3 and can be easily confirmed using elementaryphase-plane analysis and monotonicity in c . Returning to the motivation by striped patterns, we now study the formation of striped patterns in a directionallyquenched Swift-Hohenberg equation. The phase-diffusion approximation with nonlinear boundary fluxes given by thestrain-displacement relation was shown to be a correct approximation in the case c x = 0 in [31], for y -independentpatterns. Considering patterns in two spatial dimensions, one notices that patterns selected for c x ≪ k y ≪ k < Figure 5.2: Left, top: nonlocal; Left, bottom: local. Right: sketch of phase portraits with heteroclinics in local case. along stripes and higher-order corrections as in the Cross-Newell equation are necessary to fully capture dynamics;see for instance [28].We therefore focus on the quenched anisotropic Swift-Hohenberg equation, u t = − (1 + ∆ x,y ) u + β∂ yy u + µu − u , (6.1)used in [7, 20, 29, 30, 21, 26, 16, 23, 10] to describe nematic liquid crystals, electroconvection, ion bombardment,surface catalysis, or vegetation patterns; see also [17] for an analysis of dislocations in this model. For β >
0, theanisotropic term suppresses the zig-zag instability in stripes with wavenumbers k .
1. for sufficiently large β allwavenumbers within the strain-displacement relation, k ∈ ( k min , k max ), with k max = max g ( φ ), are stabilized. In thefollowing, we first derive a phase-diffusion approximation and nonlinear fluxes in the form studied in this paper fromthe anisotropic Swift-Hohenberg equation, and then describe a numerical approach to computing striped patternscreated in directional quenching, with the goal of comparing the numerical results to the quantitative predictionsfrom the phase-diffusion approximation. Throughout, we focus on the regime 0 < c x , k y ≪ µ = − µ tanh(( x − c x t ) /δ ) with δ = 0 . Derivation of phase diffusion in anisotropic Swift-Hohenberg.
Focusing on nearly parallel stripes withconstant parameter µ , we use the parabolic scaling µ = ǫ , x = ǫ ˜ x, y = ǫ ˜ y, t = ǫ ˜ t , and substitute the ansatz u ( x, y, t ) = εA (cid:0) ˜ x, ˜ y, ˜ t (cid:1) e i x + c . c . into (6.1) to obtain, at leading order, an anisotropic Ginzburg-Landau equation A ˜ t = 4 A ˜ x ˜ x + βA ˜ y ˜ y + A − A | A | . (6.2)Introducing polar coordinates A = R e i ˜ φ and expanding near R = 1 / √
3, ˜ φ = 0, one finds an exponentially dampedequation for R and an anisotropic diffusion equation for ˜ φ ,˜ φ ˜ t = 4 ˜ φ ˜ x ˜ x + β ˜ φ ˜ y ˜ y . (6.3)Note that this equation is again invariant under the parabolic scaling such that we may consider (6.3) in the originalcoordinates t, x, y to describe patterns in (6.1). 16 igure 5.3: Profiles of derivatives ∂ y ϕ (top) and ∂ x ψ in unscaled variables y for values of ˜ k y = k y /c x , c x = 10 − , passing the heteroclinicbifurcation. Profiles are roughly constant for large ˜ k y (left inset) but rapidly delocalize past the heteroclinic transition with long tails tothe left of the peak; amplitude of profiles rapidly decreases past heteroclinic bifurcation (right inset). Normal derivatives also delocalizebut always peak at minimal and maximal strain. We next turn to the effect of the spatial quenching. At the order of the Ginzburg-Landau equation, one doesnot capture the non-adiabatic effects of the parameter jump. We use the expression for the strain-displacementrelation from [31] for the strain-displacement relation in the one-dimensional case, unaffected by the anisotropicterm, ˜ φ x = g SH ( ˜ φ ) := 1 + µ sin 2 ˜ φ + O ( µ / ). The symmetry ˜ φ ˜ φ + π is present at higher orders, as well, andcaused by the u
7→ − u symmetry in the nonlinearity and the ensuing symmetry u per ( ξ )
7→ − u per ( ξ + π ) of periodicpatterns. We use the same boundary condition for two-dimensional patterns, neglecting in particular dependence of g SH on ˜ φ y , and also neglect dependence on c x , which gives the two-dimensional system˜ φ t = 4 ˜ φ xx + β ˜ φ yy + ˜ c x ˜ φ, x < , y ∈ R , ˜ φ x = g SH ( ˜ φ ) , x = 0 , y ∈ R . (6.4)With the additional scaling φ = 2 ˜ φ, x = ˜ x, y = ˜ y, c x = 8˜ c x , t = 16˜ t, we then obtain the phase-diffusionequation (1.3) with strain-displacement relation φ x = g SH ( φ/
2) at x = 0. We remark that by setting κ = µ / g SH agrees to leading order with the relation φ x = g ( φ ) employed in previous sections. Through these scalings, wecan compare the heteroclinic prediction of Section 5 with moduli curves of quenched patterned solutions u (˜ x, ˜ y, t ) = u ( k x (˜ x − ˜ c x t ) , k y (˜ y − c y ˜ t )) of the full equation (6.1). In our comparisons below, we use a value for κ slightlydifferent from µ /
16, computed directly from the one-dimensional Swift-Hohenberg equation as described in [25, 31],accounting for both error terms O ( µ / ) and corrections due to the fact that we use a smoothed version of the stepfunction for the spatially dependent parameter µ . Oblique stripe formation in the full Swift-Hohenberg equation.
In the Swift-Hohenberg equation, formationof striped patterns is described by traveling-wave solutions [15, 2] with speed vector ( c x , c y ), again requiring c y = k x ˜ c x /k y , 0 = − (1 + k x ∂ ξ + k y ∂ ζ ) u + µu − u + ˜ c x k x ( ∂ ξ + ∂ ζ ) u, ξ < , ζ ∈ R , (6.5)0 = u ( ξ, ζ + 2 π ) − u ( ξ, ζ ) , ξ ≤ , ζ ∈ R , (6.6)0 = lim ξ →∞ u ( ξ, ζ ) , ξ →−∞ | u ( ξ, ζ ) − u per (cid:16) ξ + ζ ; q k x + k y (cid:17) | , ζ ∈ R . (6.7)We numerically solve (6.5) - (6.7) using a farfield-core approach similar to [22, 2], which decomposes u = w + χu per (cid:16)q k x + k y (cid:17) , where w is localized near the quenching interface, and χ is a cutoff function supported in the17 -farfield. Here, we solve for w and k x with parameter k y , using a spectral discretization in both ξ and ζ so thatfunctions can be evaluated with the fast Fourier transform. Each Newton step of the pseudo-arclength continuationalgorithm was once again performed using gmres to solve the associated linear problem. The nonlinear system wasconjugated with exponentially localized weights and pre-conditioned with the principal symbol of the linear equation.Discretization and domain size were controlled adaptively ensuring both small tails at the end of the (periodic)domain and small amplitudes in highest Fourier modes. Typical domain sizes near the origin were x ∈ ( − , × ξ, y ). Code was again implemented in matlab with computations carried outusing an Nvidia GV100 GPU. Further details of this numerical approach are left for a companion work. For valuesof ˜ c x and k y smaller than the ones shown, gmres would usually not converge due to constraints on the number ofinner iterations caused by limited memory. Comparisons between phase-diffusion and Swift-Hohenberg.
Figure 6.1 gives slices of the moduli space for(6.1) with ˜ c x fixed and shows that the surface is a graph k x = k x ( k y , ˜ c x ) for ( k y , ˜ c x ) ∼
0. Curves, which are plottedover the scaled wavenumber ˜ k y = k y / ˜ c x , show good agreement with the heteroclinic asymptotics of Section 5, witha transition around k y / ˜ c x ∼ Figure 6.1: Wavenumber selection curves for anisotropic Swift-Hohenberg (6.5)–(6.7) for k y / ˜ c x ∼ c x fixed. Left: comparison fora range of ˜ c x values with the heteroclinic curve (black) of §
5; here, β = 1 and µ = 3 / κ = µ /
16 = 3 /
64. The heterocliniccurve (black) is obtained using numerically derived strain-displacement relation to account for higher-order corrections in µ . Right: plotof selected wavenumber k for k y / ˜ c x ∼ β values with µ = 3 / µ values with β = 1 fixed(dot-dashed), ˜ c x = 0 . Varying the anisotropy coefficient β and the parameter µ , we also show how this phase transition depends onsystem parameters. As expected, the strength of non-adiabatic effects increases with µ as small- µ averaging is lesseffective, and the strain 1 − k on the stripes created at small k y increases, roughly proportional to µ as predictedby the amplitude µ /
16 of the strain-displacement relation. The location of the transition appears to be roughlyindependent of µ , in agreement with our derivation above. Varying the strength of anisotropy does affect thetransition. Stronger anisotropy narrows the plateau where delocalized defects determine wavenumber selection. Veryweak and in particular vanishing anisotropy lead to non-monotone dependence of k on k y which is beyond the scopeof this paper. We investigated directional growth of striped phases in the absence of instabilities and for weakly oblique orientationof stripes relative to the boundary. In a reduced phase-diffusion approximation, we established existence of simple,18 igure 6.2: Plots of solutions of (6.5)-(6.7) near quenching interface in original coordinates for ˜ c x = 10 − fixed for a range of ˜ k y values:˜ k y = 118 . ... (top left), ˜ k y = 25 . ... (top right ). Bottom row illustrates delocalization of dislocation defect both in x and y for small˜ k y , with a zoom-in near a defect for ˜ k y = 25 . ... (left) ˜ k y = 4 . ... (right). Note that the odd symmetry in Swift-Hohenberg creates twoantisymmetric dislocation-type defects, a covering symmetry visible also in the phase-diffusion approximation through the dependence ofthe strain-displacement relation on 2 ˜ φ , only. resonant growth mechanisms and derived universal asymptotics in limiting regimes. Our results compare well withcomputations in a Swift-Hohenberg equation where instabilities are suppressed by weak anisotropy.Many of our results can be rephrased in coarse terms. For parallel stripes, we had earlier found that very small speedscause maximal strain, given by the minimum of the strain-dispersion relation, which decreases up to a dynamicallyaveraged (harmonic average) strain for large speeds. Zero speeds and growth at larger angles yield zero strain,with selected wavenumber given by the (energy-minimizing) average of the strain-displacement relation. At smallangles, k y ∼
0, the growth process is mediated by the emergence of a point defect at the boundary, which undergoes adelocalization bifurcation at a critical value, similar in character to the codimension-two bifurcation from a hyperbolichomoclinic orbit to a saddle-node homoclinic orbit. The growth process is described well by a glide motion of thedefect along the boundary of the patterned region, adding one stripe once the defect has moved by one period alongthe boundary. In our asymptotics, we identify the glide motion in the absence of growth, c x = 0, when a non-equilibrium strain k = − R g is imposed in the far field: the nonequilibrium strain k drives the defect at a finite speed c y ( k ), c ′ y = 0. Then, for a growth process with given speed c x and angle k y , the selected wavenumber k adjusts suchthat the induced glide speed c y ( k ) corresponds to compatible defect motion by one y -period 2 π/k y while one stripe isgrown across the interface, in time 2 π/ ( c x k x ). The effective wavenumber used in the scaling, k y /c x = k x /c y ∼ /c y ,is at leading order simply the inverse glide speed. From this perspective, the ˜ k y -dependent contribution to the strainstems from drag in the glide motion of the defect. The c x -dependence can be understood as in [14] as an interactionbetween dislocation over the finite distance 2 π/k y , effectively leading to an effective deceleration of the glide motionand reduced strain.Looking forward, we hope that this glimpse into the role of point defects in growth of crystalline phases can beextended, including for instance the effect of zigzag instabilities associated with wrinkling. More mathematically, ofthe many phenomena described here, it would be interesting to analyze the heteroclinic bifurcation at the origin,finding in particular better asymptotics near the critical value of ˜ k y . One may also hope to better understandsome of the asymptotic expansions derived here, adding mathematical rigor, or relating them more directly to ourunderstanding of dislocations, their farfield, and interaction properties.19 eferences [1] S. Akamatsu, S. Bottin-Rousseau, and G. Faivre. Experimental evidence for a zigzag bifurcation in bulk lamellareutectic growth. Phys. Rev. Lett. , 93:175701, Oct 2004.[2] M. Avery, R. Goh, O. Goodloe, A. Milewski, and A. Scheel. Growing stripes, with and without wrinkles.
SIAMJ. Appl. Dyn. Syst. , 18(2):1078–1117, 2019.[3] B. Bakker and A. Scheel. Spatial Hamiltonian identities for nonlocally coupled systems.
Forum Math. Sigma ,6:Paper No. e22, 55, 2018.[4] J. Bellay and A. Scheel. Coherent structures near the boundary between excitable and oscillatory media.
Dyn.Syst. , 25(1):111–132, 2010.[5] E. Ben-Jacob, I. Cohen, and H. Levine. Cooperative self-organization of microorganisms.
Advances in Physics ,49(4):395–554, 2000.[6] E. Bodenschatz, W. Pesch, and G. Ahlers. Recent developments in Rayleigh-B´enard convection.
Annual Reviewof Fluid Mechanics , 32(1):709–778, 2000.[7] D. Boyer. Numerical study of domain coarsening in anisotropic stripe patterns.
Phys. Rev. E , 69:066111, Jun2004.[8] X. Cabr´e, N. C´onsul, and J. V. Mand´e. Traveling wave solutions in a half-space for boundary reactions.
Anal.PDE , 8(2):333–364, 2015.[9] A. S.-H. Chen and S. W. Morris. On the origin and evolution of icicle ripples.
New Journal of Physics ,15(10):103012, 2013.[10] Y. Cheng, M. Stieglitz, G. Turk, and V. Engel. Effects of anisotropy on pattern formation in wetland ecosystems.
Geophysical Research Letters , 38(4), 2011.[11] S.-N. Chow and X.-B. Lin. Bifurcation of a homoclinic orbit with a saddle-node equilibrium.
Differential IntegralEquations , 3(3):435–466, 1990.[12] A. Doelman, B. Sandstede, A. Scheel, and G. Schneider. The dynamics of modulated wave trains.
Mem. Amer.Math. Soc. , 199(934):viii+105, 2009.[13] D. Double. Imperfections in lamellar eutectic crystals.
Materials Science and Engineering , 11(6):325 – 335,1973.[14] R. Goh, R. Beekie, D. Matthias, J. Nunley, and A. Scheel. Universal wave-number selection laws in apicalgrowth.
Phys. Rev. E , 94:022219, Aug 2016.[15] R. Goh and A. Scheel. Pattern-forming fronts in a Swift-Hohenberg equation with directional quenching—paralleland oblique stripes.
J. Lond. Math. Soc. (2) , 98(1):104–128, 2018.[16] D. Handwerk, G. Dangelmayr, I. Oprea, and P. D. Shipman. Phase-diffusion equations for the anisotropiccomplex ginzburg-landau equation, 2020.[17] M. Haragus and A. Scheel. Dislocations in an anisotropic Swift-Hohenberg equation.
Comm. Math. Phys. ,315(2):311–335, 2012.[18] M. Josien. Some mathematical properties of the Weertman equation.
Commun. Math. Sci. , 16(6):1729–1748,2018. 2019] M. Josien, Y.-P. Pellegrini, F. Legoll, and C. Le Bris. Fourier-based numerical approximation of the Weertmanequation for moving dislocations.
Internat. J. Numer. Methods Engrg. , 113(12):1827–1850, 2018.[20] C. Kamaga, F. Ibrahim, and M. Dennin. Dislocation dynamics in an anisotropic stripe pattern.
Phys. Rev. E ,69:066213, Jun 2004.[21] L. Kramer and W. Pesch. Convection instabilities in nematic liquid crystals.
Annual Review of Fluid Mechanics ,27(1):515–539, 1995.[22] D. J. B. Lloyd and A. Scheel. Continuation and bifurcation of grain boundaries in the Swift-Hohenberg equation.
SIAM J. Appl. Dyn. Syst. , 16(1):252–293, 2017.[23] A. S. Mikhailov and K. Showalter. Control of waves, patterns and turbulence in chemical systems.
PhysicsReports , 425(2):79 – 194, 2006.[24] S. A. Mollick, D. Ghose, P. D. Shipman, and R. Mark Bradley. Anomalous patterns and nearly defect-freeripples produced by bombarding silicon and germanium with a beam of gold ions.
Applied Physics Letters ,104(4):043103, 2014.[25] D. Morrissey and A. Scheel. Characterizing the effect of boundary conditions on striped phases.
SIAM J. Appl.Dyn. Syst. , 14(3):1387–1417, 2015.[26] J. Mu˜noz-Garc´ıa, L. V´azquez, M. Castro, R. Gago, A. Redondo-Cubero, A. Moreno-Barrado, and R. Cuerno.Self-organized nanopatterning of silicon surfaces by ion beam sputtering.
Materials Science and Engineering:R: Reports , 86:1 – 44, 2014.[27] A. Pauthier and A. Scheel. Advection-diffusion dynamics with nonlinear boundary flux as a model for crystalgrowth.
Math. Nachr. , 293(8):1565–1590, 2020.[28] L. M. Pismen.
Patterns and interfaces in dissipative dynamics . Springer Series in Synergetics. Springer-Verlag,Berlin, 2006. With a foreword by Y. Pomeau.[29] H. Qian and G. F. Mazenko. Growth of order in an anisotropic swift-hohenberg model.
Phys. Rev. E , 73:036117,Mar 2006.[30] A. Roxin and H. Riecke. Rotating convection in an anisotropic system.
Phys. Rev. E , 65:046219, Apr 2002.[31] A. Scheel and J. Weinburd. Wavenumber selection via spatial parameter jump.
Philos. Trans. Roy. Soc. A ,376(2117):20170191, 20, 2018.[32] R. Sheth, L. Marcon, M. F. Bastida, M. Junco, L. Quintana, R. Dahn, M. Kmita, J. Sharpe, and M. A.Ros. Hox genes regulate digit patterning by controlling the wavelength of a Turing-type mechanism.
Science ,338(6113):1476–1480, 2012.[33] A. Stegner and J. E. Wesfreid. Dynamical evolution of sand ripples under water.
Phys. Rev. E , 60:R3487–R3490,Oct 1999.[34] S. Thomas, I. Lagzi, F. Moln´ar Jr, and Z. R´acz. Probability of the emergence of helical precipitation patternsin the wake of reaction-diffusion fronts.
Physical Review Letters , 110(7):078303, 2013.[35] Wilczek, M., Tewes, W. B.H., Gurevich, S. V., K¨opf, M. H., Chi, L. F., and Thiele, U. Modelling patternformation in dip-coating experiments.